OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
trplanrot3d.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/PlaneStress/trplanrot3d.h"
36 #include "../sm/CrossSections/structuralcrosssection.h"
37 #include "../sm/Materials/structuralms.h"
38 #include "material.h"
39 #include "node.h"
40 #include "load.h"
41 #include "mathfem.h"
42 #include "classfactory.h"
43 #include "gaussintegrationrule.h"
44 #include "gausspoint.h"
45 #include "fei2dtrlin.h"
46 
47 namespace oofem {
48 REGISTER_Element(TrPlaneStrRot3d);
49 
51 {
52 }
53 
54 
55 void
57 // Returns global coordinates given in global vector
58 // transformed into local coordinate system of the
59 // receiver
60 {
61  // test the parameter
62  if ( global.giveSize() != 3 ) {
63  OOFEM_ERROR("cannot transform coordinates - size mismatch");
64  exit(1);
65  }
66 
67  // first ensure that receiver's GtoLRotationMatrix[3,3] is defined
68  if ( !GtoLRotationMatrix.isNotEmpty() ) {
70  }
71 
72  FloatArray offset = global;
73  offset.subtract( * this->giveNode(1)->giveCoordinates() );
74  answer.beProductOf(GtoLRotationMatrix, offset);
75 }
76 
77 
78 double
80 {
81  double detJ, weight;
82 
83  FloatArray x, y;
84  this->giveNodeCoordinates(x, y);
85  std :: vector< FloatArray > lc = {{x[0], y[0]}, {x[1], y[1]}, {x[2], y[2]}};
86 
87  weight = gp->giveWeight();
89  return detJ * weight * this->giveStructuralCrossSection()->give(CS_Thickness, gp);
90 }
91 
92 
93 void
95 {
96  FloatArray nc1(3), nc2(3), nc3(3);
97 
98  this->giveLocalCoordinates( nc1, * ( this->giveNode(1)->giveCoordinates() ) );
99  this->giveLocalCoordinates( nc2, * ( this->giveNode(2)->giveCoordinates() ) );
100  this->giveLocalCoordinates( nc3, * ( this->giveNode(3)->giveCoordinates() ) );
101 
102  x.resize(3);
103  x.at(1) = nc1.at(1);
104  x.at(2) = nc2.at(1);
105  x.at(3) = nc3.at(1);
106 
107  y.resize(3);
108  y.at(1) = nc1.at(2);
109  y.at(2) = nc2.at(2);
110  y.at(3) = nc3.at(2);
111 
112  //if (z) {
113  // z[0] = nc1->at(3);
114  // z[1] = nc2->at(3);
115  // z[2] = nc3->at(3);
116  //}
117 }
118 
119 
120 void
122 {
123  answer = {D_u, D_v, D_w, R_u, R_v, R_w};
124 }
125 
126 
127 const FloatMatrix *
129 // Returns the rotation matrix of the receiver of the size [3,3]
130 // coords(local) = T * coords(global)
131 //
132 // local coordinate (described by vector triplet e1',e2',e3') is defined as follows:
133 //
134 // e1' : [N2-N1] Ni - means i - th node
135 // help : [N3-N1]
136 // e3' : e1' x help
137 // e2' : e3' x e1'
138 {
139  if ( !GtoLRotationMatrix.isNotEmpty() ) {
140  FloatArray e1, e2, e3, help;
141 
142  // compute e1' = [N2-N1] and help = [N3-N1]
143  e1.beDifferenceOf( * this->giveNode(2)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
144  help.beDifferenceOf( * this->giveNode(3)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
145 
146  // let us normalize e1'
147  e1.normalize();
148 
149  // compute e3' : vector product of e1' x help
150  e3.beVectorProductOf(e1, help);
151  // let us normalize
152  e3.normalize();
153 
154  // now from e3' x e1' compute e2'
155  e2.beVectorProductOf(e3, e1);
156 
157  //
159 
160  for ( int i = 1; i <= 3; i++ ) {
161  GtoLRotationMatrix.at(1, i) = e1.at(i);
162  GtoLRotationMatrix.at(2, i) = e2.at(i);
163  GtoLRotationMatrix.at(3, i) = e3.at(i);
164  }
165  }
166 
167  return &GtoLRotationMatrix;
168 }
169 
170 
171 bool
173 // Returns the rotation matrix of the receiver of the size [9,18]
174 // r(local) = T * r(global)
175 // for one node (r written transposed): {u,v,r3} = T * {u,v,w,r1,r2,r3}
176 {
177  // test if pereviously computed
178  if ( !GtoLRotationMatrix.isNotEmpty() ) {
180  }
181 
182  answer.resize(9, 18);
183  answer.zero();
184 
185  for ( int i = 1; i <= 3; i++ ) {
186  answer.at(1, i) = answer.at(1 + 3, i + 6) = answer.at(1 + 6, i + 12) = GtoLRotationMatrix.at(1, i);
187  answer.at(2, i) = answer.at(2 + 3, i + 6) = answer.at(2 + 6, i + 12) = GtoLRotationMatrix.at(2, i);
188  answer.at(3, i + 3) = answer.at(3 + 3, i + 3 + 6) = answer.at(3 + 6, i + 3 + 12) = GtoLRotationMatrix.at(3, i);
189  }
190 
191  return 1;
192 }
193 
194 
195 void
197 // returns characteristic tensor of the receiver at given gp and tStep
198 // strain vector = (Eps_X, Eps_y, Gamma_xy, Kappa_z)
199 {
200  FloatArray charVect;
202 
203  answer.resize(3, 3);
204  answer.zero();
205 
206  if ( type == LocalForceTensor || type == GlobalForceTensor ) {
207  //this->computeStressVector(charVect, gp, tStep);
208  charVect = ms->giveStressVector();
209 
210  double h = this->giveStructuralCrossSection()->give(CS_Thickness, gp);
211  answer.at(1, 1) = charVect.at(1) * h;
212  answer.at(2, 2) = charVect.at(2) * h;
213  answer.at(1, 2) = charVect.at(3) * h;
214  answer.at(2, 1) = charVect.at(3) * h;
215  } else if ( type == LocalMomentTensor || type == GlobalMomentTensor ) {
216  //this->computeStressVector(charVect, gp, tStep);
217  charVect = ms->giveStressVector();
218 
219  answer.at(3, 3) = charVect.at(4);
220  } else if ( type == LocalStrainTensor || type == GlobalStrainTensor ) {
221  //this->computeStrainVector(charVect, gp, tStep);
222  charVect = ms->giveStrainVector();
223 
224  answer.at(1, 1) = charVect.at(1);
225  answer.at(2, 2) = charVect.at(2);
226  answer.at(1, 2) = charVect.at(3) / 2.;
227  answer.at(2, 1) = charVect.at(3) / 2.;
228  } else if ( type == LocalCurvatureTensor || type == GlobalCurvatureTensor ) {
229  //this->computeStrainVector(charVect, gp, tStep);
230  charVect = ms->giveStrainVector();
231 
232  answer.at(3, 3) = charVect.at(4);
233  } else {
234  OOFEM_ERROR("unsupported tensor mode");
235  exit(1);
236  }
237 
238  if ( type == GlobalForceTensor || type == GlobalMomentTensor ||
239  type == GlobalStrainTensor || type == GlobalCurvatureTensor ) {
242  }
243 }
244 
245 
246 int
248 {
249  FloatMatrix globTensor;
250  CharTensor cht;
251 
252  answer.resize(6);
253 
254  if ( type == IST_CurvatureTensor || type == IST_ShellStrainTensor ) {
255  if ( type == IST_CurvatureTensor ) {
256  cht = GlobalCurvatureTensor;
257  } else {
258  cht = GlobalStrainTensor;
259  }
260 
261  this->giveCharacteristicTensor(globTensor, cht, gp, tStep);
262 
263  answer.at(1) = globTensor.at(1, 1); //xx
264  answer.at(2) = globTensor.at(2, 2); //yy
265  answer.at(3) = globTensor.at(3, 3); //zz
266  answer.at(4) = 2 * globTensor.at(2, 3); //yz
267  answer.at(5) = 2 * globTensor.at(1, 3); //xz
268  answer.at(6) = 2 * globTensor.at(1, 2); //xy
269 
270  return 1;
271  } else if ( type == IST_ShellMomentTensor || type == IST_ShellForceTensor ) {
272  if ( type == IST_ShellMomentTensor ) {
273  cht = GlobalMomentTensor;
274  } else {
275  cht = GlobalForceTensor;
276  }
277 
278  this->giveCharacteristicTensor(globTensor, cht, gp, tStep);
279 
280  answer.at(1) = globTensor.at(1, 1); //xx
281  answer.at(2) = globTensor.at(2, 2); //yy
282  answer.at(3) = globTensor.at(3, 3); //zz
283  answer.at(4) = globTensor.at(2, 3); //yz
284  answer.at(5) = globTensor.at(1, 3); //xz
285  answer.at(6) = globTensor.at(1, 2); //xy
286 
287  return 1;
288  } else {
289  return TrPlaneStrRot :: giveIPValue(answer, gp, type, tStep);
290  }
291 }
292 
293 int
295 // Returns the rotation matrix of the receiver of the size [6,6]
296 // f(local) = T * f(global)
297 {
298  // test if previously computed
299  if ( !GtoLRotationMatrix.isNotEmpty() ) {
301  }
302 
303  answer.resize(6, 6);
304  answer.zero();
305 
306  for ( int i = 1; i <= 3; i++ ) {
307  answer.at(1, i) = answer.at(4, i + 3) = GtoLRotationMatrix.at(1, i);
308  answer.at(2, i) = answer.at(5, i + 3) = GtoLRotationMatrix.at(2, i);
309  answer.at(3, i) = answer.at(6, i + 3) = GtoLRotationMatrix.at(3, i);
310  }
311 
312  return 1;
313 }
314 
315 
316 void
318 // Computes numerically the load vector of the receiver due to the body loads, at tStep.
319 // load is assumed to be in global cs.
320 // load vector is then transformed to coordinate system in each node.
321 // (should be global coordinate system, but there may be defined
322 // different coordinate system in each node)
323 {
324  double dens, dV, load;
325  FloatArray force;
326  FloatMatrix T;
327 
328  if ( ( forLoad->giveBCGeoType() != BodyLoadBGT ) || ( forLoad->giveBCValType() != ForceLoadBVT ) ) {
329  OOFEM_ERROR("unknown load type");
330  }
331 
332  // note: force is assumed to be in global coordinate system; 6 components
333  forLoad->computeComponentArrayAt(force, tStep, mode);
334  // get it in local c.s.
336  force.rotatedWith(T, 'n');
337 
338  if ( force.giveSize() ) {
339  GaussIntegrationRule iRule (1, this, 1, 1);
340  iRule.SetUpPointsOnTriangle(1, _Unknown);
341  GaussPoint *gp = iRule.getIntegrationPoint(0);
342 
343  dens = this->giveStructuralCrossSection()->give('d', gp);
344  dV = this->computeVolumeAround(gp);// * this->giveCrossSection()->give(CS_Thickness, gp);
345 
346  answer.resize(9);
347  answer.zero();
348 
349  load = force.at(1) * dens * dV / 3.0;
350  answer.at(1) = load;
351  answer.at(4) = load;
352  answer.at(7) = load;
353 
354  load = force.at(2) * dens * dV / 3.0;
355  answer.at(2) = load;
356  answer.at(5) = load;
357  answer.at(8) = load;
358 
359  load = force.at(6) * dens * dV / 3.0;
360  answer.at(3) = load;
361  answer.at(6) = load;
362  answer.at(9) = load;
363 
364  // transform result from global cs to local element cs.
365  //if ( this->computeGtoLRotationMatrix(T) ) {
366  // answer.rotatedWith(T, 'n');
367  //}
368  } else {
369  answer.clear(); // nil resultant
370  }
371 }
372 
373 void
375 {
376  FloatMatrix ne;
377  this->computeNmatrixAt(sgp->giveNaturalCoordinates(), ne);
378 
379  answer.resize(6, 18);
380  answer.zero();
381  int ri[] = {
382  0, 1, 5
383  };
384  int ci[] = {
385  0, 1, 5, 6, 7, 11, 12, 13, 17
386  };
387 
388  for ( int i = 0; i < 3; i++ ) {
389  for ( int j = 0; j < 9; j++ ) {
390  answer(ri [ i ], ci [ j ]) = ne(i, j);
391  }
392  }
393 }
394 
395 void
397 {
398  answer.resize(18);
399  answer.zero();
400  if ( iSurf == 1 ) {
401  answer.at(1) = 1; // node 1
402  answer.at(2) = 2;
403  answer.at(6) = 3;
404 
405  answer.at(7) = 4; // node 2
406  answer.at(8) = 5;
407  answer.at(12) = 6;
408 
409  answer.at(13) = 7; // node 3
410  answer.at(14) = 8;
411  answer.at(18) = 9;
412  } else {
413  OOFEM_ERROR("wrong surface number");
414  }
415 }
416 
419 {
420  IntegrationRule *iRule = new GaussIntegrationRule(1, this, 1, 1);
421  int npoints = iRule->getRequiredNumberOfIntegrationPoints(_Triangle, approxOrder);
422  iRule->SetUpPointsOnTriangle(npoints, _Unknown);
423  return iRule;
424 }
425 
426 double
428 {
429  return this->computeVolumeAround(gp);
430 }
431 
432 
433 int
435 {
436  return 0;
437 }
438 
439 
440 void
442 // Performs end-of-step operations.
443 {
444  FloatArray v;
445 
446  fprintf( file, "element %d (%8d) :\n", this->giveLabel(), this->giveNumber() );
447 
448  for ( int i = 0; i < (int)integrationRulesArray.size(); i++ ) {
449  for ( GaussPoint *gp: *integrationRulesArray [ i ] ) {
450 
451  fprintf( file, " GP %2d.%-2d :", i + 1, gp->giveNumber() );
452 
453  this->giveIPValue(v, gp, IST_ShellStrainTensor, tStep);
454  fprintf(file, " strains ");
455  for ( auto &val : v ) fprintf(file, " %.4e", val);
456 
457  this->giveIPValue(v, gp, IST_CurvatureTensor, tStep);
458  fprintf(file, "\n curvatures ");
459  for ( auto &val : v ) fprintf(file, " %.4e", val);
460 
461  // Forces - Moments
462  this->giveIPValue(v, gp, IST_ShellForceTensor, tStep);
463  fprintf(file, "\n stresses ");
464  for ( auto &val : v ) fprintf(file, " %.4e", val);
465 
466  this->giveIPValue(v, gp, IST_ShellMomentTensor, tStep);
467  fprintf(file, "\n moments ");
468  for ( auto &val : v ) fprintf(file, " %.4e", val);
469 
470  fprintf(file, "\n");
471  }
472  }
473 }
474 } // end namespace oofem
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
virtual void computeBodyLoadVectorAt(FloatArray &answer, Load *forLoad, TimeStep *tStep, ValueModeType mode)
Computes the load vector due to body load acting on receiver, at given time step. ...
Definition: trplanrot3d.C:317
void subtract(const FloatArray &src)
Subtracts array src to receiver.
Definition: floatarray.C:258
const FloatMatrix * computeGtoLRotationMatrix()
Definition: trplanrot3d.C:128
TrPlaneStrRot3d(int n, Domain *d)
Definition: trplanrot3d.C:50
Class implements an triangular three-node plane- stress elasticity finite element with independent ro...
Definition: trplanrot.h:56
void beVectorProductOf(const FloatArray &v1, const FloatArray &v2)
Computes vector product (or cross product) of vectors given as parameters, , and stores the result in...
Definition: floatarray.C:415
Class and object Domain.
Definition: domain.h:115
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: trplanrot3d.C:79
Wrapper around cell with vertex coordinates stored in FloatArray**.
Definition: feinterpol.h:115
virtual IntegrationRule * GetSurfaceIntegrationRule(int iSurf)
Definition: trplanrot3d.C:418
void zero()
Sets all component to zero.
Definition: intarray.C:52
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
virtual double computeSurfaceVolumeAround(GaussPoint *gp, int iSurf)
Computes volume related to integration point on local surface.
Definition: trplanrot3d.C:427
FloatMatrix GtoLRotationMatrix
Transformation Matrix form GtoL(3,3) is stored at the element level for computation efficiency...
Definition: trplanrot3d.h:75
virtual void computeComponentArrayAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
Computes boundary condition value - its components values at given time.
Definition: load.C:82
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
void rotatedWith(FloatMatrix &r, char mode)
Returns the receiver a rotated according the change-of-base matrix r.
Definition: floatarray.C:799
Abstract base class representing integration rule.
Distributed body load.
Definition: bcgeomtype.h:43
void beDifferenceOf(const FloatArray &a, const FloatArray &b)
Sets receiver to be a - b.
Definition: floatarray.C:341
virtual void giveSurfaceDofMapping(IntArray &answer, int iSurf) const
Assembles surface dof mapping mask, which provides mapping between surface local DOFs and "global" el...
Definition: trplanrot3d.C:396
bool isNotEmpty() const
Tests for empty matrix.
Definition: floatmatrix.h:162
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: trplanrot3d.C:121
#define OOFEM_ERROR(...)
Definition: error.h:61
void giveLocalCoordinates(FloatArray &answer, FloatArray &global)
Definition: trplanrot3d.C:56
REGISTER_Element(LSpace)
virtual double giveWeight()
Returns integration weight of receiver.
Definition: gausspoint.h:181
StructuralCrossSection * giveStructuralCrossSection()
Helper function which returns the structural cross-section for the element.
void rotatedWith(const FloatMatrix &r, char mode= 'n')
Returns the receiver &#39;a&#39; transformed using give transformation matrix r.
Definition: floatmatrix.C:1557
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: trplanrot3d.C:247
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to stream, for given time step.
Definition: trplanrot3d.C:441
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
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
virtual void giveNodeCoordinates(FloatArray &x, FloatArray &y)
Definition: trplanrot3d.C:94
virtual int computeLoadGToLRotationMtrx(FloatMatrix &answer)
Returns transformation matrix from global coordinate system to local element coordinate system for el...
Definition: trplanrot3d.C:294
void giveCharacteristicTensor(FloatMatrix &answer, CharTensor type, GaussPoint *gp, TimeStep *tStep)
Definition: trplanrot3d.C:196
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
GaussPoint * getIntegrationPoint(int n)
Access particular integration point of receiver.
virtual double giveTransformationJacobian(const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation.
Definition: fei2dtrlin.C:147
const FloatArray & giveStressVector() const
Returns the const pointer to receiver&#39;s stress vector.
Definition: structuralms.h:107
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
Definition: gausspoint.h:205
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: trplanrot.C:241
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
virtual void computeSurfaceNMatrixAt(FloatMatrix &answer, int iSurf, GaussPoint *gp)
Definition: trplanrot3d.C:374
virtual int getRequiredNumberOfIntegrationPoints(integrationDomain dType, int approxOrder)
Abstract service.
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual bcGeomType giveBCGeoType() const
Returns geometry character of boundary condition.
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 int SetUpPointsOnTriangle(int, MaterialMode mode)
Sets up receiver&#39;s integration points on triangular (area coords) integration domain.
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
Load is base abstract class for all loads.
Definition: load.h:61
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
int giveLabel() const
Definition: element.h:1055
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual bcValType giveBCValType() const
Returns receiver load type.
int giveNumber() const
Definition: femcmpnn.h:107
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
const FloatArray & giveStrainVector() const
Returns the const pointer to receiver&#39;s strain vector.
Definition: structuralms.h:105
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: trplanrot.C:674
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
virtual int SetUpPointsOnTriangle(int nPoints, MaterialMode mode)
Sets up receiver&#39;s integration points on triangular (area coords) integration domain.
Class representing Gaussian-quadrature integration rule.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
static FEI2dTrLin interp
Definition: trplanstrss.h:67
virtual int computeLoadLSToLRotationMatrix(FloatMatrix &answer, int iSurf, GaussPoint *gp)
Returns transformation matrix from local surface c.s to element local coordinate system of load vecto...
Definition: trplanrot3d.C:434

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