OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
lattice2d.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  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #include "../sm/Elements/lattice2d.h"
36 #include "../sm/Materials/latticematstatus.h"
37 #include "../sm/Elements/latticestructuralelement.h"
39 #include "domain.h"
40 #include "node.h"
41 #include "material.h"
42 #include "gausspoint.h"
43 #include "gaussintegrationrule.h"
44 #include "floatmatrix.h"
45 #include "intarray.h"
46 #include "floatarray.h"
47 #include "mathfem.h"
48 #include "datastream.h"
49 #include "contextioerr.h"
50 #include "classfactory.h"
51 
52 #ifdef __OOFEG
53  #include "oofeggraphiccontext.h"
54  #include "../sm/Materials/structuralmaterial.h"
55 #endif
56 
57 namespace oofem {
58 REGISTER_Element(Lattice2d);
59 
61 {
62  numberOfDofMans = 2;
63 
64  length = 0.;
65  pitch = 10.; // a dummy value
67 }
68 
70 { }
71 
72 
73 int
75 {
77  LatticeMaterialStatus *status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
78  return status->giveCrackFlag();
79 }
80 
81 
82 double
84 {
86  LatticeMaterialStatus *status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
87  return status->giveCrackWidth();
88 }
89 
90 double
92 {
93  LatticeMaterialStatus *status;
94 
96  GaussPoint *gp = iRule->getIntegrationPoint(0);
97  status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
98  return status->giveOldCrackWidth();
99 }
100 
101 
102 double
104 {
106  LatticeMaterialStatus *status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
107  return status->giveDissipation();
108 }
109 
110 
111 double
113 {
115  LatticeMaterialStatus *status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
116  return status->giveDeltaDissipation();
117 }
118 
119 void
121 // Returns the strain matrix of the receiver.
122 {
123  double l = this->giveLength();
124  double x1, x2, y1, y2, xp, yp;
125 
126  //Coordinates of the nodes
127  x1 = this->giveNode(1)->giveCoordinate(1);
128  y1 = this->giveNode(1)->giveCoordinate(2);
129  x2 = this->giveNode(2)->giveCoordinate(1);
130  y2 = this->giveNode(2)->giveCoordinate(2);
131 
132  xp = this->gpCoords.at(1);
133  yp = this->gpCoords.at(2);
134 
135  double ecc;
136 
137  double areaHelp = 0.5 * ( x1 * y2 + x2 * yp + xp * y1 - ( xp * y2 + yp * x1 + x2 * y1 ) );
138 
139  ecc = 2 * areaHelp / l;
140 
141  // Assemble Bmatrix (used to compute strains and rotations
142  answer.resize(3, 6);
143  answer.zero();
144  answer.at(1, 1) = -1.;
145  answer.at(1, 2) = 0.;
146  answer.at(1, 3) = ecc;
147  answer.at(1, 4) = 1.;
148  answer.at(1, 5) = 0.;
149  answer.at(1, 6) = -ecc;
150 
151  answer.at(2, 1) = 0.;
152  answer.at(2, 2) = -1.;
153  answer.at(2, 3) = -l / 2.;
154  answer.at(2, 4) = 0.;
155  answer.at(2, 5) = 1.;
156  answer.at(2, 6) = -l / 2.;
157 
158  answer.at(3, 1) = 0.;
159  answer.at(3, 2) = 0.;
160  answer.at(3, 3) = -this->width / sqrt(12.);
161  answer.at(3, 4) = 0.;
162  answer.at(3, 5) = 0.;
163  answer.at(3, 6) = this->width / sqrt(12.);
164 
165  answer.times(1. / l);
166 }
167 
168 void
170 {
171  this->giveStructuralCrossSection()->giveCharMaterialStiffnessMatrix(answer, rMode, gp, tStep);
172  //this->giveStructuralCrossSection()->give2dLatticeStiffMtrx(answer, rMode, gp, tStep);
173 }
174 
175 void
177 {
178  this->giveStructuralCrossSection()->giveRealStresses(answer, gp, strain, tStep);
179  //this->giveStructuralCrossSection()->giveLatticeStress(answer, rMode, gp, tStep);
180 }
181 
182 void
184  TimeStep *tStep)
185 // Computes numerically the stiffness matrix of the receiver.
186 {
187  double dV;
188  FloatMatrix d, bj, dbj;
189  answer.resize(6, 6);
190  answer.zero();
191  this->computeBmatrixAt(integrationRulesArray [ 0 ]->getIntegrationPoint(0), bj);
192  this->computeConstitutiveMatrixAt(d, rMode, integrationRulesArray [ 0 ]->getIntegrationPoint(0), tStep);
193  dV = this->computeVolumeAround( integrationRulesArray [ 0 ]->getIntegrationPoint(0) );
194  dbj.beProductOf(d, bj);
195  answer.plusProductUnsym(bj, dbj, dV);
196 }
197 
198 
200 // Sets up the array of Gauss Points of the receiver.
201 {
202  // the gauss point is used only when methods from crosssection and/or material
203  // classes are requested
204  integrationRulesArray.resize(1);
205  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 3) );
206  integrationRulesArray [ 0 ]->SetUpPointsOnLine(1, _2dLattice);
207 }
208 
209 
210 bool
212 {
213  double sine, cosine;
214  answer.resize(6, 6);
215  answer.zero();
216 
217  sine = sin( this->givePitch() );
218  cosine = cos(pitch);
219  answer.at(1, 1) = cosine;
220  answer.at(1, 2) = sine;
221  answer.at(2, 1) = -sine;
222  answer.at(2, 2) = cosine;
223  answer.at(3, 3) = 1.;
224  answer.at(4, 4) = cosine;
225  answer.at(4, 5) = sine;
226  answer.at(5, 4) = -sine;
227  answer.at(5, 5) = cosine;
228  answer.at(6, 6) = 1.;
229  return 1;
230 }
231 
232 
233 double
235 {
236  double area = this->width * this->thickness;
237  double weight = gp->giveWeight();
238  return weight * 0.5 * this->giveLength() * area;
239 }
240 
241 
242 void
244 {
245  answer = {
246  D_u, D_v, R_w
247  };
248 }
249 
250 
252 // Returns the length of the receiver.
253 {
254  double dx, dy;
255  Node *nodeA, *nodeB;
256 
257  if ( length == 0. ) {
258  nodeA = this->giveNode(1);
259  nodeB = this->giveNode(2);
260  dx = nodeB->giveCoordinate(1) - nodeA->giveCoordinate(1);
261  dy = nodeB->giveCoordinate(2) - nodeA->giveCoordinate(2);
262  length = sqrt(dx * dx + dy * dy);
263  }
264 
265  return length;
266 }
267 
268 
270 // Returns the pitch of the receiver.
271 {
272  double xA, xB, yA, yB;
273  Node *nodeA, *nodeB;
274 
275  if ( pitch == 10. ) { // 10. : dummy initialization value
276  nodeA = this->giveNode(1);
277  nodeB = this->giveNode(2);
278  xA = nodeA->giveCoordinate(1);
279  xB = nodeB->giveCoordinate(1);
280  yA = nodeA->giveCoordinate(2);
281  yB = nodeB->giveCoordinate(2);
282  pitch = atan2(yB - yA, xB - xA);
283  }
284 
285  return pitch;
286 }
287 
288 double
290 {
292  LatticeMaterialStatus *status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
293  double normalStress = 0;
294  normalStress = status->giveNormalStress();
295 
296  return normalStress;
297 }
298 
299 
300 double
302 {
303  LatticeMaterialStatus *status;
304 
306  GaussPoint *gp = iRule->getIntegrationPoint(0);
307  status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
308  double normalStress = 0;
309  normalStress = status->giveOldNormalStress();
310 
311  return normalStress;
312 }
313 
314 int
316 {
317  LatticeMaterialStatus *status;
318 
320  GaussPoint *gp = iRule->getIntegrationPoint(0);
321  status = static_cast< LatticeMaterialStatus * >( gp->giveMaterialStatus() );
322  int updateFlag = 0;
323  updateFlag = status->hasBeenUpdated();
324 
325  return updateFlag;
326 }
327 
328 
329 int
331 //
332 // returns a unit vectors of local coordinate system at element
333 // stored rowwise (mainly used by some materials with ortho and anisotrophy)
334 //
335 {
336  double sine, cosine;
337 
338  answer.resize(3, 3);
339  answer.zero();
340 
341  sine = sin( this->givePitch() );
342  cosine = cos(pitch);
343 
344  answer.at(1, 1) = cosine;
345  answer.at(1, 2) = sine;
346  answer.at(2, 1) = -sine;
347  answer.at(2, 2) = cosine;
348  answer.at(3, 3) = 1.0;
349 
350  return 1;
351 }
352 
355 {
356  IRResultType result; // Required by IR_GIVE_FIELD macro
357 
359 
361 
363 
364  couplingFlag = 0;
366 
369  if ( couplingFlag == 1 ) {
371  }
372 
373  // first call parent
375 }
376 
377 
378 int
380 {
381  answer.resize(3);
382  answer.at(1) = this->gpCoords.at(1);
383  answer.at(2) = this->gpCoords.at(2);
384 
385  return 1;
386 }
387 
388 
389 void
391 {
392  answer.resize(3);
393  answer.at(1) = this->gpCoords.at(1);
394  answer.at(2) = this->gpCoords.at(2);
395 }
396 
397 
399 {
400  contextIOResultType iores;
401 
402  if ( ( iores = LatticeStructuralElement :: saveContext(stream, mode, obj) ) != CIO_OK ) {
403  THROW_CIOERR(iores);
404  }
405 
406  if ( ( mode & CM_Definition ) ) {
407 
408  if ( !stream.write(width) ) {
410  }
411 
412  if ( !stream.write(thickness) ) {
414  }
415 
416  if ( !stream.write(couplingFlag) ) {
418  }
419 
420  if ( ( iores = couplingNumbers.storeYourself(stream) ) != CIO_OK ) {
421  THROW_CIOERR(iores);
422  }
423 
424  if ( ( iores = gpCoords.storeYourself(stream) ) != CIO_OK ) {
425  THROW_CIOERR(iores);
426  }
427  }
428 
429  return CIO_OK;
430 }
431 
432 
433 
435 {
436  contextIOResultType iores;
437 
438  if ( ( iores = LatticeStructuralElement :: restoreContext(stream, mode, obj) ) != CIO_OK ) {
439  THROW_CIOERR(iores);
440  }
441 
442  if ( mode & CM_Definition ) {
443 
444  if ( !stream.read(width) ) {
446  }
447 
448  if ( !stream.read(thickness) ) {
450  }
451 
452  if ( !stream.read(couplingFlag) ) {
454  }
455 
456  if ( ( iores = couplingNumbers.restoreYourself(stream) ) != CIO_OK ) {
457  THROW_CIOERR(iores);
458  }
459 
460  if ( ( iores = gpCoords.restoreYourself(stream) ) != CIO_OK ) {
461  THROW_CIOERR(iores);
462  }
463 
464  }
465 
466  return CIO_OK;
467 }
468 
469 #ifdef __OOFEG
470 
471 void
473 {
475 
476  if ( mode == OGC_rawGeometry ) {
477  this->drawRawGeometry(gc, tStep);
478  this->drawRawCrossSections(gc, tStep);
479  } else if ( mode == OGC_deformedGeometry ) {
480  this->drawDeformedGeometry(gc, tStep, DisplacementVector);
481  } else if ( mode == OGC_eigenVectorGeometry ) {
482  this->drawDeformedGeometry(gc, tStep, EigenVector);
483  } else if ( mode == OGC_scalarPlot ) {
484  this->drawScalar(gc, tStep);
485  } else if ( mode == OGC_elemSpecial ) {
486  this->drawSpecial(gc, tStep);
487  } else {
488  OOFEM_ERROR("unsupported mode");
489  }
490 }
491 
492 void
494 {
495  GraphicObj *go;
496 
497  WCRec p [ 2 ]; /* points */
498  if ( !gc.testElementGraphicActivity(this) ) {
499  return;
500  }
501 
502  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
503  EASValsSetColor( gc.getElementColor() );
504  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
505 
506  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
507  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveCoordinate(2);
508  p [ 0 ].z = 0.;
509  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
510  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveCoordinate(2);
511  p [ 1 ].z = 0.;
512 
513  go = CreateLine3D(p);
514  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
515  EGAttachObject(go, ( EObjectP ) this);
516  EMAddGraphicsToModel(ESIModel(), go);
517 }
518 
519 void
521 {
522  GraphicObj *go;
523 
524  if ( !gc.testElementGraphicActivity(this) ) {
525  return;
526  }
527 
528  double defScale = gc.getDefScale();
529 
530  WCRec p [ 2 ]; /* poin */
531  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
532  EASValsSetColor( gc.getDeformedElementColor() );
533  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
534 
535  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
536  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale);
537  p [ 0 ].z = 0.;
538 
539  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
540  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale);
541  p [ 1 ].z = 0.;
542 
543  go = CreateLine3D(p);
544  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
545  EGAttachObject(go, ( EObjectP ) this);
546  EMAddGraphicsToModel(ESIModel(), go);
547 }
548 
549 void
551 {
552  WCRec p [ 2 ];
553  GraphicObj *tr;
554  GaussPoint *gp;
555  FloatArray crackStatuses, cf;
556 
557  if ( !gc.testElementGraphicActivity(this) ) {
558  return;
559  }
560 
561  if ( gc.giveIntVarType() == IST_CrackState ) {
562  gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
563  this->giveIPValue(crackStatuses, gp, IST_CrackStatuses, tStep);
564  if ( crackStatuses(0) == 1. || crackStatuses(0) == 2. || crackStatuses(0) == 3 || crackStatuses(0) == 4 ) {
565  FloatArray coords;
566  this->giveCrossSectionCoordinates(coords);
567 
568  p [ 0 ].x = ( FPNum ) coords.at(1);
569  p [ 0 ].y = ( FPNum ) coords.at(2);
570  p [ 0 ].z = ( FPNum ) coords.at(3);
571  p [ 1 ].x = ( FPNum ) coords.at(4);
572  p [ 1 ].y = ( FPNum ) coords.at(5);
573  p [ 1 ].z = ( FPNum ) coords.at(6);
574 
575 
576  EASValsSetLayer(OOFEG_CRACK_PATTERN_LAYER);
577  EASValsSetLineWidth(OOFEG_CRACK_PATTERN_WIDTH);
578  if ( ( crackStatuses(0) == 1. ) ) {
579  EASValsSetColor( gc.getActiveCrackColor() );
580  } else if ( crackStatuses(0) == 2. ) {
581  EASValsSetColor( gc.getCrackPatternColor() );
582  } else if ( crackStatuses(0) == 3. ) {
583  EASValsSetColor( gc.getActiveCrackColor() );
584  } else if ( crackStatuses(0) == 4. ) {
585  EASValsSetColor( gc.getActiveCrackColor() );
586  }
587 
588 
589  tr = CreateLine3D(p);
590  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, tr);
591  EGAttachObject(tr, ( EObjectP ) this);
592  EMAddGraphicsToModel(ESIModel(), tr);
593  }
594  }
595 }
596 
597 void
599 {
600  GraphicObj *go;
601 
602  // if (!go) { // create new one
603  WCRec p [ 2 ]; /* poin */
604  if ( !gc.testElementGraphicActivity(this) ) {
605  return;
606  }
607 
608  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
609  EASValsSetColor( gc.getCrossSectionColor() );
610  EASValsSetLayer(OOFEG_RAW_CROSSSECTION_LAYER);
611 
612  FloatArray coords;
613  this->giveCrossSectionCoordinates(coords);
614 
615  p [ 0 ].x = ( FPNum ) coords.at(1);
616  p [ 0 ].y = ( FPNum ) coords.at(2);
617  p [ 0 ].z = ( FPNum ) coords.at(3);
618  p [ 1 ].x = ( FPNum ) coords.at(4);
619  p [ 1 ].y = ( FPNum ) coords.at(5);
620  p [ 1 ].z = ( FPNum ) coords.at(6);
621 
622  go = CreateLine3D(p);
623  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
624  EGAttachObject(go, ( EObjectP ) this);
625  EMAddGraphicsToModel(ESIModel(), go);
626 }
627 
628 void
630 {
631  double x1, y1, x2, y2;
632  x1 = this->giveNode(1)->giveCoordinate(1);
633  y1 = this->giveNode(1)->giveCoordinate(2);
634  x2 = this->giveNode(2)->giveCoordinate(1);
635  y2 = this->giveNode(2)->giveCoordinate(2);
636 
637  //Compute normal and shear direction
638  FloatArray normalDirection;
639  FloatArray shearDirection;
640  normalDirection.resize(2);
641  normalDirection.zero();
642  shearDirection.resize(2);
643  shearDirection.zero();
644  normalDirection.at(1) = x2 - x1;
645  normalDirection.at(2) = y2 - y1;
646  normalDirection.normalize();
647  if ( normalDirection.at(2) == 0. ) {
648  shearDirection.at(1) = 0.;
649  shearDirection.at(2) = 1.;
650  } else {
651  shearDirection.at(1) = 1.0;
652  shearDirection.at(2) =
653  -normalDirection.at(1) / normalDirection.at(2);
654  }
655 
656  shearDirection.normalize();
657 
658  coords.resize(6);
659  coords.at(1) = this->gpCoords.at(1) - shearDirection.at(1) * this->width / 2.;
660  coords.at(2) = this->gpCoords.at(2) - shearDirection.at(2) * this->width / 2.;
661  coords.at(3) = 0.;
662  coords.at(4) = this->gpCoords.at(1) + shearDirection.at(1) * this->width / 2.;
663  coords.at(5) = this->gpCoords.at(2) + shearDirection.at(2) * this->width / 2.;
664  coords.at(6) = 0.;
665 }
666 
667 #endif
668 } // end namespace oofem
virtual ~Lattice2d()
Definition: lattice2d.C:69
int testElementGraphicActivity(Element *)
Test if particular element passed fulfills various filtering criteria for its graphics output...
contextIOResultType storeYourself(DataStream &stream) const
Stores array to output stream.
Definition: intarray.C:289
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
virtual int giveLocalCoordinateSystem(FloatMatrix &answer)
Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy...
Definition: lattice2d.C:330
Class and object Domain.
Definition: domain.h:115
virtual IntegrationRule * giveDefaultIntegrationRulePtr()
Access method for default integration rule.
Definition: element.h:822
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Restores the receiver state previously written in stream.
Definition: lattice2d.C:434
virtual int giveCrackFlag()
Returns the crack flag.
virtual double giveNormalStress()
Gives the last equilibrated normal stress.
The purpose of DataStream abstract class is to allow to store/restore context to different streams...
Definition: datastream.h:54
contextIOResultType storeYourself(DataStream &stream) const
Definition: floatarray.C:872
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Stores receiver state to output stream.
Definition: lattice2d.C:398
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: lattice2d.C:493
void zero()
Sets all component to zero.
Definition: intarray.C:52
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual void giveCrossSectionCoordinates(FloatArray &coords)
Definition: lattice2d.C:629
#define OOFEG_RAW_CROSSSECTION_LAYER
virtual int giveCrackFlag()
Returns the crack flag.
Definition: lattice2d.C:74
#define OOFEG_RAW_GEOMETRY_LAYER
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: lattice2d.C:234
#define _IFT_Lattice2d_couplingflag
Definition: lattice2d.h:46
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: lattice2d.C:199
virtual int hasBeenUpdated()
Returns flag indicating if status has been updated.
Definition: lattice2d.C:315
virtual void giveGpCoordinates(FloatArray &coords)
Gives the GP coordinates.
Definition: lattice2d.C:390
virtual bool computeGtoLRotationMatrix(FloatMatrix &)
Returns transformation matrix from global c.s.
Definition: lattice2d.C:211
General IO error.
This class implements a base lattice material status.
virtual double giveDissipation()
Returns the energy dissipation computed at the GaussPoint of the element.
Definition: lattice2d.C:103
virtual double giveCoordinate(int i)
Definition: node.C:82
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
virtual int read(int *data, int count)=0
Reads count integer values into array pointed by data.
MatResponseMode
Describes the character of characteristic material matrix.
#define THROW_CIOERR(e)
Definition: contextioerr.h:61
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...
#define OOFEG_DEFORMED_GEOMETRY_LAYER
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
Definition: element.h:1007
Abstract base class representing integration rule.
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Restores the receiver state previously written in stream.
Definition: element.C:970
virtual int hasBeenUpdated()
Gives the last equilibrated normal stress.
virtual int write(const int *data, int count)=0
Writes count integer values from array pointed by data.
void plusProductUnsym(const FloatMatrix &a, const FloatMatrix &b, double dV)
Adds to the receiver the product .
Definition: floatmatrix.C:780
InternalStateType giveIntVarType()
#define OOFEG_CRACK_PATTERN_LAYER
FloatArray gpCoords
Definition: lattice2d.h:60
virtual double giveDeltaDissipation()
Returns the increment of dissipation computed at the GaussPoint of the element.
Definition: lattice2d.C:112
virtual double giveDeltaDissipation()
Returns the increment of dissipation computed at the GaussPoint of the element.
#define OOFEM_ERROR(...)
Definition: error.h:61
REGISTER_Element(LSpace)
virtual void computeBmatrixAt(GaussPoint *, FloatMatrix &, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: lattice2d.C:120
#define OOFEG_RAW_GEOMETRY_WIDTH
virtual double giveWeight()
Returns integration weight of receiver.
Definition: gausspoint.h:181
UnknownType
Type representing particular unknown (its physical meaning).
Definition: unknowntype.h:55
StructuralCrossSection * giveStructuralCrossSection()
Helper function which returns the structural cross-section for the element.
void times(double f)
Multiplies receiver by factor f.
Definition: floatmatrix.C:1594
#define OOFEG_CRACK_PATTERN_WIDTH
virtual void drawSpecial(oofegGraphicContext &gc, TimeStep *tStep)
Definition: lattice2d.C:550
contextIOResultType restoreYourself(DataStream &stream)
Definition: floatarray.C:895
#define _IFT_Lattice2d_gpcoords
Definition: lattice2d.h:45
virtual double giveUpdatedCoordinate(int ic, TimeStep *tStep, double scale=1.)
Returns updated ic-th coordinate of receiver.
Definition: node.C:245
virtual double giveOldCrackWidth()
Definition: lattice2d.C:91
double givePitch()
Definition: lattice2d.C:269
double at(int i, int j) const
Coefficient access function.
Definition: floatmatrix.h:176
void resize(int n)
Checks size of receiver towards requested bounds.
Definition: intarray.C:124
void giveRealStresses(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
Computes the real stress vector for given strain and integration point.
double thickness
Definition: lattice2d.h:59
virtual double giveCrackWidth()
Definition: lattice2d.C:83
contextIOResultType restoreYourself(DataStream &stream)
Restores array from image on stream.
Definition: intarray.C:305
#define _IFT_Lattice2d_width
Definition: lattice2d.h:44
Class representing vector of real numbers.
Definition: floatarray.h:82
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 double giveLength()
Returns the element length.
Definition: lattice2d.C:251
GaussPoint * getIntegrationPoint(int n)
Access particular integration point of receiver.
virtual double giveDissipation()
Returns the energy dissipation computed at the GaussPoint of the element.
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
Definition: gausspoint.h:205
#define _IFT_Lattice2d_thick
Definition: lattice2d.h:43
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes numerically stiffness matrix of receiver.
Definition: lattice2d.C:183
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
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: lattice2d.C:243
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Stores receiver state to output stream.
Definition: element.C:885
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: lattice2d.C:169
virtual void drawRawCrossSections(oofegGraphicContext &gc, TimeStep *tStep)
Definition: lattice2d.C:598
long ContextMode
Context mode (mask), defining the type of information written/read to/from context.
Definition: contextmode.h:43
std::vector< std::unique_ptr< IntegrationRule > > integrationRulesArray
List of integration rules of receiver (each integration rule contains associated integration points a...
Definition: element.h:170
#define CM_Definition
Definition: contextmode.h:47
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: lattice2d.C:354
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
#define _IFT_Lattice2d_couplingnumber
Definition: lattice2d.h:47
virtual double giveOldCrackWidth()
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
Definition: lattice2d.C:520
virtual double giveNormalStress()
Returns the normal stress.
Definition: lattice2d.C:289
virtual double giveOldNormalStress()
Returns the old normal stress.
Definition: lattice2d.C:301
Lattice2d(int n, Domain *d)
Definition: lattice2d.C:60
This class implements the base of a special lattice element following the concepts orginally develope...
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:337
virtual void drawYourself(oofegGraphicContext &gc, TimeStep *tStep)
Definition: lattice2d.C:472
the oofem namespace is to define a context or scope in which all oofem names are defined.
Class implementing node in finite element mesh.
Definition: node.h:87
IntArray couplingNumbers
Definition: lattice2d.h:62
double normalize()
Normalizes receiver.
Definition: floatarray.C:828
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords)
This function is different from the standard computeGlobalCorrdinates function as it returns the glob...
Definition: lattice2d.C:379
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
int numberOfDofMans
Number of dofmanagers.
Definition: element.h:149
virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)
Computes the stress vector of receiver at given integration point, at time step tStep.
Definition: lattice2d.C:176
Class representing Gaussian-quadrature integration rule.
virtual double giveOldNormalStress()
Gives the last equilibrated normal stress.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
OGC_PlotModeType giveIntVarPlotMode()

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:29 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011