OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
phasefieldelement.C
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  *
14  * Copyright (C) 1993 - 2013 Borek Patzak
15  *
16  *
17  *
18  * Czech Technical University, Faculty of Civil Engineering,
19  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
20  *
21  * This library is free software; you can redistribute it and/or
22  * modify it under the terms of the GNU Lesser General Public
23  * License as published by the Free Software Foundation; either
24  * version 2.1 of the License, or (at your option) any later version.
25  *
26  * This program is distributed in the hope that it will be useful,
27  * but WITHOUT ANY WARRANTY; without even the implied warranty of
28  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
29  * Lesser General Public License for more details.
30  *
31  * You should have received a copy of the GNU Lesser General Public
32  * License along with this library; if not, write to the Free Software
33  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
34  */
35 
36 #include "../sm/Elements/phasefieldelement.h"
37 #include "../sm/CrossSections/structuralcrosssection.h"
38 #include "../sm/Materials/structuralms.h"
39 #include "node.h"
40 #include "material.h"
41 #include "gausspoint.h"
42 #include "gaussintegrationrule.h"
43 #include "domain.h"
44 #include "mathfem.h"
45 #include "timestep.h"
46 #include <cstdio>
47 
48 
49 namespace oofem {
50 
52 {
54  internalLength = 6.0;
55  criticalEnergy = 1.0e0;
56  relaxationTime = 1.0;
57 }
58 
59 void
61 {
62  // Routine to extract compute the location array an element given an dofid array.
63  answer.clear();
64  NLStructuralElement *el = this->giveElement();
65  int k = 0;
66  for ( int i = 1; i <= el->giveNumberOfDofManagers(); i++ ) {
67  DofManager *dMan = el->giveDofManager( i );
68  for ( int j = 1; j <= dofIdArray.giveSize( ); j++ ) {
69 
70  if ( dMan->hasDofID( (DofIDItem) dofIdArray.at( j ) ) ) {
71  Dof *d = dMan->giveDofWithID( dofIdArray.at( j ) );
72  answer.followedBy( k + d->giveNumber( ) );
73  }
74  }
75  k += dMan->giveNumberOfDofs( );
76  }
77 }
78 
79 void
82 {
83  IntArray dofIdArray;
84  this->giveDofManDofIDMask_u(dofIdArray);
85  Element :: computeVectorOf(dofIdArray, valueMode, stepN, answer);
86 }
87 
88 void
90 {
91  IntArray dofIdArray;
92  this->giveDofManDofIDMask_d(dofIdArray);
93  Element :: computeVectorOf(dofIdArray, valueMode, stepN, answer);
94 }
95 
96 void
97 PhaseFieldElement :: giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
98 {
99  // Computes the internal forces corresponding to the two fields u & d
100  IntArray IdMask_u, IdMask_d;
101  this->giveDofManDofIDMask_u( IdMask_u );
102  this->giveDofManDofIDMask_d( IdMask_d );
103  this->computeLocationArrayOfDofIDs( IdMask_u, loc_u );
104  this->computeLocationArrayOfDofIDs( IdMask_d, loc_d );
105 
106  int ndofs = this->giveElement()->computeNumberOfDofs();
107  answer.resize(ndofs);
108  answer.zero();
109 
110  FloatArray answer_u(0);
111  FloatArray answer_d(0);
112 
113  this->giveInternalForcesVector_u(answer_u, tStep, useUpdatedGpRecord);
114  this->giveInternalForcesVector_d(answer_d, tStep, useUpdatedGpRecord);
115  answer.assemble(answer_u, loc_u);
116  answer.assemble(answer_d, loc_d);
117 }
118 
119 void
121 {
122  // computes int_V ( B_u^t * BSgima_u ) * dV
123  FloatArray NStress, BStress, vGenStress, NS;
124  FloatMatrix N, B;
125  NLStructuralElement *el = this->giveElement();
126 
127  answer.clear();
128  for ( auto &gp: el->giveIntegrationRule(0) ) {
129  double dV = el->computeVolumeAround(gp);
130 
131  // compute generalized stress measure
132  el->computeBmatrixAt(gp, B);
133  this->computeBStress_u(BStress, gp, tStep, useUpdatedGpRecord);
134  answer.plusProduct(B, BStress, dV);
135  }
136 }
137 
138 void
140 {
141  // Computes int_V ( N^t *Nstress_d + B^t * g_c*l*grad(d) ) dV
142  FloatArray NStress, BStress, NS, a_d, grad_d;
143  FloatMatrix N, B;
144  NLStructuralElement *el = this->giveElement();
145  this->computeDamageUnknowns( a_d, VM_Total, tStep );
146 
147  for ( auto &gp: el->giveIntegrationRule(0) ) {
148  double dV = el->computeVolumeAround(gp);
149 
150  // compute generalized stress measures
151  this->computeNd_matrixAt( *gp->giveNaturalCoordinates(), N);
152  computeNStress_d(NStress, gp, tStep, useUpdatedGpRecord);
153  NS.beTProductOf(N, NStress);
154  answer.add(dV, NS);
155 
156  this->computeBd_matrixAt(gp, B);
157  grad_d.beProductOf(B, a_d);
158  double l = this->giveInternalLength();
159  double g_c = this->giveCriticalEnergy();
160  BStress = grad_d * l * g_c;
161  answer.plusProduct(B, BStress, dV);
162  }
163 }
164 
165 void
166 PhaseFieldElement :: computeBStress_u(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, int useUpdatedGpRecord)
167 {
168  // computes G(d)*sig(u)
169  NLStructuralElement *el = this->giveElement();
170  FloatArray strain, a_u;
171  FloatMatrix B_u;
172 
173  el->computeBmatrixAt(gp, B_u);
174  this->computeDisplacementUnknowns(a_u, VM_Total, tStep);
175  strain.beProductOf(B_u, a_u);
176  el->computeStressVector(answer, strain, gp, tStep);
177  answer.times( this->computeG(gp, VM_Total, tStep) );
178 
179 }
180 
181 double
183 {
184  StructuralMaterialStatus *matStat = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() );
185  FloatArray strain, stress;
186  stress = matStat->giveTempStressVector();
187  strain = matStat->giveTempStrainVector();
188  return 0.5 * stress.dotProduct( strain );
189 }
190 
191 void
192 PhaseFieldElement :: computeNStress_d(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, int useUpdatedGpRecord)
193 {
194  // computes (t*/Delta t) *(d - d_old) + g_c/l*d + G'*Psibar
195 
196  //PhaseFieldCrossSection *cs = static_cast...
197  double Delta_t = tStep->giveTimeIncrement();
198  double t_star = this->giveRelaxationTime();
199  double d = computeDamageAt( gp, VM_Total, tStep );
200  double Delta_d = computeDamageAt( gp, VM_Incremental, tStep );
201 
202  double l = this->giveInternalLength();
203  double g_c = this->giveCriticalEnergy();
204  double Gprim = this->computeGPrim(gp, VM_Total, tStep);
205  double Psibar = this->computeFreeEnergy( gp, tStep );
206  answer.resize( 1 );
207  answer.at( 1 ) = t_star / Delta_t * Delta_d + g_c / l *d + Gprim * Psibar;
208 }
209 
210 double
212 {
213  // d = N_d * a_d
214  NLStructuralElement *el = this->giveElement();
215  FloatArray dVec;
216  computeDamageUnknowns(dVec, valueMode, stepN);
217  FloatArray Nvec;
219  return Nvec.dotProduct(dVec);
220 }
221 
222 double
224 {
225  // computes Dg/Dd = (1-d)^2 + r0
226  double d = this->computeDamageAt(gp, valueMode, stepN);
227  double r0 = 1.0e-10;
228  return (1.0 - d) * (1.0 - d) + r0;
229 }
230 
231 double
233 {
234  // compute -2*(1-d)
235  double d = this->computeDamageAt(gp, valueMode, stepN);
236  return -2.0 * (1.0 - d);
237 }
238 
239 void
241 {
242  NLStructuralElement *el = this->giveElement();
243  FloatArray Nvec;
244  el->giveInterpolation( )->evalN( Nvec, lCoords, FEIElementGeometryWrapper( el ) );
245  N.resize(1, Nvec.giveSize());
246  N.beNMatrixOf(Nvec,1);
247 
248 }
249 
250 void
252 {
253  // Returns the [numSpaceDim x nDofs] gradient matrix {B_d} of the receiver,
254  // evaluated at gp.
255 
256  NLStructuralElement *el = this->giveElement();
257  FloatMatrix dNdx;
259  answer.beTranspositionOf( dNdx );
260 }
261 
262 void
264 {
265  //set displacement and nonlocal location array
267  IntArray IdMask_u, IdMask_d;
268  this->giveDofManDofIDMask_u( IdMask_u );
269  this->giveDofManDofIDMask_d( IdMask_d );
270  this->computeLocationArrayOfDofIDs( IdMask_u, loc_u );
271  this->computeLocationArrayOfDofIDs( IdMask_d, loc_d );
272 
273  int nDofs = this->giveElement()->computeNumberOfDofs();
274  answer.resize( nDofs, nDofs );
275  answer.zero();
276 
277  FloatMatrix answer1, answer2, answer3, answer4;
278  this->computeStiffnessMatrix_uu(answer1, rMode, tStep);
279  this->computeStiffnessMatrix_ud(answer2, rMode, tStep);
280  //this->computeStiffnessMatrix_du(answer3, rMode, tStep); //symmetric
281  answer3.beTranspositionOf( answer2 );
282  this->computeStiffnessMatrix_dd(answer4, rMode, tStep);
283 
284  answer.assemble( answer1, loc_u, loc_u );
285  answer.assemble( answer2, loc_u, loc_d );
286  answer.assemble( answer3, loc_d, loc_u );
287  answer.assemble( answer4, loc_d, loc_d );
288 }
289 
290 void
292 {
293  // This is the regular stiffness matrix times G
294  FloatMatrix B, DB, N, D_B;
295  NLStructuralElement *el = this->giveElement();
296  StructuralCrossSection *cs = dynamic_cast<StructuralCrossSection* > (el->giveCrossSection() );
297 
298  bool matStiffSymmFlag = cs->isCharacteristicMtrxSymmetric(rMode);
299  answer.clear();
300 
301  for ( auto gp: el->giveIntegrationRule(0) ) {
302  double dV = el->computeVolumeAround(gp);
303  // compute int_V ( B^t * D_B * B )dV
304  el->computeBmatrixAt(gp, B );
305  cs->giveCharMaterialStiffnessMatrix(D_B, rMode, gp, tStep);
306  D_B.times( computeG(gp, VM_Total, tStep) );
307  DB.beProductOf(D_B, B);
308 
309  if ( matStiffSymmFlag ) {
310  answer.plusProductSymmUpper(B, DB, dV);
311  } else {
312  answer.plusProductUnsym(B, DB, dV);
313  }
314  }
315 
316  if ( matStiffSymmFlag ) {
317  answer.symmetrized();
318  }
319 }
320 
321 void
323 {
324  FloatMatrix B, DB, N_d, DN, S(3,1);
325  FloatArray stress;
326  NLStructuralElement *el = this->giveElement();
327  StructuralCrossSection *cs = dynamic_cast<StructuralCrossSection* > (el->giveCrossSection() );
328 
329  answer.clear();
330 
331  for ( auto &gp: el->giveIntegrationRule(0) ) {
332  StructuralMaterialStatus *matStat = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() );
333 
334  double dV = el->computeVolumeAround(gp);
335  // compute int_V ( B^t * D_B * B )dV
336 
337  this->computeNd_matrixAt(*gp->giveNaturalCoordinates(), N_d);
338 
339  // stress
340  FloatArray reducedStrain, a_u;
341  FloatMatrix B_u;
342  el->computeBmatrixAt( gp, B_u );
343  stress = matStat->giveTempStressVector();
344  stress.times( this->computeGPrim( gp, VM_Total, tStep ) );
345 
346  S.setColumn(stress,1);
347  DN.beProductOf(S, N_d);
348  answer.plusProductUnsym(B_u, DN, dV);
349  }
350 
351 }
352 
353 void
355 {
356  double Delta_t = tStep->giveTimeIncrement();
357  double t_star = this->giveRelaxationTime();
358  double l = this->giveInternalLength();
359  double g_c = this->giveCriticalEnergy();
360 
361  FloatMatrix B_d, N_d;
362  //StructuralCrossSection *cs = dynamic_cast<StructuralCrossSection* > (this->giveCrossSection() );
363  NLStructuralElement *el = this->giveElement();
364  answer.clear();
365 
366  for ( auto &gp: el->giveIntegrationRule(0) ) {
367  double dV = el->computeVolumeAround(gp);
368 
369  this->computeNd_matrixAt(*gp->giveNaturalCoordinates(), N_d);
370  this->computeBd_matrixAt(gp, B_d, 1, 3);
371 
372  double Gprim = this->computeGPrim(gp, VM_Total, tStep);
373  double psiBar = this->computeFreeEnergy( gp, tStep );
374  //double factorN = t_star/Delta_t + g_c/l + psiBar*Gbis;
375  double factorN = t_star / Delta_t + g_c / l + psiBar*(-2.0);
376  double factorB = g_c*l;
377 
378  answer.plusProductSymmUpper(N_d, N_d, factorN*dV);
379  answer.plusProductSymmUpper(B_d, B_d, factorB*dV);
380  }
381 
382  answer.symmetrized();
383 }
384 
385 
388 {
389  //IRResultType result; // Required by IR_GIVE_FIELD macro
390  //nlGeo = 0;
391 
392  return IRRT_OK;
393 }
394 
395 
396 #if 0
397 void
398 PhaseFieldElement :: computeStressVectorAndLocalCumulatedStrain(FloatArray &answer, double localCumulatedStrain, GaussPoint *gp, TimeStep *stepN)
399 //PhaseFieldElement :: computeStressVector(FloatArray &answer, double localCumulatedStrain, GaussPoint *gp, TimeStep *stepN)
400 {
401  NLStructuralElement *elem = this->giveNLStructuralElement();
402 
403  nlGeo = elem->giveGeometryMode();
404  StructuralCrossSection *cs = this->giveNLStructuralElement()->giveStructuralCrossSection();
405 
406  //this->computeDamage(alpha, gp, stepN);
407  if ( nlGeo == 0 ) {
408  FloatArray Epsilon;
409  this->computeLocalStrainVector(Epsilon, gp, stepN);
410  dpmat->giveRealStressVector(answer, gp, Epsilon, nlCumulatedStrain, stepN);
411  } else if ( nlGeo == 1 ) {
412  if ( elem->giveDomain()->giveEngngModel()->giveFormulation() == TL ) {
413  FloatArray vF;
414  this->computeDeformationGradientVector(vF, gp, stepN);
415  dpmat->giveFirstPKStressVector(answer, gp, vF, nlCumulatedStrain, stepN);
416  }
417 
418  }
419 }
420 
421 
422 void
423 PhaseFieldElement :: giveInternalForcesVector_u(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
424 {
425  // computes int_V ( B_u^t * BSgima_u ) * dV
426  this->giveInternalForcesVectorGen(answer, tStep, useUpdatedGpRecord,
427  NULL,
428  &this->computeBmatrixAt,
429  NULL,
430  &this->computeBStress_u,
431  &this->computeVolumeAround
432  );
433 }
434 
435 void
436 PhaseFieldElement :: giveInternalForcesVector_d(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
437 {
438  // computes int_V ( N_d^t * NSgima_d + B_d^t * BSgima_d ) * dV
439  this->giveInternalForcesVectorGen(answer, tStep, useUpdatedGpRecord,
440  &this->computeNmatrixAt,
441  &this->computeBmatrixAt,
442  &this->computeNStress_d,
443  &this->computeBStress_d,
444  &this->computeVolumeAround
445  );
446 }
447 
448 void
449 PhaseFieldElement :: computeBStress_d(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, int useUpdatedGpRecord)
450 {
451  // computes g_c*l
452  // PhaseFieldCrossSection *cs = static_cast...
453  double l = this->giveInternalLength();
454  double g_c = this->giveCriticalEnergy();
455  answer.resize(1);
456  answer.at(1) = g_c*l;
457 }
458 
459 void
461 {
462 
463 }
464 void
466 {
467  StructuralCrossSection *cs = this->giveStructuralCrossSection();
468 
469  this->computeStiffnessMatrixGen(answer, rMode, tStep,
470  NULL,
471  &this->computeBmatrixAt,
472  NULL,
473  &this->Duu_B,
474  this->computeVolumeAround
475  );
476 }
477 
478 void
480 {
481  // This is the regular stiffness matrix times G
482 
483  this->computeStiffnessMatrixGen(answer, rMode, tStep,
484  NULL,
485  &this->computeBmatrixAt,
486  NULL,
487  &this->Duu_B,
488  this->computeVolumeAround
489  );
490 }
491 
492 void
493 PhaseFieldElement :: Duu_B(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
494 {
495  // G*dSig/dEps
496  //StructuralCrossSection *cs = this->giveStructuralCrossSection();
497  //cs->giveCharMaterialStiffnessMatrix(answer, rMode, gp, tStep);
498 
499  //answer.times( this->computeG(gp) );
500 }
501 #endif
502 
503 
504 } // end namespace oofem
CrossSection * giveCrossSection()
Definition: element.C:495
virtual void evalN(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)=0
Evaluates the array of interpolation functions (shape functions) at given point.
void giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
Class and object Domain.
Definition: domain.h:115
void computeVectorOf(ValueModeType u, TimeStep *tStep, FloatArray &answer)
Returns local vector of unknowns.
Definition: element.C:86
virtual int computeNumberOfDofs()
Computes or simply returns total number of element&#39;s local DOFs.
Definition: element.h:424
virtual IntegrationRule * giveIntegrationRule(int i)
Definition: element.h:835
Abstract base class for "structural" finite elements with geometrical nonlinearities.
Total Lagrange.
Definition: fmode.h:44
double computeGPrim(GaussPoint *gp, ValueModeType valueMode, TimeStep *stepN)
void computeStiffnessMatrix_dd(FloatMatrix &, MatResponseMode, TimeStep *)
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
This class implements a structural material status information.
Definition: structuralms.h:65
void clear()
Clears receiver (zero size).
Definition: floatarray.h:206
void giveInternalForcesVector_u(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
EngngModel * giveEngngModel()
Returns engineering model to which receiver is associated.
Definition: domain.C:433
void computeLocationArrayOfDofIDs(const IntArray &dofIdArray, IntArray &answer)
Base class for dof managers.
Definition: dofmanager.h:113
#define S(p)
Definition: mdm.C:481
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
MatResponseMode
Describes the character of characteristic material matrix.
virtual int giveNumberOfDofManagers() const
Definition: element.h:656
virtual void giveCharMaterialStiffnessMatrix(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)=0
Computes the stiffness matrix of receiver in given integration point, respecting its history...
virtual FEInterpolation * giveInterpolation() const
Definition: element.h:629
void computeDamageUnknowns(FloatArray &answer, ValueModeType valueMode, TimeStep *stepN)
double giveTimeIncrement()
Returns solution step associated time increment.
Definition: timestep.h:150
void plusProductUnsym(const FloatMatrix &a, const FloatMatrix &b, double dV)
Adds to the receiver the product .
Definition: floatmatrix.C:780
virtual IRResultType initializeFrom(InputRecord *ir)
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: element.h:518
virtual bool isCharacteristicMtrxSymmetric(MatResponseMode mode)=0
Check for symmetry of stiffness matrix.
double computeG(GaussPoint *gp, ValueModeType valueMode, TimeStep *stepN)
int giveNumberOfDofs() const
Definition: dofmanager.C:279
double dotProduct(const FloatArray &x) const
Computes the dot product (or inner product) of receiver and argument.
Definition: floatarray.C:463
void computeStiffnessMatrix_uu(FloatMatrix &, MatResponseMode, TimeStep *)
void clear()
Clears the array (zero size).
Definition: intarray.h:177
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
void times(double f)
Multiplies receiver by factor f.
Definition: floatmatrix.C:1594
void computeDisplacementUnknowns(FloatArray &answer, ValueModeType valueMode, TimeStep *stepN)
virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)=0
Computes the stress vector of receiver at given integration point, at time step tStep.
Wrapper around element definition to provide FEICellGeometry interface.
Definition: feinterpol.h:95
virtual fMode giveFormulation()
Indicates type of non linear computation (total or updated formulation).
Definition: engngm.h:1069
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
#define N(p, q)
Definition: mdm.C:367
const FloatArray & giveTempStrainVector() const
Returns the const pointer to receiver&#39;s temporary strain vector.
Definition: structuralms.h:115
void plusProductSymmUpper(const FloatMatrix &a, const FloatMatrix &b, double dV)
Adds to the receiver the product .
Definition: floatmatrix.C:698
virtual void giveDofManDofIDMask_u(IntArray &answer)=0
void beTProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix^T and anArray.
Definition: floatarray.C:708
void beNMatrixOf(const FloatArray &n, int nsd)
Assigns the receiver to be a repeated diagonal matrix.
Definition: floatmatrix.C:505
int giveGeometryMode()
Returns the geometry mode describing the formulation used in the internal work 0 - Engineering (small...
void computeStiffnessMatrix_du(FloatMatrix &, MatResponseMode, TimeStep *)
Class representing vector of real numbers.
Definition: floatarray.h:82
double computeDamageAt(GaussPoint *gp, ValueModeType valueMode, TimeStep *stepN)
bool hasDofID(DofIDItem id) const
Checks if receiver contains dof with given ID.
Definition: dofmanager.C:166
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
virtual NLStructuralElement * giveElement()=0
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
Definition: gausspoint.h:205
void computeStiffnessMatrix_ud(FloatMatrix &, MatResponseMode, TimeStep *)
void assemble(const FloatArray &fe, const IntArray &loc)
Assembles the array fe (typically, the load vector of a finite element) into the receiver, using loc as location array.
Definition: floatarray.C:551
void computeBStress_u(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, int useUpdatedGpRecord)
void computeNStress_d(FloatArray &answer, GaussPoint *gp, TimeStep *tStep, int useUpdatedGpRecord)
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
Class representing the general Input Record.
Definition: inputrecord.h:101
void giveInternalForcesVector_d(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
Dof * giveDofWithID(int dofID) const
Returns DOF with given dofID; issues error if not present.
Definition: dofmanager.C:119
void followedBy(const IntArray &b, int allocChunk=0)
Appends array b at the end of receiver.
Definition: intarray.C:145
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual void giveDofManDofIDMask_d(IntArray &answer)=0
void setColumn(const FloatArray &src, int c)
Sets the values of the matrix in specified column.
Definition: floatmatrix.C:648
void times(double s)
Multiplies receiver with scalar.
Definition: floatarray.C:818
void beTranspositionOf(const FloatMatrix &src)
Assigns to the receiver the transposition of parameter.
Definition: floatmatrix.C:323
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int lowerIndx=1, int upperIndx=ALL_STRAINS)=0
Computes the geometrical matrix of receiver in given integration point.
double computeFreeEnergy(GaussPoint *gp, TimeStep *tStep)
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
Domain * giveDomain() const
Definition: femcmpnn.h:100
virtual void computeBd_matrixAt(GaussPoint *, FloatMatrix &, int=1, int=ALL_STRAINS)
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:337
int giveSize() const
Definition: intarray.h:203
const FloatArray & giveTempStressVector() const
Returns the const pointer to receiver&#39;s temporary stress vector.
Definition: structuralms.h:117
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
Abstract base class for all structural cross section models.
the oofem namespace is to define a context or scope in which all oofem names are defined.
DofManager * giveDofManager(int i) const
Definition: element.C:514
void assemble(const FloatMatrix &src, const IntArray &loc)
Assembles the contribution using localization array into receiver.
Definition: floatmatrix.C:251
void clear()
Sets size of receiver to be an empty matrix. It will have zero rows and zero columns size...
Definition: floatmatrix.h:516
PhaseFieldElement(int i, Domain *aDomain)
Abstract class Dof represents Degree Of Freedom in finite element mesh.
Definition: dof.h:93
void symmetrized()
Initializes the lower half of the receiver according to the upper half.
Definition: floatmatrix.C:1576
virtual void computeStiffnessMatrix(FloatMatrix &, MatResponseMode, TimeStep *)
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
void add(const FloatArray &src)
Adds array src to receiver.
Definition: floatarray.C:156
virtual double evaldNdx(FloatMatrix &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)=0
Evaluates the matrix of derivatives of interpolation functions (shape functions) at given point...
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
virtual void computeNd_matrixAt(const FloatArray &lCoords, FloatMatrix &N)
void plusProduct(const FloatMatrix &b, const FloatArray &s, double dV)
Adds the product .
Definition: floatarray.C:226

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:30 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011