OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
pfemelement.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 "pfemelement.h"
36 #include "domain.h"
37 #include "timestep.h"
38 #include "node.h"
39 #include "dof.h"
40 #include "load.h"
41 #include "boundaryload.h"
42 #include "gausspoint.h"
43 #include "gaussintegrationrule.h"
44 #include "intarray.h"
45 #include "floatarray.h"
46 #include "floatmatrix.h"
47 #include "verbose.h"
48 #include "material.h"
49 #include "elementside.h"
50 #include "mathfem.h"
51 
52 #ifdef __OOFEG
53  #include "oofeggraphiccontext.h"
54  #include "connectivitytable.h"
55 #endif
56 
57 namespace oofem {
59  FMElement(n, aDomain)
60  // Constructor. Creates an element with number n, belonging to aDomain.
61 { }
62 
63 
65 // Destructor.
66 { }
67 
70 {
71  //const char *__proc = "initializeFrom"; // Required by IR_GIVE_FIELD macro
72  //IRResultType result; // Required by IR_GIVE_FIELD macro
73 
75 
76  this->computeGaussPoints();
77  return IRRT_OK;
78 }
79 
80 
81 void
83  CharType mtrx, TimeStep *tStep)
84 //
85 // returns characteristics matrix of receiver according to mtrx
86 //
87 {
88  FloatMatrix reducedAnswer;
89  int size = this->computeNumberOfDofs();
90  answer.resize(size, size);
91  answer.zero();
92 
93  if ( mtrx == LumpedMassMatrix ) {
94  this->computeDiagonalMassMtrx(reducedAnswer, tStep);
95  answer.assemble(reducedAnswer, this->giveVelocityDofMask(), this->giveVelocityDofMask());
96  } else if ( mtrx == StiffnessMatrix ) { //Tangent stiffness matrix set as default
97  this->computeStiffnessMatrix(reducedAnswer, TangentStiffness, tStep);
98  answer.assemble(reducedAnswer, this->giveVelocityDofMask(), this->giveVelocityDofMask());
99  } else if ( mtrx == TangentStiffnessMatrix ) {
100  this->computeStiffnessMatrix(reducedAnswer, TangentStiffness, tStep);
101  answer.assemble(reducedAnswer, this->giveVelocityDofMask(), this->giveVelocityDofMask());
102  } else if ( mtrx == PressureGradientMatrix ) {
103  this->computeGradientMatrix(reducedAnswer, tStep);
104  answer.assemble(reducedAnswer, this->giveVelocityDofMask(), this->givePressureDofMask());
105  } else if ( mtrx == DivergenceMatrix ) {
106  this->computeDivergenceMatrix(reducedAnswer, tStep);
107  answer.assemble(reducedAnswer, this->givePressureDofMask(), this->giveVelocityDofMask());
108  } else if ( mtrx == PressureLaplacianMatrix ) {
109  this->computePressureLaplacianMatrix(reducedAnswer, tStep);
110  answer.assemble(reducedAnswer, this->givePressureDofMask(), this->givePressureDofMask());
111  } else {
112  OOFEM_ERROR("giveCharacteristicMatrix: Unknown Type of characteristic mtrx.");
113  }
114 }
115 
116 
117 void
119  TimeStep *tStep)
120 //
121 // returns characteristics vector of receiver according to requested type
122 //
123 {
124  FloatArray reducedVector;
125  answer.resize(this->computeNumberOfDofs());
126  answer.zero();
127 
128  if ( mtrx == ExternalForcesVector ) {
129  //this->computeForceVector(reducedVector, tStep);
130  //answer.assemble(reducedVector, this->giveVelocityDofMask());
131  } else if ( mtrx == PressureGradientVector ) {
132  FloatMatrix g;
133  this->giveCharacteristicMatrix(g, PressureGradientMatrix, tStep);
134  FloatArray p;
135  this->computeVectorOfPressures(VM_Total, tStep, reducedVector);
136  p.resize(this->computeNumberOfDofs());
137  p.assemble(reducedVector, this->givePressureDofMask());
138  answer.beProductOf(g, p);
139  } else if ( mtrx == MassVelocityVector ) {
140  FloatMatrix m;
141  FloatArray u;
142  this->giveCharacteristicMatrix(m, LumpedMassMatrix, tStep);
143  this->computeVectorOfVelocities(VM_Total, tStep, reducedVector);
144  u.resize(this->computeNumberOfDofs());
145  u.assemble(reducedVector, this->giveVelocityDofMask());
146  answer.beProductOf(m, u);
147  } else if ( mtrx == LaplaceVelocityVector ) {
148  FloatMatrix l;
149  FloatArray u;
150  this->giveCharacteristicMatrix(l, StiffnessMatrix, tStep);
151  this->computeVectorOfVelocities(VM_Total, tStep, reducedVector);
152  u.resize(this->computeNumberOfDofs());
153  u.assemble(reducedVector, this->giveVelocityDofMask());
154  answer.beProductOf(l, u);
155  } else if ( mtrx == MassAuxVelocityVector ) {
156  FloatMatrix m;
157  FloatArray u;
158  this->giveCharacteristicMatrix(m, LumpedMassMatrix, tStep);
159  this->computeVectorOf(VM_Intermediate, tStep, u);
160  answer.beProductOf(m, u);
161  } else if ( mtrx == DivergenceAuxVelocityVector ) {
162  FloatMatrix d;
163  this->giveCharacteristicMatrix(d, DivergenceMatrix, tStep);
164  FloatArray u_star;
165  this->computeVectorOf(VM_Intermediate, tStep, u_star);
166  answer.beProductOf(d, u_star);
167  } else if ( mtrx == DivergenceVelocityVector ) {
168  FloatMatrix d;
169  this->giveCharacteristicMatrix(d, DivergenceMatrix, tStep);
170  FloatArray u;
171  this->computeVectorOfVelocities(VM_Total, tStep, u);
172  answer.beProductOf(d, u);
173  } else if ( mtrx == LumpedMassMatrix ) {
174  this->computeDiagonalMassMtrx(reducedVector, tStep);
175  answer.assemble(reducedVector, this->giveVelocityDofMask());
176  } else {
177  OOFEM_ERROR("giveCharacteristicVector: Unknown Type of characteristic mtrx.");
178  }
179 }
180 
181 
183 {
184  if ( type != ExternalForcesVector ) {
185  answer.clear();
186  return;
187  }
188  answer.resize(this->computeNumberOfDofs());
189  answer.zero();
190  FloatArray reducedVector;
191  this->computeBodyLoadVectorAt(reducedVector, load, tStep, mode);
192  answer.assemble(reducedVector, this->giveVelocityDofMask());
193 }
194 
195 
196 
197 int
199 // no check at the moment
200 {
201  int result = 1;
202 
203  return result;
204 }
205 
206 void
208 {
209  FloatArray stress;
210  // force updating strains & stresses
211  for ( auto &iRule: integrationRulesArray ) {
212  for ( GaussPoint *gp: *iRule ) {
213  computeDeviatoricStress(stress, gp, stepN);
214  }
215  }
216 }
217 
218 void
220 // Performs end-of-step operations.
221 {
222 #ifdef __PARALLEL_MODE
223  fprintf( file, "element %d [%8d] :\n", this->giveNumber(), this->giveGlobalNumber() );
224 #else
225  fprintf(file, "element %d :\n", number);
226 #endif
227  for ( auto &iRule: integrationRulesArray ) {
228  iRule->printOutputAt(file, tStep);
229  }
230 }
231 
232 
233 #ifdef __OOFEG
234 int
236  int node, TimeStep *atTime)
237 {
238  int indx = 1;
239  Node *n = this->giveNode(node);
240 
241  if ( type == IST_Velocity ) {
242  answer.resize( this->giveSpatialDimension() );
243  int dofindx;
244  if ( ( dofindx = n->findDofWithDofId(V_u) ) ) {
245  answer.at(indx++) = n->giveDof(dofindx)->giveUnknown(VM_Total, atTime);
246  }
247 
248  if ( ( dofindx = n->findDofWithDofId(V_v) ) ) {
249  answer.at(indx++) = n->giveDof(dofindx)->giveUnknown(VM_Total, atTime);
250  }
251 
252  if ( ( dofindx = n->findDofWithDofId(V_w) ) ) {
253  answer.at(indx++) = n->giveDof(dofindx)->giveUnknown(VM_Total, atTime);
254  }
255 
256  return 1;
257  } else if ( type == IST_Pressure ) {
258  int dofindx;
259  if ( ( dofindx = n->findDofWithDofId(P_f) ) ) {
260  answer.resize(1);
261  answer.at(1) = n->giveDof(dofindx)->giveUnknown(VM_Total, atTime);
262  return 1;
263  } else {
264  return 0;
265  }
266  } else {
267  return Element :: giveInternalStateAtNode(answer, type, mode, node, atTime);
268  }
269 }
270 
271 #endif
272 
273 // might be needed for OOFEG-plot
274 // int
275 // PFEMElement :: giveIntVarCompFullIndx(IntArray &answer, InternalStateType type)
276 // {
277 // if ( ( type == IST_Velocity ) ) {
278 // IntArray mask;
279 // int indx = 1;
280 // answer.resize(3);
281 // this->giveElementDofIDMask(EID_MomentBalance, mask);
282 // if ( mask.findFirstIndexOf(V_u) ) {
283 // answer.at(1) = indx++;
284 // }
285 //
286 // if ( mask.findFirstIndexOf(V_v) ) {
287 // answer.at(2) = indx++;
288 // }
289 //
290 // if ( mask.findFirstIndexOf(V_w) ) {
291 // answer.at(3) = indx++;
292 // }
293 //
294 // return 1;
295 // } else if ( type == IST_Pressure ) {
296 // answer.resize(1);
297 // answer.at(1) = 1;
298 // return 1;
299 // } else {
300 // return Element :: giveIntVarCompFullIndx(answer, type);
301 // }
302 // }
303 } // end namespace oofem
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
virtual void computeBodyLoadVectorAt(FloatArray &answer, BodyLoad *load, TimeStep *tStep, ValueModeType mode)=0
Computes the load vector due to body load acting on receiver, at given time step. ...
int number
Component number.
Definition: femcmpnn.h:80
virtual void computeDeviatoricStress(FloatArray &answer, GaussPoint *gp, TimeStep *)=0
Computes deviatoric stress vector in given integration point and solution step from given total strai...
virtual void computeDiagonalMassMtrx(FloatArray &answer, TimeStep *)=0
Calculates diagonal mass matrix as vector.
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's local DOFs.
Definition: element.h:424
virtual int checkConsistency()
Performs consistency check.
Definition: pfemelement.C:198
IRResultType initializeFrom(InputRecord *ir)
Initializes receiver acording to object description stored in input record.
Definition: pfemelement.C:69
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
int giveGlobalNumber() const
Definition: element.h:1059
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
Class implementing element body load, acting over whole element volume (e.g., the dead weight)...
Definition: bodyload.h:49
void clear()
Clears receiver (zero size).
Definition: floatarray.h:206
~PFEMElement()
Destructor.
Definition: pfemelement.C:64
void computeLoadVector(FloatArray &answer, BodyLoad *load, CharType type, ValueModeType mode, TimeStep *tStep)
Computes the contribution of the given body load (volumetric).
Definition: pfemelement.C:182
virtual void updateInternalState(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
Definition: pfemelement.C:207
virtual void computeDivergenceMatrix(FloatMatrix &answer, TimeStep *atTime)=0
Calculates the velocity divergence matrix.
virtual int giveInternalStateAtNode(FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *atTime)
Returns internal state variable (like stress,strain) at node of element in Reduced form...
Definition: pfemelement.C:235
virtual void computePressureLaplacianMatrix(FloatMatrix &answer, TimeStep *atTime)=0
Calculates the pressure laplacian matrix.
void giveCharacteristicMatrix(FloatMatrix &answer, CharType, TimeStep *)
Computes characteristic matrix of receiver of requested type in given time step.
Definition: pfemelement.C:82
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: element.h:1158
virtual void computeGradientMatrix(FloatMatrix &answer, TimeStep *atTime)=0
Calculates the pressure gradient matrix.
virtual void computeVectorOfPressures(ValueModeType mode, TimeStep *tStep, FloatArray &pressures)
Definition: fmelement.C:55
PFEMElement(int, Domain *)
Constructor.
Definition: pfemelement.C:58
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: element.C:638
virtual void computeVectorOfVelocities(ValueModeType mode, TimeStep *tStep, FloatArray &velocities)
Definition: fmelement.C:48
#define OOFEM_ERROR(...)
Definition: error.h:61
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
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...
Definition: element.C:1661
Initializes the variable VERBOSE, in order to get a few intermediate messages on screen: beginning an...
Class representing vector of real numbers.
Definition: floatarray.h:82
This abstract class represent a general base element class for fluid dynamic problems.
Definition: fmelement.h:54
void giveCharacteristicVector(FloatArray &answer, CharType, ValueModeType, TimeStep *)
Computes characteristic vector of receiver of requested type in given time step.
Definition: pfemelement.C:118
virtual int giveSpatialDimension()
Returns the element spatial dimension (1, 2, or 3).
Definition: element.C:1347
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
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
CharType
Definition: chartype.h:87
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
std::vector< Dof * >::const_iterator findDofWithDofId(DofIDItem dofID) const
Finds index of DOF with required physical meaning of receiver.
Definition: dofmanager.C:266
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode mode, TimeStep *atTime)=0
Calculates the stiffness matrix.
virtual const IntArray & givePressureDofMask() const =0
Returns mask of pressure Dofs.
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
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
the oofem namespace is to define a context or scope in which all oofem names are defined.
void assemble(const FloatMatrix &src, const IntArray &loc)
Assembles the contribution using localization array into receiver.
Definition: floatmatrix.C:251
Class implementing node in finite element mesh.
Definition: node.h:87
int giveNumber() const
Definition: femcmpnn.h:107
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
virtual const IntArray & giveVelocityDofMask() const =0
Returns mask of velocity Dofs.
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to stream, for given time step.
Definition: pfemelement.C:219
InternalStateMode
Determines the mode of internal variable.
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:30 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011