OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
cct.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/Plates/cct.h"
36 #include "../sm/Materials/structuralms.h"
37 #include "../sm/CrossSections/structuralcrosssection.h"
38 #include "fei2dtrlin.h"
39 #include "node.h"
40 #include "material.h"
41 #include "crosssection.h"
42 #include "gausspoint.h"
43 #include "gaussintegrationrule.h"
44 #include "floatmatrix.h"
45 #include "floatarray.h"
46 #include "intarray.h"
47 #include "load.h"
48 #include "mathfem.h"
49 #include "classfactory.h"
50 
51 #ifdef __OOFEG
52  #include "oofeggraphiccontext.h"
53 #endif
54 
55 namespace oofem {
56 REGISTER_Element(CCTPlate);
57 
58 FEI2dTrLin CCTPlate :: interp_lin(1, 2);
59 //FEI2dTrRot CCTPlate :: interp_rot(1, 2);
60 
61 CCTPlate :: CCTPlate(int n, Domain *aDomain) :
62  NLStructuralElement(n, aDomain),
65 {
66  numberOfDofMans = 3;
68  area = 0;
69 }
70 
71 
74 
75 
78 {
79  if ( id == D_w ) {
80  return & interp_lin;
81  } else {
82  return NULL; //&interp_rot;
83  }
84 }
85 
86 
87 void
89 // Sets up the array containing the four Gauss points of the receiver.
90 {
91  if ( integrationRulesArray.size() == 0 ) {
92  integrationRulesArray.resize( 1 );
93  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 5) );
95  }
96 }
97 
98 
99 void
101 // Computes numerically the load vector of the receiver due to the body loads, at tStep.
102 // load is assumed to be in global cs.
103 // load vector is then transformed to coordinate system in each node.
104 // (should be global coordinate system, but there may be defined
105 // different coordinate system in each node)
106 {
107  FloatArray force;
108  FloatMatrix T;
109 
110  if ( ( forLoad->giveBCGeoType() != BodyLoadBGT ) || ( forLoad->giveBCValType() != ForceLoadBVT ) ) {
111  OOFEM_ERROR("unknown load type");
112  }
113 
114  GaussIntegrationRule irule(1, this, 1, 5);
115  irule.SetUpPointsOnTriangle(1, _2dPlate);
116 
117  // note: force is assumed to be in global coordinate system.
118  forLoad->computeComponentArrayAt(force, tStep, mode);
119 
120  if ( force.giveSize() ) {
121  GaussPoint *gp = irule.getIntegrationPoint(0);
122  // constant density and thickness assumed
123  double dens = this->giveStructuralCrossSection()->give('d', gp);
124  double dV = this->computeVolumeAround(gp) * this->giveCrossSection()->give(CS_Thickness, gp);
125 
126  answer.resize(9);
127  answer.zero();
128 
129  double load = force.at(1) * dens * dV / 3.0;
130  answer.at(1) = load;
131  answer.at(4) = load;
132  answer.at(7) = load;
133 
134  // transform result from global cs to local element cs.
135  if ( this->computeGtoLRotationMatrix(T) ) {
136  answer.rotatedWith(T, 'n');
137  }
138  } else {
139  answer.clear(); // nil resultant
140  }
141 }
142 
143 
144 void
146 // Returns the [5x9] strain-displacement matrix {B} of the receiver,
147 // evaluated at gp.
148 {
149  // get node coordinates
150  double x1, x2, x3, y1, y2, y3, z1, z2, z3;
151  this->giveNodeCoordinates(x1, x2, x3, y1, y2, y3, z1, z2, z3);
152 
153  //
154  double area, b1, b2, b3, c1, c2, c3, l1, l2, l3;
155 
156  b1 = y2 - y3;
157  b2 = y3 - y1;
158  b3 = y1 - y2;
159 
160  c1 = x3 - x2;
161  c2 = x1 - x3;
162  c3 = x2 - x1;
163 
164  l1 = 1. / 3.;
165  l2 = 1. / 3.;
166  l3 = 1. / 3.;
167 
168  area = this->computeArea();
169 
170  answer.resize(5, 9);
171  answer.zero();
172 
173  answer.at(1, 3) = b1;
174  answer.at(1, 6) = b2;
175  answer.at(1, 9) = b3;
176 
177  answer.at(2, 2) = -c1;
178  answer.at(2, 5) = -c2;
179  answer.at(2, 8) = -c3;
180 
181  answer.at(3, 2) = -b1;
182  answer.at(3, 3) = c1;
183  answer.at(3, 5) = -b2;
184  answer.at(3, 6) = c2;
185  answer.at(3, 8) = -b3;
186  answer.at(3, 9) = c3;
187 
188  answer.at(4, 1) = b1;
189  answer.at(4, 2) = ( -b1 * b3 * l2 + b1 * b2 * l3 ) * 0.5;
190  answer.at(4, 3) = ( -b1 * c3 * l2 - b2 * c3 * l1 + b3 * c2 * l1 + b1 * c2 * l3 ) * 0.5 + l1 * 2. * area;
191  answer.at(4, 4) = b2;
192  answer.at(4, 5) = ( -b2 * b1 * l3 + b2 * b3 * l1 ) * 0.5;
193  answer.at(4, 6) = ( -b2 * c1 * l3 - b3 * c1 * l2 + b1 * c3 * l2 + b2 * c3 * l1 ) * 0.5 + l2 * 2. * area;
194  answer.at(4, 7) = b3;
195  answer.at(4, 8) = ( -b3 * b2 * l1 + b3 * b1 * l2 ) * 0.5;
196  answer.at(4, 9) = ( -b3 * c2 * l1 - b1 * c2 * l3 + b2 * c1 * l3 + b3 * c1 * l2 ) * 0.5 + l3 * 2. * area;
197 
198  answer.at(5, 1) = c1;
199  answer.at(5, 2) = ( -b3 * c1 * l2 - b3 * c2 * l1 + b2 * c3 * l1 + b2 * c1 * l3 ) * 0.5 - l1 * 2. * area;
200  answer.at(5, 3) = ( -c1 * c3 * l2 + c1 * c2 * l3 ) * 0.5;
201  answer.at(5, 4) = c2;
202  answer.at(5, 5) = ( -b1 * c2 * l3 - b1 * c3 * l2 + b3 * c1 * l2 + b3 * c2 * l1 ) * 0.5 - l2 * 2. * area;
203  answer.at(5, 6) = ( -c2 * c1 * l3 + c2 * c3 * l1 ) * 0.5;
204  answer.at(5, 7) = c3;
205  answer.at(5, 8) = ( -b2 * c3 * l1 - b2 * c1 * l3 + b1 * c2 * l3 + b1 * c3 * l2 ) * 0.5 - l3 * 2. * area;
206  answer.at(5, 9) = ( -c3 * c2 * l1 + c3 * c1 * l2 ) * 0.5;
207 
208  answer.times( 1. / ( 2. * area ) );
209 }
210 
211 
212 void
214 // Returns the [3x9] displacement interpolation matrix {N} of the receiver,
215 // evaluated at gp.
216 {
217  // get node coordinates
218  double x1, x2, x3, y1, y2, y3, z1, z2, z3;
219  this->giveNodeCoordinates(x1, x2, x3, y1, y2, y3, z1, z2, z3);
220 
221  //
222  double l1, l2, l3, b1, b2, b3, c1, c2, c3;
223 
224  b1 = y2 - y3;
225  b2 = y3 - y1;
226  b3 = y1 - y2;
227 
228  c1 = x3 - x2;
229  c2 = x1 - x3;
230  c3 = x2 - x1;
231 
232  l1 = iLocCoord.at(1);
233  l2 = iLocCoord.at(2);
234  l3 = 1.0 - l1 - l2;
235 
236  //
237  answer.resize(3, 9);
238  answer.zero();
239 
240  answer.at(1, 1) = l1;
241  answer.at(1, 2) = l1 * ( l2 * b3 - l3 * b2 ) * 0.5;
242  answer.at(1, 3) = l1 * ( l2 * c3 - l3 * c2 ) * 0.5;
243  answer.at(1, 4) = l2;
244  answer.at(1, 5) = l2 * ( l3 * b1 - l1 * b3 ) * 0.5;
245  answer.at(1, 6) = l2 * ( l3 * c1 - l1 * c3 ) * 0.5;
246  answer.at(1, 7) = l3;
247  answer.at(1, 8) = l3 * ( l1 * b2 - l2 * b1 ) * 0.5;
248  answer.at(1, 9) = l3 * ( l1 * c2 - l2 * c1 ) * 0.5;
249 
250  answer.at(2, 2) = l1;
251  answer.at(2, 5) = l2;
252  answer.at(2, 8) = l3;
253 
254  answer.at(3, 3) = l1;
255  answer.at(3, 6) = l2;
256  answer.at(3, 9) = l3;
257 }
258 
259 
260 void
262 {
263  this->giveStructuralCrossSection()->giveGeneralizedStress_Plate(answer, gp, strain, tStep);
264 }
265 
266 
267 void
269 {
270  this->giveStructuralCrossSection()->give2dPlateStiffMtrx(answer, rMode, gp, tStep);
271 }
272 
273 
274 void
275 CCTPlate :: giveNodeCoordinates(double &x1, double &x2, double &x3,
276  double &y1, double &y2, double &y3,
277  double &z1, double &z2, double &z3)
278 {
279  FloatArray *nc1, *nc2, *nc3;
280  nc1 = this->giveNode(1)->giveCoordinates();
281  nc2 = this->giveNode(2)->giveCoordinates();
282  nc3 = this->giveNode(3)->giveCoordinates();
283 
284  x1 = nc1->at(1);
285  x2 = nc2->at(1);
286  x3 = nc3->at(1);
287 
288  y1 = nc1->at(2);
289  y2 = nc2->at(2);
290  y3 = nc3->at(2);
291 
292  z1 = nc1->at(3);
293  z2 = nc2->at(3);
294  z3 = nc3->at(3);
295 
296 }
297 
298 double
300 // returns the area occupied by the receiver
301 {
302  // get node coordinates
303  double x1, x2, x3, y1, y2, y3, z1, z2, z3;
304  this->giveNodeCoordinates(x1, x2, x3, y1, y2, y3, z1, z2, z3);
305 
306  if (area > 0) return area; // check if previously computed
307 
308  return (area = 0.5*(x2*y3+x1*y2+y1*x3-x2*y1-x3*y2-x1*y3)) ;
309 
310 }
311 
314 {
317 }
318 
319 
320 void
321 CCTPlate :: giveDofManDofIDMask(int inode, IntArray &answer) const
322 {
323  answer = {D_w, R_u, R_v};
324 }
325 
326 
327 void
329 // returns normal vector to midPlane in GaussPoinr gp of receiver
330 {
331  FloatArray u, v;
332  u.beDifferenceOf( * this->giveNode(2)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
333  v.beDifferenceOf( * this->giveNode(3)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
334 
335  answer.beVectorProductOf(u, v);
336  answer.normalize();
337 }
338 
339 
340 double
342 // returns receiver's characteristic length (for crack band models)
343 // for a crack formed in the plane with normal normalToCrackPlane.
344 {
345  return this->giveCharacteristicLengthForPlaneElements(normalToCrackPlane);
346 }
347 
348 
349 double
351 // Returns the portion of the receiver which is attached to gp.
352 {
353  double detJ, weight;
354 
355  std :: vector< FloatArray > lc = {FloatArray(3), FloatArray(3), FloatArray(3)};
356  this->giveNodeCoordinates(lc[0].at(1), lc[1].at(1), lc[2].at(1),
357  lc[0].at(2), lc[1].at(2), lc[2].at(2),
358  lc[0].at(3), lc[1].at(3), lc[2].at(3));
359 
360  weight = gp->giveWeight();
362  return detJ * weight;
363 }
364 
365 
366 void
368 // Returns the lumped mass matrix of the receiver.
369 {
370  answer.resize(9, 9);
371  answer.zero();
372 
373  double mass = 0.0;
375  if ( iRule ) {
376  for ( GaussPoint *gp: *iRule ) {
377  mass += this->computeVolumeAround(gp) * this->giveStructuralCrossSection()->give('d', gp) * this->giveCrossSection()->give(CS_Thickness, gp);
378  }
379  }
380  double mass3 = mass/3.0;
381 
382  answer.at(1, 1) = mass3;
383  answer.at(4, 4) = mass3;
384  answer.at(7, 7) = mass3;
385 }
386 
387 
388 Interface *
390 {
391  if ( interface == LayeredCrossSectionInterfaceType ) {
392  return static_cast< LayeredCrossSectionInterface * >(this);
393  } else if ( interface == ZZNodalRecoveryModelInterfaceType ) {
394  return static_cast< ZZNodalRecoveryModelInterface * >(this);
395  } else if ( interface == NodalAveragingRecoveryModelInterfaceType ) {
396  return static_cast< NodalAveragingRecoveryModelInterface * >(this);
397  } else if ( interface == SPRNodalRecoveryModelInterfaceType ) {
398  return static_cast< SPRNodalRecoveryModelInterface * >(this);
399  } else if ( interface == ZZErrorEstimatorInterfaceType ) {
400  return static_cast< ZZErrorEstimatorInterface * >(this);
401  }
402 
403 
404  return NULL;
405 }
406 
407 
408 #define POINT_TOL 1.e-3
409 
410 bool
412 //converts global coordinates to local planar area coordinates,
413 //does not return a coordinate in the thickness direction, but
414 //does check that the point is in the element thickness
415 {
416  // get node coordinates
417  double x1, x2, x3, y1, y2, y3, z1, z2, z3;
418  this->giveNodeCoordinates(x1, x2, x3, y1, y2, y3, z1, z2, z3);
419 
420  // Fetch local coordinates.
421  bool ok = this->interp_lin.global2local( answer, coords, FEIElementGeometryWrapper(this) ) > 0;
422 
423  //get midplane location at this point
424  double midplZ;
425  midplZ = z1 * answer.at(1) + z2 * answer.at(2) + z3 * answer.at(3);
426 
427  //check that the z is within the element
428  GaussPoint _gp(NULL, 1, answer, 1.0, _2dPlate);
430  double elthick = cs->give(CS_Thickness, & _gp);
431 
432  if ( elthick / 2.0 + midplZ - fabs( coords.at(3) ) < -POINT_TOL ) {
433  answer.zero();
434  return false;
435  }
436 
437  //check that the point is in the element and set flag
438  for ( int i = 1; i <= 3; i++ ) {
439  if ( answer.at(i) < ( 0. - POINT_TOL ) ) {
440  return false;
441  }
442 
443  if ( answer.at(i) > ( 1. + POINT_TOL ) ) {
444  return false;
445  }
446  }
447 
448  return ok;
449 }
450 
451 
452 int
454 {
455  FloatArray help;
456  answer.resize(6);
457  if ( type == IST_ShellForceTensor || type == IST_ShellStrainTensor ) {
458  if ( type == IST_ShellForceTensor ) {
459  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
460  } else {
461  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
462  }
463  answer.at(1) = 0.0; // nx
464  answer.at(2) = 0.0; // ny
465  answer.at(3) = 0.0; // nz
466  answer.at(4) = help.at(5); // vyz
467  answer.at(5) = help.at(4); // vxz
468  answer.at(6) = 0.0; // vxy
469  return 1;
470  } else if ( type == IST_ShellMomentTensor || type == IST_CurvatureTensor ) {
471  if ( type == IST_ShellMomentTensor ) {
472  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
473  } else {
474  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
475  }
476  answer.at(1) = help.at(1); // mx
477  answer.at(2) = help.at(2); // my
478  answer.at(3) = 0.0; // mz
479  answer.at(4) = 0.0; // myz
480  answer.at(5) = 0.0; // mxz
481  answer.at(6) = help.at(3); // mxy
482  return 1;
483  } else {
484  return NLStructuralElement :: giveIPValue(answer, gp, type, tStep);
485  }
486 }
487 
488 void
490  InternalStateType type, TimeStep *tStep)
491 {
492  GaussPoint *gp;
493  if ( ( type == IST_ShellForceTensor ) || ( type == IST_ShellMomentTensor ) ||
494  ( type == IST_ShellStrainTensor ) || ( type == IST_CurvatureTensor ) ) {
495  gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
496  this->giveIPValue(answer, gp, type, tStep);
497  } else {
498  answer.clear();
499  }
500 }
501 
502 
503 void
505 {
506  pap.resize(3);
507  pap.at(1) = this->giveNode(1)->giveNumber();
508  pap.at(2) = this->giveNode(2)->giveNumber();
509  pap.at(3) = this->giveNode(3)->giveNumber();
510 }
511 
512 
513 void
515 {
516  answer.resize(1);
517  if ( ( pap == this->giveNode(1)->giveNumber() ) ||
518  ( pap == this->giveNode(2)->giveNumber() ) ||
519  ( pap == this->giveNode(3)->giveNumber() ) ) {
520  answer.at(1) = pap;
521  } else {
522  OOFEM_ERROR("node unknown");
523  }
524 }
525 
526 
529 {
530  return SPRPatchType_2dxy;
531 }
532 
533 
534 //
535 // layered cross section support functions
536 //
537 void
539  GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep)
540 // returns full 3d strain vector of given layer (whose z-coordinate from center-line is
541 // stored in slaveGp) for given tStep
542 {
543  double layerZeta, layerZCoord, top, bottom;
544 
545  top = this->giveCrossSection()->give(CS_TopZCoord, masterGp);
546  bottom = this->giveCrossSection()->give(CS_BottomZCoord, masterGp);
547  layerZeta = slaveGp->giveNaturalCoordinate(3);
548  layerZCoord = 0.5 * ( ( 1. - layerZeta ) * bottom + ( 1. + layerZeta ) * top );
549 
550  answer.resize(5); // {Exx,Eyy,GMyz,GMzx,GMxy}
551 
552  answer.at(1) = masterGpStrain.at(1) * layerZCoord;
553  answer.at(2) = masterGpStrain.at(2) * layerZCoord;
554  answer.at(5) = masterGpStrain.at(3) * layerZCoord;
555  answer.at(3) = masterGpStrain.at(5);
556  answer.at(4) = masterGpStrain.at(4);
557 }
558 
559 // Edge load support
560 /*
561 void
562 CCTPlate :: computeEgdeNMatrixAt(FloatMatrix &answer, int iedge, GaussPoint *gp)
563 {
564  IntArray edgeNodes;
565  FloatArray n;
566  double b, c, n12;
567 
568  this->interp_lin.edgeEvalN( n, iedge, gp->giveNaturalCoordinates(), FEIElementGeometryWrapper(this) );
569  this->interp_lin.computeLocalEdgeMapping(edgeNodes, iedge);
570 
571  n12 = 0.5 * n.at(1) * n.at(2);
572  b = this->giveNode( edgeNodes.at(1) )->giveCoordinate(2) - this->giveNode( edgeNodes.at(2) )->giveCoordinate(2);
573  c = this->giveNode( edgeNodes.at(2) )->giveCoordinate(1) - this->giveNode( edgeNodes.at(1) )->giveCoordinate(1);
574 
575 
576  answer.resize(3, 6);
577  answer.at(1, 1) = n.at(1);
578  answer.at(1, 2) = n12 * b;
579  answer.at(1, 3) = n12 * c;
580  answer.at(1, 4) = n.at(2);
581  answer.at(1, 5) = -n12 * b;
582  answer.at(1, 6) = -n12 * c;
583  //
584  answer.at(2, 2) = answer.at(3, 3) = n.at(1);
585  answer.at(2, 5) = answer.at(3, 6) = n.at(2);
586 }
587 */
588 
589 void
590 CCTPlate :: giveEdgeDofMapping(IntArray &answer, int iEdge) const
591 {
592  /*
593  * provides dof mapping of local edge dofs (only nonzero are taken into account)
594  * to global element dofs
595  */
596 
597  answer.resize(6);
598  if ( iEdge == 1 ) { // edge between nodes 1,2
599  answer.at(1) = 1;
600  answer.at(2) = 2;
601  answer.at(3) = 3;
602  answer.at(4) = 4;
603  answer.at(5) = 5;
604  answer.at(6) = 6;
605  } else if ( iEdge == 2 ) { // edge between nodes 2 3
606  answer.at(1) = 4;
607  answer.at(2) = 5;
608  answer.at(3) = 6;
609  answer.at(4) = 7;
610  answer.at(5) = 8;
611  answer.at(6) = 9;
612  } else if ( iEdge == 3 ) { // edge between nodes 3 1
613  answer.at(1) = 7;
614  answer.at(2) = 8;
615  answer.at(3) = 9;
616  answer.at(4) = 1;
617  answer.at(5) = 2;
618  answer.at(6) = 3;
619  } else {
620  OOFEM_ERROR("wrong edge number");
621  }
622 }
623 
624 double
626 {
628  return detJ *gp->giveWeight();
629 }
630 
631 
632 int
634 {
635  // returns transformation matrix from
636  // edge local coordinate system
637  // to element local c.s
638  // (same as global c.s in this case)
639  //
640  // i.e. f(element local) = T * f(edge local)
641  //
642  double dx, dy, length;
643  IntArray edgeNodes;
644  Node *nodeA, *nodeB;
645 
646  answer.resize(3, 3);
647  answer.zero();
648 
649  this->interp_lin.computeLocalEdgeMapping(edgeNodes, iEdge);
650 
651  nodeA = this->giveNode( edgeNodes.at(1) );
652  nodeB = this->giveNode( edgeNodes.at(2) );
653 
654  dx = nodeB->giveCoordinate(1) - nodeA->giveCoordinate(1);
655  dy = nodeB->giveCoordinate(2) - nodeA->giveCoordinate(2);
656  length = sqrt(dx * dx + dy * dy);
657 
658  answer.at(1, 1) = 1.0;
659  answer.at(2, 2) = dx / length;
660  answer.at(2, 3) = -dy / length;
661  answer.at(3, 2) = dy / length;
662  answer.at(3, 3) = dx / length;
663 
664  return 1;
665 }
666 
667 
668 void
669 CCTPlate::computeSurfaceNMatrix(FloatMatrix &answer, int boundaryID, const FloatArray &lcoords)
670 {
671  if (boundaryID == 1) {
672  this->computeNmatrixAt(lcoords, answer);
673  } else {
674  OOFEM_ERROR("computeSurfaceNMatrix: Unsupported surface id");
675  }
676 }
677 
678 
679 
680 //
681 // io routines
682 //
683 #ifdef __OOFEG
684 void
686 {
687  WCRec p [ 3 ];
688  GraphicObj *go;
689 
690  if ( !gc.testElementGraphicActivity(this) ) {
691  return;
692  }
693 
694  if ( this->giveMaterial()->isActivated(tStep) ) {
695  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
696  EASValsSetColor( gc.getElementColor() );
697  EASValsSetEdgeColor( gc.getElementEdgeColor() );
698  EASValsSetEdgeFlag(true);
699  EASValsSetFillStyle(FILL_SOLID);
700  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
701  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
702  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveCoordinate(2);
703  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveCoordinate(3);
704  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
705  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveCoordinate(2);
706  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveCoordinate(3);
707  p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveCoordinate(1);
708  p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveCoordinate(2);
709  p [ 2 ].z = ( FPNum ) this->giveNode(3)->giveCoordinate(3);
710 
711  go = CreateTriangle3D(p);
712  EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
713  EGAttachObject(go, ( EObjectP ) this);
714  EMAddGraphicsToModel(ESIModel(), go);
715  }
716 }
717 
718 
719 void
721 {
722  WCRec p [ 3 ];
723  GraphicObj *go;
724  double defScale = gc.getDefScale();
725 
726  if ( !gc.testElementGraphicActivity(this) ) {
727  return;
728  }
729 
730  if ( this->giveMaterial()->isActivated(tStep) ) {
731  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
732  EASValsSetColor( gc.getDeformedElementColor() );
733  EASValsSetEdgeColor( gc.getElementEdgeColor() );
734  EASValsSetEdgeFlag(true);
735  EASValsSetFillStyle(FILL_SOLID);
736  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
737  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
738  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale);
739  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(3, tStep, defScale);
740  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
741  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale);
742  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(3, tStep, defScale);
743  p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(1, tStep, defScale);
744  p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(2, tStep, defScale);
745  p [ 2 ].z = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(3, tStep, defScale);
746 
747  go = CreateTriangle3D(p);
748  EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
749  EMAddGraphicsToModel(ESIModel(), go);
750  }
751 }
752 
753 
754 void
756 {
757  int i, indx, result = 0;
758  WCRec p [ 3 ];
759  GraphicObj *tr;
760  FloatArray v1, v2, v3;
761  double s [ 3 ], defScale;
762 
763  if ( !gc.testElementGraphicActivity(this) ) {
764  return;
765  }
766 
767  if ( !this->giveMaterial()->isActivated(tStep) ) {
768  return;
769  }
770 
771  if ( gc.giveIntVarMode() == ISM_recovered ) {
772  result += this->giveInternalStateAtNode(v1, gc.giveIntVarType(), gc.giveIntVarMode(), 1, tStep);
773  result += this->giveInternalStateAtNode(v2, gc.giveIntVarType(), gc.giveIntVarMode(), 2, tStep);
774  result += this->giveInternalStateAtNode(v3, gc.giveIntVarType(), gc.giveIntVarMode(), 3, tStep);
775  } else if ( gc.giveIntVarMode() == ISM_local ) {
776  GaussPoint *gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
777  result += giveIPValue(v1, gp, gc.giveIntVarType(), tStep);
778  v2 = v1;
779  v3 = v1;
780  result *= 3;
781  }
782 
783  if ( result != 3 ) {
784  return;
785  }
786 
787  indx = gc.giveIntVarIndx();
788 
789  s [ 0 ] = v1.at(indx);
790  s [ 1 ] = v2.at(indx);
791  s [ 2 ] = v3.at(indx);
792 
793  EASValsSetLayer(OOFEG_VARPLOT_PATTERN_LAYER);
794 
795  if ( gc.getScalarAlgo() == SA_ISO_SURF ) {
796  for ( i = 0; i < 3; i++ ) {
797  if ( gc.getInternalVarsDefGeoFlag() ) {
798  // use deformed geometry
799  defScale = gc.getDefScale();
800  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
801  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
802  p [ i ].z = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(3, tStep, defScale);
803  } else {
804  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
805  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
806  p [ i ].z = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(3);
807  }
808  }
809 
810  //EASValsSetColor(gc.getYieldPlotColor(ratio));
811  gc.updateFringeTableMinMax(s, 3);
812  tr = CreateTriangleWD3D(p, s [ 0 ], s [ 1 ], s [ 2 ]);
813  EGWithMaskChangeAttributes(LAYER_MASK, tr);
814  EMAddGraphicsToModel(ESIModel(), tr);
815  }
816 }
817 
818 #endif
819 } // 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.
virtual void giveNodeCoordinates(double &x1, double &x2, double &x3, double &y1, double &y2, double &y3, double &z1, double &z2, double &z3)
Definition: cct.C:275
The element interface required by NodalAvergagingRecoveryModel.
virtual SPRPatchType SPRNodalRecoveryMI_givePatchType()
Definition: cct.C:528
virtual bool computeLocalCoordinates(FloatArray &answer, const FloatArray &gcoords)
Computes the element local coordinates from given global coordinates.
Definition: cct.C:411
The element interface required by ZZNodalRecoveryModel.
virtual FEInterpolation * giveInterpolation() const
Definition: cct.C:73
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: cct.C:268
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
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
ScalarAlgorithmType getScalarAlgo()
Bottom z coordinate.
Definition: crosssection.h:72
Wrapper around cell with vertex coordinates stored in FloatArray**.
Definition: feinterpol.h:115
virtual void computeLocalEdgeMapping(IntArray &edgeNodes, int iedge)
Definition: fei2dtrlin.C:230
The element interface required by ZZNodalRecoveryModel.
Abstract base class for "structural" finite elements with geometrical nonlinearities.
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...
virtual Interface * giveInterface(InterfaceType it)
Interface requesting service.
Definition: cct.C:389
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
#define OOFEG_RAW_GEOMETRY_LAYER
This class implements a structural material status information.
Definition: structuralms.h:65
void clear()
Clears receiver (zero size).
Definition: floatarray.h:206
virtual void SPRNodalRecoveryMI_giveSPRAssemblyPoints(IntArray &pap)
Definition: cct.C:504
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual double edgeGiveTransformationJacobian(int iedge, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the edge Jacobian of transformation between local and global coordinates.
Definition: feinterpol2d.C:175
Top z coordinate.
Definition: crosssection.h:71
virtual void computeComponentArrayAt(FloatArray &answer, TimeStep *tStep, ValueModeType mode)
Computes boundary condition value - its components values at given time.
Definition: load.C:82
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: cct.C:453
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
MatResponseMode
Describes the character of characteristic material matrix.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: cct.C:313
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
Definition: cct.C:755
#define OOFEG_DEFORMED_GEOMETRY_LAYER
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 give2dPlateStiffMtrx(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)=0
Method for computing 2d plate stiffness matrix.
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
virtual int computeLoadLEToLRotationMatrix(FloatMatrix &answer, int iEdge, GaussPoint *gp)
Returns transformation matrix from local edge c.s to element local coordinate system of load vector c...
Definition: cct.C:633
InternalStateType giveIntVarType()
double giveNaturalCoordinate(int i) const
Returns i-th natural element coordinate of receiver.
Definition: gausspoint.h:136
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType type)
Definition: cct.C:720
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: cct.C:213
The element interface corresponding to ZZErrorEstimator.
virtual double computeArea()
Computes the area (zero for all but 2d geometries).
Definition: cct.C:299
#define OOFEM_ERROR(...)
Definition: error.h:61
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: cct.C:321
REGISTER_Element(LSpace)
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
#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
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
Definition: cct.C:367
virtual void giveGeneralizedStress_Plate(FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)=0
virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge)
Computes volume related to integration point on local edge.
Definition: cct.C:625
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: cct.C:145
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
virtual int global2local(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Default implementation using Newton&#39;s method to find the local coordinates.
Definition: fei2dtrlin.C:103
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 computeBodyLoadVectorAt(FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode)
Computes the load vector due to body load acting on receiver, at given time step. ...
Definition: cct.C:100
virtual void computeStrainVectorInLayer(FloatArray &answer, const FloatArray &masterGpStrain, GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep)
Computes full 3D strain vector in element layer.
Definition: cct.C:538
virtual int setupIntegrationPoints(IntegrationRule &irule, int npoints, Element *element)
Sets up integration rule for the given element.
Definition: crosssection.C:54
int numberOfGaussPoints
Number of integration points as specified by nip.
Definition: element.h:188
InternalStateMode giveIntVarMode()
#define POINT_TOL
Definition: cct.C:408
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 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
double giveCharacteristicLengthForPlaneElements(const FloatArray &normalToCrackPlane)
Returns the size of element in the given direction if the direction is in the XY plane, otherwise gives the mean size defined as the square root of the element area.
Definition: element.C:1170
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
Definition: gausspoint.h:205
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 drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: cct.C:685
Class Interface.
Definition: interface.h:82
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
virtual bcGeomType giveBCGeoType() const
Returns geometry character of boundary condition.
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: cct.C:88
virtual void SPRNodalRecoveryMI_giveDofMansDeterminedByPatch(IntArray &answer, int pap)
Definition: cct.C:514
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
CCTPlate(int n, Domain *d)
Definition: cct.C:61
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
virtual bool isActivated(TimeStep *tStep)
Definition: material.h:161
virtual FloatArray * giveCoordinates()
Definition: node.h:114
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
double area
Definition: cct.h:67
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
void updateFringeTableMinMax(double *s, int size)
Load is base abstract class for all loads.
Definition: load.h:61
The element interface required by LayeredCrossSection.
virtual void computeMidPlaneNormal(FloatArray &answer, const GaussPoint *gp)
Computes mid-plane normal of receiver at integration point.
Definition: cct.C:328
virtual bool computeGtoLRotationMatrix(FloatMatrix &answer)
Returns transformation matrix from global c.s.
Definition: element.C:1537
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
Abstract base class for all structural cross section models.
virtual void NodalAveragingRecoveryMI_computeNodalValue(FloatArray &answer, int node, InternalStateType type, TimeStep *tStep)
Computes the element value in given node.
Definition: cct.C:489
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: cct.C:350
virtual bcValType giveBCValType() const
Returns receiver load type.
Class implementing node in finite element mesh.
Definition: node.h:87
int giveNumber() const
Definition: femcmpnn.h:107
virtual double giveCharacteristicLength(const FloatArray &normalToCrackPlane)
Returns the size of element in the given direction, in some cases adjusted (e.g.
Definition: cct.C:341
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
static FEI2dTrLin interp_lin
Definition: cct.h:64
virtual void giveEdgeDofMapping(IntArray &answer, int iEdge) const
Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element ...
Definition: cct.C:590
#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
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: cct.C:261
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.
virtual Material * giveMaterial()
Definition: element.C:484
Class representing Gaussian-quadrature integration rule.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631
virtual void computeSurfaceNMatrix(FloatMatrix &answer, int boundaryID, const FloatArray &lcoords)
Computes surface interpolation matrix.
Definition: cct.C:669

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