OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
truss1d.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/Bars/truss1d.h"
36 #include "../sm/CrossSections/structuralcrosssection.h"
37 #include "fei1dlin.h"
38 #include "node.h"
39 #include "material.h"
40 #include "gausspoint.h"
41 #include "gaussintegrationrule.h"
42 #include "floatmatrix.h"
43 #include "floatarray.h"
44 #include "intarray.h"
45 #include "mathfem.h"
46 #include "classfactory.h"
47 
48 #ifdef __OOFEG
49  #include "oofeggraphiccontext.h"
50 #endif
51 
52 namespace oofem {
53 REGISTER_Element(Truss1d);
54 
55 FEI1dLin Truss1d :: interp(1); // Initiates the static interpolator
56 
57 
58 Truss1d :: Truss1d(int n, Domain *aDomain) :
59  StructuralElement(n, aDomain),
64 {
65  numberOfDofMans = 2;
66 }
67 
68 
70 // Sets up the array of Gauss Points of the receiver.
71 {
72  if ( integrationRulesArray.size() == 0 ) {
73  integrationRulesArray.resize( 1 );
74  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 2) );
76  }
77 }
78 
80 Truss1d :: giveInterpolation() const { return & interp; }
81 
82 void
84 // Returns the lumped mass matrix of the receiver. This expression is
85 // valid in both local and global axes.
86 {
87  GaussPoint *gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
88  double density = this->giveStructuralCrossSection()->give('d', gp);
89  double halfMass = density * this->giveCrossSection()->give(CS_Area, gp) * this->computeLength() * 0.5;
90  answer.resize(2, 2);
91  answer.zero();
92  answer.at(1, 1) = halfMass;
93  answer.at(2, 2) = halfMass;
94 }
95 
96 
97 void
98 Truss1d :: computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int li, int ui)
99 //
100 // Returns linear part of geometrical equations of the receiver at gp.
101 // Returns the linear part of the B matrix
102 //
103 {
104  FloatMatrix dN;
106  answer.beTranspositionOf(dN);
107 }
108 
109 
110 void
112 //
113 // Returns the [1x2] displacement gradient matrix {BH} of the receiver,
114 // evaluated at gp.
115 // @todo not checked if correct
116 {
117  this->computeBmatrixAt(gp, answer);
118 }
119 
120 
121 
122 void
124 // Returns the displacement interpolation matrix {N} of the receiver,
125 // evaluated at gp.
126 {
127  FloatArray n;
128  this->interp.evalN( n, iLocCoord, FEIElementGeometryWrapper(this) );
129  // Reshape
130  answer.resize(1, 2);
131  answer.at(1, 1) = n.at(1);
132  answer.at(1, 2) = n.at(2);
133 }
134 
135 
136 double
138 // Returns the length of the receiver. This method is valid only if 1
139 // Gauss point is used.
140 {
141  double detJ = fabs( this->interp.giveTransformationJacobian( gp->giveNaturalCoordinates(), FEIElementGeometryWrapper(this) ) );
142  return detJ *gp->giveWeight() * this->giveCrossSection()->give(CS_Area, gp);
143 }
144 
145 
146 void
148 {
149  this->giveStructuralCrossSection()->giveRealStress_1d(answer, gp, strain, tStep);
150 }
151 
152 void
154 {
155  this->giveStructuralCrossSection()->giveStiffnessMatrix_1d(answer, rMode, gp, tStep);
156 }
157 
158 void
159 Truss1d :: giveDofManDofIDMask(int inode, IntArray &answer) const
160 {
161  answer = {D_u};
162 }
163 
164 
165 void
167  IntArray &localNodeIdArray, IntArray &globalNodeIdArray,
169  int &localNodeId, int &localElemId, int &localBcId,
170  IntArray &controlNode, IntArray &controlDof,
172 {
173  int inode, nodes = 2;
174  FloatArray *corner [ 2 ], midNode, cor [ 2 ];
175  double x = 0.0;
176 
179  for ( inode = 0; inode < nodes; inode++ ) {
180  corner [ inode ] = this->giveNode(inode + 1)->giveCoordinates();
181  if ( corner [ inode ]->giveSize() != 3 ) {
182  cor [ inode ].resize(3);
183  cor [ inode ].at(1) = corner [ inode ]->at(1);
184  cor [ inode ].at(2) = 0.0;
185  cor [ inode ].at(3) = 0.0;
186 
187  corner [ inode ] = & ( cor [ inode ] );
188  }
189 
190  x += corner [ inode ]->at(1);
191  }
192 
193  midNode.resize(3);
194 
195  midNode.at(1) = x / nodes;
196  midNode.at(2) = 0.0;
197  midNode.at(3) = 0.0;
198  }
199 
200  this->setupRefinedElementProblem1D(this, refinedElement, level, nodeId, localNodeIdArray, globalNodeIdArray,
201  sMode, tStep, nodes, corner, midNode,
202  localNodeId, localElemId, localBcId,
203  controlNode, controlDof, aMode, "Truss1d");
204 }
205 
206 
208 {
210 }
211 
212 
213 #ifdef __OOFEG
215 {
216  GraphicObj *go;
217  // if (!go) { // create new one
218  WCRec p [ 2 ]; /* point */
219  if ( !gc.testElementGraphicActivity(this) ) {
220  return;
221  }
222 
223  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
224  EASValsSetColor( gc.getElementColor() );
225  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
226  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
227  p [ 0 ].y = 0.;
228  p [ 0 ].z = 0.0;
229  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
230  p [ 1 ].y = 0.;
231  p [ 1 ].z = 0.0;
232  go = CreateLine3D(p);
233  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
234  EGAttachObject(go, ( EObjectP ) this);
235  EMAddGraphicsToModel(ESIModel(), go);
236 }
237 
238 
240 {
241  GraphicObj *go;
242  double defScale = gc.getDefScale();
243  // if (!go) { // create new one
244  WCRec p [ 2 ]; /* point */
245  if ( !gc.testElementGraphicActivity(this) ) {
246  return;
247  }
248 
249  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
250  EASValsSetColor( gc.getDeformedElementColor() );
251  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
252  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
253  p [ 0 ].y = 0.;
254  p [ 0 ].z = 0.;
255 
256  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
257  p [ 1 ].y = 0.;
258  p [ 1 ].z = 0.;
259  go = CreateLine3D(p);
260  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
261  EMAddGraphicsToModel(ESIModel(), go);
262 }
263 
264 
266 {
267  int i, indx, result = 0;
268  WCRec p [ 2 ];
269  GraphicObj *tr;
270  FloatArray v1, v2;
271  double s [ 2 ], defScale;
272 
273  if ( !gc.testElementGraphicActivity(this) ) {
274  return;
275  }
276 
277  if ( gc.giveIntVarMode() == ISM_recovered ) {
278  result += this->giveInternalStateAtNode(v1, gc.giveIntVarType(), gc.giveIntVarMode(), 1, tStep);
279  result += this->giveInternalStateAtNode(v2, gc.giveIntVarType(), gc.giveIntVarMode(), 2, tStep);
280  } else if ( gc.giveIntVarMode() == ISM_local ) {
281  GaussPoint *gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
282  result += giveIPValue(v1, gp, gc.giveIntVarType(), tStep);
283  v2 = v1;
284  result *= 2;
285  }
286 
287  if ( result != 2 ) {
288  return;
289  }
290 
291  indx = gc.giveIntVarIndx();
292 
293  s [ 0 ] = v1.at(indx);
294  s [ 1 ] = v2.at(indx);
295 
296  EASValsSetLayer(OOFEG_VARPLOT_PATTERN_LAYER);
297 
298  if ( ( gc.getScalarAlgo() == SA_ISO_SURF ) || ( gc.getScalarAlgo() == SA_ISO_LINE ) ) {
299  for ( i = 0; i < 2; i++ ) {
300  if ( gc.getInternalVarsDefGeoFlag() ) {
301  // use deformed geometry
302  defScale = gc.getDefScale();
303  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
304  p [ i ].y = 0.;
305  p [ i ].z = 0.;
306  } else {
307  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
308  p [ i ].y = 0.;
309  p [ i ].z = 0.;
310  }
311  }
312 
313  //EASValsSetColor(gc.getYieldPlotColor(ratio));
314  tr = CreateLine3D(p);
315  EGWithMaskChangeAttributes(LAYER_MASK, tr);
316  EMAddGraphicsToModel(ESIModel(), tr);
317  } else if ( ( gc.getScalarAlgo() == SA_ZPROFILE ) || ( gc.getScalarAlgo() == SA_COLORZPROFILE ) ) {
318  double landScale = gc.getLandScale();
319 
320  for ( i = 0; i < 2; i++ ) {
321  if ( gc.getInternalVarsDefGeoFlag() ) {
322  // use deformed geometry
323  defScale = gc.getDefScale();
324  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
325  p [ i ].y = 0.0;
326  p [ i ].z = s [ i ] * landScale;
327  } else {
328  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
329  p [ i ].y = 0.0;
330  p [ i ].z = s [ i ] * landScale;
331  }
332  }
333 
334  if ( gc.getScalarAlgo() == SA_ZPROFILE ) {
335  /*
336  * EASValsSetColor(gc.getDeformedElementColor());
337  * EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
338  * tr = CreateLine3D(p);
339  * EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, tr);
340  */
341  WCRec pp [ 4 ];
342  pp [ 0 ].x = p [ 0 ].x;
343  pp [ 0 ].y = 0.0;
344  pp [ 0 ].z = 0.0;
345  pp [ 1 ].x = p [ 0 ].x;
346  pp [ 1 ].y = 0.0;
347  pp [ 1 ].z = p [ 0 ].z;
348  pp [ 2 ].x = p [ 1 ].x;
349  pp [ 2 ].y = 0.0;
350  pp [ 2 ].z = p [ 1 ].z;
351  pp [ 3 ].x = p [ 1 ].x;
352  pp [ 3 ].y = 0.0;
353  pp [ 3 ].z = 0.0;
354  tr = CreateQuad3D(pp);
355  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
356  EASValsSetColor( gc.getDeformedElementColor() );
357  //EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
358  EASValsSetFillStyle(FILL_HOLLOW);
359  EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | LAYER_MASK, tr);
360  EMAddGraphicsToModel(ESIModel(), tr);
361  } else {
362  //tr = CreateTriangleWD3D(p, s[0], s[1], s[2]);
363  EASValsSetColor( gc.getDeformedElementColor() );
364  tr = CreateLine3D(p);
365  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, tr);
366  EMAddGraphicsToModel(ESIModel(), tr);
367  }
368  }
369 }
370 
371 #endif
372 
373 
374 Interface *
376 {
377  if ( interface == ZZNodalRecoveryModelInterfaceType ) {
378  return static_cast< ZZNodalRecoveryModelInterface * >(this);
379  } else if ( interface == NodalAveragingRecoveryModelInterfaceType ) {
380  return static_cast< NodalAveragingRecoveryModelInterface * >(this);
381  } else if ( interface == SpatialLocalizerInterfaceType ) {
382  return static_cast< SpatialLocalizerInterface * >(this);
383  } else if ( interface == ZZErrorEstimatorInterfaceType ) {
384  return static_cast< ZZErrorEstimatorInterface * >(this);
385  } else if ( interface == HuertaErrorEstimatorInterfaceType ) {
386  return static_cast< HuertaErrorEstimatorInterface * >(this);
387  }
388 
389  return NULL;
390 }
391 
392 
393 void
395  InternalStateType type, TimeStep *tStep)
396 {
397  GaussPoint *gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
398  this->giveIPValue(answer, gp, type, tStep);
399 }
400 
401 } // end namespace oofem
CrossSection * giveCrossSection()
Definition: element.C:495
int testElementGraphicActivity(Element *)
Test if particular element passed fulfills various filtering criteria for its graphics output...
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
The element interface required by NodalAvergagingRecoveryModel.
virtual void NodalAveragingRecoveryMI_computeNodalValue(FloatArray &answer, int node, InternalStateType type, TimeStep *tStep)
Computes the element value in given node.
Definition: truss1d.C:394
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Class and object Domain.
Definition: domain.h:115
virtual double computeLength()
Computes the length (zero for all but 1D geometries)
Definition: element.C:1129
ScalarAlgorithmType getScalarAlgo()
virtual FEInterpolation * giveInterpolation() const
Definition: truss1d.C:80
virtual void HuertaErrorEstimatorI_computeNmatrixAt(GaussPoint *gp, FloatMatrix &answer)
Definition: truss1d.C:207
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: truss1d.C:159
The element interface required by ZZNodalRecoveryModel.
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
Definition: truss1d.C:265
const FloatArray & giveSubPatchCoordinates()
Returns local sub-patch coordinates of the receiver.
Definition: gausspoint.h:144
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: truss1d.C:153
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual int giveInternalStateAtNode(FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *tStep)
Returns internal state variable (like stress,strain) at node of element in Reduced form...
#define OOFEG_RAW_GEOMETRY_LAYER
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual void evalN(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the array of interpolation functions (shape functions) at given point.
Definition: fei1dlin.C:48
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: truss1d.C:214
virtual Interface * giveInterface(InterfaceType it)
Interface requesting service.
Definition: truss1d.C:375
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: truss1d.C:137
virtual double giveCoordinate(int i)
Definition: node.C:82
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: truss1d.C:147
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
#define OOFEG_DEFORMED_GEOMETRY_LAYER
virtual void HuertaErrorEstimatorI_setupRefinedElementProblem(RefinedElement *refinedElement, int level, int nodeId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray, HuertaErrorEstimatorInterface::SetupMode sMode, TimeStep *tStep, int &localNodeId, int &localElemId, int &localBcId, IntArray &controlNode, IntArray &controlDof, HuertaErrorEstimator::AnalysisMode aMode)
Definition: truss1d.C:166
virtual void computeBHmatrixAt(GaussPoint *gp, FloatMatrix &answer)
Definition: truss1d.C:111
virtual double giveTransformationJacobian(const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation.
Definition: fei1dlin.C:94
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
virtual double evaldNdx(FloatMatrix &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the matrix of derivatives of interpolation functions (shape functions) at given point...
Definition: fei1dlin.C:58
Truss1d(int n, Domain *d)
Definition: truss1d.C:58
InternalStateType giveIntVarType()
Abstract base class for all "structural" finite elements.
The element interface corresponding to ZZErrorEstimator.
The element interface corresponding to HuertaErrorEstimator.
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: truss1d.C:98
REGISTER_Element(LSpace)
#define OOFEG_RAW_GEOMETRY_WIDTH
static FEI1dLin interp
Definition: truss1d.h:62
virtual double giveWeight()
Returns integration weight of receiver.
Definition: gausspoint.h:181
void setupRefinedElementProblem1D(Element *element, RefinedElement *refinedElement, int level, int nodeId, IntArray &localNodeIdArray, IntArray &globalNodeIdArray, HuertaErrorEstimatorInterface::SetupMode mode, TimeStep *tStep, int nodes, FloatArray **corner, FloatArray &midNode, int &localNodeId, int &localElemId, int &localBcId, IntArray &controlNode, IntArray &controlDof, HuertaErrorEstimator::AnalysisMode aMode, const char *edgetype)
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.
Wrapper around element definition to provide FEICellGeometry interface.
Definition: feinterpol.h:95
virtual double giveUpdatedCoordinate(int ic, TimeStep *tStep, double scale=1.)
Returns updated ic-th coordinate of receiver.
Definition: node.C:245
double at(int i, int j) const
Coefficient access function.
Definition: floatmatrix.h:176
virtual int setupIntegrationPoints(IntegrationRule &irule, int npoints, Element *element)
Sets up integration rule for the given element.
Definition: crosssection.C:54
InternalStateMode giveIntVarMode()
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
Definition: truss1d.C:83
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
virtual double give(CrossSectionProperty a, GaussPoint *gp)
Returns the value of cross section property at given point.
Definition: crosssection.C:151
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
Class Interface.
Definition: interface.h:82
virtual void giveRealStress_1d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)=0
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
void beTranspositionOf(const FloatMatrix &src)
Assigns to the receiver the transposition of parameter.
Definition: floatmatrix.C:323
The spatial localizer element interface associated to spatial localizer.
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
virtual FloatArray * giveCoordinates()
Definition: node.h:114
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
virtual void giveStiffnessMatrix_1d(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)=0
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
Definition: truss1d.C:239
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: truss1d.C:123
the oofem namespace is to define a context or scope in which all oofem names are defined.
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: truss1d.C:69
#define OOFEG_VARPLOT_PATTERN_LAYER
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
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
Class representing Gaussian-quadrature integration rule.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631

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