OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
interfaceelement1d.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 "interfaceelement1d.h"
36 #include "domain.h"
37 #include "node.h"
38 #include "gaussintegrationrule.h"
39 #include "floatmatrix.h"
40 #include "floatarray.h"
41 #include "intarray.h"
42 #include "mathfem.h"
43 #include "feinterpol.h"
44 #include "../sm/CrossSections/structuralinterfacecrosssection.h"
45 #include "classfactory.h"
46 
47 #ifdef __OOFEG
48  #include "oofeggraphiccontext.h"
49 
50  #include <Emarkwd3d.h>
51 #endif
52 
53 namespace oofem {
54 REGISTER_Element(InterfaceElem1d);
55 
57  StructuralElement(n, aDomain)
58 {
59  numberOfDofMans = 2;
60  referenceNode = 0;
61  normal.resize(3);
62  normal.zero();
63 }
64 
65 
66 void
68 {
70  case 1:
71  this->mode = ie1d_1d;
72  break;
73  case 2:
74  this->mode = ie1d_2d;
75  break;
76  case 3:
77  this->mode = ie1d_3d;
78  break;
79  default:
80  OOFEM_ERROR("Unsupported domain type")
81  }
82 }
83 
84 
85 void
87 {
88  setCoordMode();
89  switch ( mode ) {
90  case ie1d_1d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->giveEngTraction_1d(answer, gp, strain, tStep); return;
91  case ie1d_2d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->giveEngTraction_2d(answer, gp, strain, tStep); return;
92  case ie1d_3d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->giveEngTraction_3d(answer, gp, strain, tStep); return;
93  }
94 }
95 
96 
97 void
99 {
100  setCoordMode();
101  switch ( mode ) {
102  case ie1d_1d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->give1dStiffnessMatrix_Eng(answer, rMode, gp, tStep); return;
103  case ie1d_2d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->give2dStiffnessMatrix_Eng(answer, rMode, gp, tStep); return;
104  case ie1d_3d: static_cast< StructuralInterfaceCrossSection* >(this->giveCrossSection())->give3dStiffnessMatrix_Eng(answer, rMode, gp, tStep); return;
105  }
106 }
107 
108 
109 
112 {
113  setCoordMode();
114  switch ( mode ) {
115  case ie1d_1d: return _1dInterface;
116 
117  case ie1d_2d: return _2dInterface;
118 
119  case ie1d_3d: return _3dInterface;
120 
121  default: OOFEM_ERROR("Unsupported coord mode");
122  }
123  return _1dInterface; // to make the compiler happy
124 }
125 
126 
127 void
129 // Returns the lumped mass matrix (which is zero matrix) of the receiver. This expression is
130 // valid in both local and global axes.
131 {
132  int ndofs = this->computeNumberOfDofs();
133  answer.resize(ndofs, ndofs);
134  answer.zero();
135 }
136 
137 
138 void
140 //
141 // Returns linear part of geometrical equations of the receiver at gp.
142 // Returns the linear part of the B matrix
143 //
144 {
145  setCoordMode();
146  //double area = this->giveCrossSection()->give(CS_Area);
147  double area = 1.;
149  FloatMatrix bLoc, lcs;
151  switch ( mode ) {
152  case ie1d_1d:
153  bLoc.resize(1, 2);
154  bLoc.at(1, 1) = -1.0;
155  bLoc.at(1, 2) = +1.0;
156  break;
157  case ie1d_2d:
158  bLoc.resize(2, 4);
159  bLoc.zero();
160  bLoc.at(1, 1) = -1.0;
161  bLoc.at(1, 3) = +1.0;
162  bLoc.at(2, 2) = -1.0;
163  bLoc.at(2, 4) = +1.0;
164  break;
165  case ie1d_3d:
166  bLoc.resize(3, 6);
167  bLoc.zero();
168  bLoc.at(1, 1) = -1.0;
169  bLoc.at(1, 4) = +1.0;
170  bLoc.at(2, 2) = -1.0;
171  bLoc.at(2, 5) = +1.0;
172  bLoc.at(3, 3) = -1.0;
173  bLoc.at(3, 6) = +1.0;
174  break;
175  default:
176  OOFEM_ERROR("unsupported mode");
177  }
178 
179  bLoc.times(area);
180  answer.beProductOf(lcs, bLoc);
181 }
182 
183 
184 void
186 //
187 // Computes unit vectors of local coordinate system, stored by rows.
188 //
189 {
190  setCoordMode();
191  switch ( mode ) {
192  case ie1d_1d:
193  lcs.resize(1, 1);
194  lcs.at(1, 1) = 1.;
195  return;
196 
197  case ie1d_2d:
198  lcs.resize(2, 2);
199  lcs.at(1, 1) = normal.at(1);
200  lcs.at(1, 2) = normal.at(2);
201  lcs.at(2, 1) = -normal.at(2);
202  lcs.at(2, 2) = normal.at(1);
203  return;
204 
205  case ie1d_3d:
206  {
207  FloatArray ly(3), lz(3);
208  normal.normalize();
209  ly.zero();
210  if ( fabs( normal.at(1) ) > fabs( normal.at(2) ) ) {
211  ly.at(2) = 1.;
212  } else {
213  ly.at(1) = 1.;
214  }
215 
216  lz.beVectorProductOf(normal, ly);
217  lz.normalize();
218  ly.beVectorProductOf(lz, normal);
219  ly.normalize();
220 
221  lcs.resize(3, 3);
222  int i;
223  for ( i = 1; i <= 3; i++ ) {
224  lcs.at(1, i) = normal.at(i);
225  lcs.at(2, i) = ly.at(i);
226  lcs.at(3, i) = lz.at(i);
227  }
228 
229  return;
230  }
231 
232  default:
233  OOFEM_ERROR("unsupported mode");
234  }
235 }
236 
237 
238 void
240 // Sets up the array of Gauss Points of the receiver.
241 {
242  if ( integrationRulesArray.size() == 0 ) {
243  integrationRulesArray.resize( 1 );
244  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 2) );
245  integrationRulesArray [ 0 ]->SetUpPointsOnLine(1, this->giveMaterialMode() );
246  }
247 }
248 
249 
250 int
252 {
253  answer = *this->giveNode(1)->giveCoordinates();
254 
255  return 1;
256 }
257 
258 
259 bool
261 {
262  OOFEM_ERROR("Not implemented");
263  return false;
264 }
265 
266 double
268 // Returns the length of the receiver. This method is valid only if 1
269 // Gauss point is used.
270 {
271  return 1.0; //MUST be set to 1.0
272 }
273 
274 
277 {
278  IRResultType result; // Required by IR_GIVE_FIELD macro
279 
281  if ( result != IRRT_OK ) {
282  return result;
283  }
284 
287  if ( referenceNode == 0 && normal.at(1) == 0 && normal.at(2) == 0 && normal.at(1) == 0 && normal.at(3) == 0 ) {
288  OOFEM_ERROR("wrong reference node or normal specified");
289  }
292  } else {
293  switch ( domain->giveNumberOfSpatialDimensions() ) {
294  case 1:
295  this->dofids = IntArray{D_u};
296  break;
297  case 2:
298  this->dofids = {D_u, D_v};
299  break;
300  case 3:
301  this->dofids = {D_u, D_v, D_w};
302  break;
303  default:
304  OOFEM_ERROR("Unsupported domain type")
305  }
306  }
307 
308  this->computeLocalSlipDir(normal);
309  return IRRT_OK;
310 }
311 
312 
313 int
315 {
316  setCoordMode();
317  switch ( mode ) {
318  case ie1d_1d:
319  return 2;
320 
321  case ie1d_2d:
322  return 4;
323 
324  case ie1d_3d:
325  return 6;
326 
327  default:
328  OOFEM_ERROR("unsupported mode");
329  }
330 
331  return 0; // to suppress compiler warning
332 }
333 
334 
335 void
337 {
338  answer = this->dofids;
339 }
340 
341 
342 void
344 {
345  normal.resizeWithValues(3);
346  if ( this->referenceNode ) {
347  // normal
348  normal.at(1) = domain->giveNode(this->referenceNode)->giveCoordinate(1) - this->giveNode(1)->giveCoordinate(1);
349  normal.at(2) = domain->giveNode(this->referenceNode)->giveCoordinate(2) - this->giveNode(1)->giveCoordinate(2);
350  normal.at(3) = domain->giveNode(this->referenceNode)->giveCoordinate(3) - this->giveNode(1)->giveCoordinate(3);
351  } else {
352  if ( normal.at(1) == 0 && normal.at(2) == 0 && normal.at(3) == 0 ) {
353  OOFEM_ERROR("normal is not defined (referenceNode=0,normal=(0,0,0))");
354  }
355  }
356 
357  normal.normalize();
358 }
359 
360 
361 #ifdef __OOFEG
363 {
364  GraphicObj *go;
365  // if (!go) { // create new one
366  WCRec p [ 1 ]; /* poin */
367  if ( !gc.testElementGraphicActivity(this) ) {
368  return;
369  }
370 
371  EASValsSetColor( gc.getElementColor() );
372  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
373  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
374  EASValsSetColor( gc.getDeformedElementColor() );
375  p [ 0 ].x = ( FPNum ) ( this->giveNode(1)->giveCoordinate(1) );
376  p [ 0 ].y = ( FPNum ) ( this->giveNode(1)->giveCoordinate(2) );
377  p [ 0 ].z = ( FPNum ) ( this->giveNode(1)->giveCoordinate(3) );
378 
379  EASValsSetMType(CIRCLE_MARKER);
380  go = CreateMarker3D(p);
381  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
382  EMAddGraphicsToModel(ESIModel(), go);
383 }
384 
385 
387 {
388  GraphicObj *go;
389  // if (!go) { // create new one
390  WCRec p [ 1 ]; /* poin */
391  if ( !gc.testElementGraphicActivity(this) ) {
392  return;
393  }
394 
395  double defScale = gc.getDefScale();
396 
397  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
398  EASValsSetColor( gc.getDeformedElementColor() );
399  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
400  p [ 0 ].x = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale) +
401  this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale) );
402  p [ 0 ].y = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale) +
403  this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale) );
404  p [ 0 ].z = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(3, tStep, defScale) +
405  this->giveNode(2)->giveUpdatedCoordinate(3, tStep, defScale) );
406 
407  EASValsSetMType(CIRCLE_MARKER);
408  go = CreateMarker3D(p);
409  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
410  EMAddGraphicsToModel(ESIModel(), go);
411 }
412 
413 
415 {
416  int indx, result = 0;
417  FloatArray gcoord(3), v1;
418  WCRec p [ 1 ];
419  GraphicObj *go;
420  double val [ 1 ];
421 
422  if ( !gc.testElementGraphicActivity(this) ) {
423  return;
424  }
425 
426  if ( gc.getInternalVarsDefGeoFlag() ) {
427  double defScale = gc.getDefScale();
428  p [ 0 ].x = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale) +
429  this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale) );
430  p [ 0 ].y = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale) +
431  this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale) );
432  p [ 0 ].z = ( FPNum ) 0.5 * ( this->giveNode(1)->giveUpdatedCoordinate(3, tStep, defScale) +
433  this->giveNode(2)->giveUpdatedCoordinate(3, tStep, defScale) );
434  } else {
435  p [ 0 ].x = ( FPNum ) ( this->giveNode(1)->giveCoordinate(1) );
436  p [ 0 ].y = ( FPNum ) ( this->giveNode(1)->giveCoordinate(2) );
437  p [ 0 ].z = ( FPNum ) ( this->giveNode(1)->giveCoordinate(3) );
438  }
439 
441  //result += giveIPValue(v1, iRule->getIntegrationPoint(0), gc.giveIntVarType(), tStep);
442 
443 
444  for ( GaussPoint *gp: *this->giveDefaultIntegrationRulePtr() ) {
445  result = 0;
446  result += giveIPValue(v1, gp, gc.giveIntVarType(), tStep);
447  if ( result != 1 ) {
448  continue;
449  }
450 
451  indx = gc.giveIntVarIndx();
452 
453  val [ 0 ] = v1.at(indx);
454  gc.updateFringeTableMinMax(val, 1);
455 
456  EASValsSetLayer(OOFEG_VARPLOT_PATTERN_LAYER);
457  EASValsSetMType(FILLED_CIRCLE_MARKER);
458  go = CreateMarkerWD3D(p, val [ 0 ]);
459  EGWithMaskChangeAttributes(LAYER_MASK | FILL_MASK | MTYPE_MASK, go);
460  EMAddGraphicsToModel(ESIModel(), go);
461  //}
462  }
463 }
464 
465 #endif
466 } // 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...
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
void computeLocalSlipDir(FloatArray &normal)
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 IntegrationRule * giveDefaultIntegrationRulePtr()
Access method for default integration rule.
Definition: element.h:822
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
#define OOFEG_RAW_GEOMETRY_LAYER
enum oofem::InterfaceElem1d::cmode mode
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
#define _IFT_InterfaceElem1d_dofIDs
virtual bool hasField(InputFieldType id)=0
Returns true if record contains field identified by idString keyword.
virtual double giveCoordinate(int i)
Definition: node.C:82
int giveNumberOfSpatialDimensions()
Returns number of spatial dimensions.
Definition: domain.C:1067
MaterialMode
Type representing material mode of integration point.
Definition: materialmode.h:89
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 computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
InternalStateType giveIntVarType()
Abstract base class for all "structural" finite elements.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
#define _IFT_InterfaceElem1d_normal
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
virtual MaterialMode giveMaterialMode()
Returns material mode for receiver integration points.
#define OOFEM_ERROR(...)
Definition: error.h:61
REGISTER_Element(LSpace)
UnknownType
Type representing particular unknown (its physical meaning).
Definition: unknowntype.h:55
#define _IFT_InterfaceElem1d_refnode
void times(double f)
Multiplies receiver by factor f.
Definition: floatmatrix.C:1594
void evaluateLocalCoordinateSystem(FloatMatrix &)
virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords)
Computes the global coordinates from given element&#39;s local coordinates.
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
virtual double giveUpdatedCoordinate(int ic, TimeStep *tStep, double scale=1.)
Returns updated ic-th coordinate of receiver.
Definition: node.C:245
void resizeWithValues(int s, int allocChunk=0)
Checks size of receiver towards requested bounds.
Definition: floatarray.C:615
Base class for all structural interface cross section models.
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
double at(int i, int j) const
Coefficient access function.
Definition: floatmatrix.h:176
virtual void giveDofManDofIDMask(int inode, IntArray &answer) const
Returns dofmanager dof mask for node.
virtual bool computeLocalCoordinates(FloatArray &answer, const FloatArray &gcoords)
Computes the element local coordinates from given global coordinates.
InterfaceElem1d(int n, Domain *d)
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
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 zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
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 computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)
Computes the stress vector of receiver at given integration point, at time step tStep.
virtual int computeNumberOfDofs()
Computes or simply returns total number of element&#39;s local DOFs.
void updateFringeTableMinMax(double *s, int size)
#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
Node * giveNode(int n)
Service for accessing particular domain node.
Definition: domain.h:371
the oofem namespace is to define a context or scope in which all oofem names are defined.
#define IR_GIVE_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:69
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
#define OOFEG_VARPLOT_PATTERN_LAYER
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
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
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:29 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011