OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
qdkt.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/qdkt.h"
36 #include "../sm/Materials/structuralms.h"
37 #include "../sm/CrossSections/structuralcrosssection.h"
38 #include "fei2dquadlin.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(QDKTPlate);
57 
58 FEI2dQuadLin QDKTPlate :: interp_lin(1, 2);
59 
60 QDKTPlate :: QDKTPlate(int n, Domain *aDomain) :
61  NLStructuralElement(n, aDomain),
65 {
66  numberOfDofMans = 4;
68 }
69 
70 
73 
74 
77 {
78  return & interp_lin;
79 }
80 
81 
82 void
84 // Sets up the array containing the four Gauss points of the receiver.
85 {
86  if ( integrationRulesArray.size() == 0 ) {
87  integrationRulesArray.resize( 1 );
88  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 5) );
90  }
91 }
92 
93 void
95 // Returns the [5x12] strain-displacement matrix {B} of the receiver,
96 // evaluated at gp.
97 {
98  // get node coordinates
99  double x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4;
100  this->giveNodeCoordinates(x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4);
101 
102  // get gp coordinates
103  double ksi = gp->giveNaturalCoordinate(1);
104  double eta = gp->giveNaturalCoordinate(2);
105 
106  // geometrical characteristics of element sides
107  double dx4 = x2-x1;
108  double dy4 = y2-y1;
109  double l4 = sqrt(dx4*dx4+dy4*dy4);
110 
111  double dx5 = x3-x2;
112  double dy5 = y3-y2;
113  double l5 = sqrt(dx5*dx5+dy5*dy5);
114 
115  double dx6 = x4-x3;
116  double dy6 = y4-y3;
117  double l6 = sqrt(dx6*dx6+dy6*dy6);
118 
119  double dx7 = x1-x4;
120  double dy7 = y1-y4;
121  double l7 = sqrt(dx7*dx7+dy7*dy7);
122 
123  double c4 = dy4/l4;
124  double s4 = -dx4/l4;
125 
126  double c5 = dy5/l5;
127  double s5 = -dx5/l5;
128 
129  double c6 = dy6/l6;
130  double s6 = -dx6/l6;
131 
132  double c7 = dy7/l7;
133  double s7 = -dx7/l7;
134 
135  // transformation matrix from vertex dofs (w,fi_x, fi_y) to quadratic rotation DOFs
136  double T101 = -3./2./l4*c4;
137  double T102 = -1./4.*c4*c4+1./2.*s4*s4;
138  double T103 = -1./4.*c4*s4-1./2.*c4*s4;
139  double T104 = 3./2./l4*c4;
140  double T105 = -1./4.*c4*c4+1./2.*s4*s4;
141  double T106 = -1./4.*c4*s4-1./2.*c4*s4;
142 
143  double T201 = -3./2./l4*s4;
144  double T202 = -1./4.*c4*s4-1./2.*c4*s4;
145  double T203 = -1./4.*s4*s4+1./2.*c4*c4;
146  double T204 = 3./2./l4*s4;
147  double T205 = -1./4.*c4*s4-1./2.*c4*s4;
148  double T206 = -1./4.*s4*s4+1./2.*c4*c4;
149 
150  double T304 = -3./2./l5*c5;
151  double T305 = -1./4.*c5*c5+1./2.*s5*s5;
152  double T306 = -1./4.*c5*s5-1./2.*c5*s5;
153  double T307 = 3./2./l5*c5;
154  double T308 = -1./4.*c5*c5+1./2.*s5*s5;
155  double T309 = -1./4.*c5*s5-1./2.*c5*s5;
156 
157  double T404 = -3./2./l5*s5;
158  double T405 = -1./4.*c5*s5-1./2.*c5*s5;
159  double T406 = -1./4.*s5*s5+1./2.*c5*c5;
160  double T407 = 3./2./l5*s5;
161  double T408 = -1./4.*c5*s5-1./2.*c5*s5;
162  double T409 = -1./4.*s5*s5+1./2.*c5*c5;
163 
164  double T507 = -3./2./l6*c6;
165  double T508 = -1./4.*c6*c6+1./2.*s6*s6;
166  double T509 = -1./4.*c6*s6-1./2.*c6*s6;
167  double T510 = 3./2./l6*c6;
168  double T511 = -1./4.*c6*c6+1./2.*s6*s6;
169  double T512 = -1./4.*c6*s6-1./2.*c6*s6;
170 
171  double T607 = -3./2./l6*s6;
172  double T608 = -1./4.*c6*s6-1./2.*c6*s6;
173  double T609 = -1./4.*s6*s6+1./2.*c6*c6;
174  double T610 = 3./2./l6*s6;
175  double T611 = -1./4.*c6*s6-1./2.*c6*s6;
176  double T612 = -1./4.*s6*s6+1./2.*c6*c6;
177 
178  double T701 = 3./2./l7*c7;
179  double T702 = -1./4.*c7*c7+1./2.*s7*s7;
180  double T703 = -1./4.*c7*s7-1./2.*c7*s7;
181  double T710 = -3./2./l7*c7;
182  double T711 = -1./4.*c7*c7+1./2.*s7*s7;
183  double T712 = -1./4.*c7*s7-1./2.*c7*s7;
184 
185  double T801 = 3./2./l7*s7;
186  double T802 = -1./4.*c7*s7-1./2.*c7*s7;
187  double T803 = -1./4.*s7*s7+1./2.*c7*c7;
188  double T810 = -3./2./l7*s7;
189  double T811 = -1./4.*c7*s7-1./2.*c7*s7;
190  double T812 = -1./4.*s7*s7+1./2.*c7*c7;
191 
192  // derivatives of quadratic interpolation functions
193  // we do not have "midside" nodes -> explicit here
194  double N1dk = 0.25*(2.0*ksi+eta)*(1.0+eta);
195  double N2dk = 0.25*(2.0*ksi-eta)*(1.0+eta);
196  double N3dk = 0.25*(2.0*ksi+eta)*(1.0-eta);
197  double N4dk = 0.25*(2.0*ksi-eta)*(1.0-eta);
198  double N7dk = -ksi*(1.0-eta);
199  double N8dk = 0.5*(1.0-eta*eta);
200  double N5dk = -ksi*(1.0+eta);
201  double N6dk = -0.5*(1.0-eta*eta);
202 
203  double N3de = 0.25*(2.0*eta+ksi)*(1.0-ksi);
204  double N4de = 0.25*(2.0*eta-ksi)*(1.0+ksi);
205  double N1de = 0.25*(2.0*eta+ksi)*(1.0+ksi);
206  double N2de = 0.25*(2.0*eta-ksi)*(1.0-ksi);
207  double N7de = -0.5*(1.0-ksi*ksi);
208  double N8de = -eta*(1.0+ksi);
209  double N5de = 0.5*(1.0-ksi*ksi);
210  double N6de = -eta*(1.0-ksi);
211 
212  double detJ = 1./8.*((y4-y2)*(x3-x1)-(y3-y1)*(x4-x2))+
213  ksi/8*((y3-y4)*(x2-x1)-(y2-y1)*(x3-x4))+
214  eta/8*((y4-y1)*(x3-x2)-(y3-y2)*(x4-x1));
215 
216  double dxdk = -1.0/detJ * ((y3-y2)+(y4-y1+ksi*(y1-y2+y3-y4)))/4.0;
217  double dxde = 1.0/detJ * ((y2-y1)+(y3-y4+eta*(y1-y2+y3-y4)))/4.0;
218  double dydk = 1.0/detJ * ((x3-x2)+(x4-x1+ksi*(x1-x2+x3-x4)))/4.0;
219  double dyde = -1.0/detJ * ((x2-x1)+(x3-x4+eta*(x1-x2+x3-x4)))/4.0;
220 
221  double dN102 = N1dk*dxdk+N1de*dxde;
222  double dN104 = N2dk*dxdk+N2de*dxde;
223  double dN106 = N3dk*dxdk+N3de*dxde;
224  double dN108 = N4dk*dxdk+N4de*dxde;
225  double dN110 = N5dk*dxdk+N5de*dxde;
226  double dN112 = N6dk*dxdk+N6de*dxde;
227  double dN114 = N7dk*dxdk+N7de*dxde;
228  double dN116 = N8dk*dxdk+N8de*dxde;
229 
230  double dN201 = -N1dk*dydk-N1de*dyde;
231  double dN203 = -N2dk*dydk-N2de*dyde;
232  double dN205 = -N3dk*dydk-N3de*dyde;
233  double dN207 = -N4dk*dydk-N4de*dyde;
234  double dN209 = -N5dk*dydk-N5de*dyde;
235  double dN211 = -N6dk*dydk-N6de*dyde;
236  double dN213 = -N7dk*dydk-N7de*dyde;
237  double dN215 = -N8dk*dydk-N8de*dyde;
238 
239  answer.resize(5, 12);
240  answer.zero();
241 
242  answer.at(1,1) = T201*dN110 + T801*dN116;
243  answer.at(1,2) = T202*dN110 + T802*dN116;
244  answer.at(1,3) = dN102 + T203*dN110 + T803*dN116;
245  answer.at(1,4) = T204*dN110 + T404*dN112;
246  answer.at(1,5) = T205*dN110 + T405*dN112;
247  answer.at(1,6) = dN104 + T206*dN110 + T406*dN112;
248  answer.at(1,7) = T407*dN112 + T607*dN114;
249  answer.at(1,8) = T408*dN112 + T608*dN114;
250  answer.at(1,9) = dN106 + T409*dN112 + T609*dN114;
251  answer.at(1,10)= T610*dN114 + T810*dN116;
252  answer.at(1,11)= T611*dN114 + T811*dN116;
253  answer.at(1,12)= dN108 + T612*dN114 + T812*dN116;
254 
255  answer.at(2,1) = T101*dN209 + T701*dN215;
256  answer.at(2,2) = dN201 + T102*dN209 + T702*dN215;
257  answer.at(2,3) = T103*dN209 + T703*dN215;
258  answer.at(2,4) = T104*dN209 + T304*dN211;
259  answer.at(2,5) = dN203 + T105*dN209 + T305*dN211;
260  answer.at(2,6) = T106*dN209 + T306*dN211;
261  answer.at(2,7) = T307*dN211 + T507*dN213;
262  answer.at(2,8) = dN205 + T308*dN211 + T508*dN213;
263  answer.at(2,9) = T309*dN211 + T509*dN213;
264  answer.at(2,10)= T510*dN213 + T710*dN215;
265  answer.at(2,11)= dN207 + T511*dN213 + T711*dN215;
266  answer.at(2,12)= T512*dN213 + T712*dN215;
267 
268  answer.at(3,1) = - T101*dN110 - T201*dN209 - T701*dN116 - T801*dN215;
269  answer.at(3,2) = - dN102 - T102*dN110 - T202*dN209 - T702*dN116 - T802*dN215;
270  answer.at(3,3) = - dN201 - T103*dN110 - T203*dN209 - T703*dN116 - T803*dN215;
271  answer.at(3,4) = - T104*dN110 - T204*dN209 - T304*dN112 - T404*dN211;
272  answer.at(3,5) = - dN104 - T105*dN110 - T205*dN209 - T305*dN112 - T405*dN211;
273  answer.at(3,6) = - dN203 - T106*dN110 - T206*dN209 - T306*dN112 - T406*dN211;
274  answer.at(3,7) = - T307*dN112 - T407*dN211 - T507*dN114 - T607*dN213;
275  answer.at(3,8) = - dN106 - T308*dN112 - T408*dN211 - T508*dN114 - T608*dN213;
276  answer.at(3,9) = - dN205 - T309*dN112 - T409*dN211 - T509*dN114 - T609*dN213;
277  answer.at(3,10)= - T510*dN114 - T610*dN213 - T710*dN116 - T810*dN215;
278  answer.at(3,11)= - dN108 - T511*dN114 - T611*dN213 - T711*dN116 - T811*dN215;
279  answer.at(3,12)= - dN207 - T512*dN114 - T612*dN213 - T712*dN116 - T812*dN215;
280 
281  // Note: no shear strains, no shear forces => the 4th and 5th rows are zero
282 }
283 
284 
285 void
287 // Returns the [3x9] displacement interpolation matrix {N} of the receiver,
288 // evaluated at gp.
289 // Note: this interpolation is not available, as the deflection is cubic along the edges,
290 // but not define in the interior of the element
291 // Note: the interpolation of rotations is quadratic
292 // NOTE: linear interpolation returned instead
293 {
294  FloatArray N;
295 
296  giveInterpolation()->evalN( N, iLocCoord, FEIElementGeometryWrapper(this) );
297 
298  answer.beNMatrixOf(N, 3);
299 
300 }
301 
302 
303 void
305 {
306  this->giveStructuralCrossSection()->giveGeneralizedStress_Plate(answer, gp, strain, tStep);
307 }
308 
309 
310 void
312 {
313  this->giveStructuralCrossSection()->give2dPlateStiffMtrx(answer, rMode, gp, tStep);
314 }
315 
316 
317 void
318 QDKTPlate :: giveNodeCoordinates(double &x1, double &x2, double &x3, double &x4,
319  double &y1, double &y2, double &y3, double &y4,
320  double &z1, double &z2, double &z3, double &z4)
321 {
322  FloatArray *nc1, *nc2, *nc3, *nc4;
323  nc1 = this->giveNode(1)->giveCoordinates();
324  nc2 = this->giveNode(2)->giveCoordinates();
325  nc3 = this->giveNode(3)->giveCoordinates();
326  nc4 = this->giveNode(4)->giveCoordinates();
327 
328  x1 = nc1->at(1);
329  x2 = nc2->at(1);
330  x3 = nc3->at(1);
331  x4 = nc4->at(1);
332 
333  y1 = nc1->at(2);
334  y2 = nc2->at(2);
335  y3 = nc3->at(2);
336  y4 = nc4->at(2);
337 
338  z1 = nc1->at(3);
339  z2 = nc2->at(3);
340  z3 = nc3->at(3);
341  z4 = nc4->at(3);
342 
343 }
344 
345 
348 {
350 }
351 
352 
353 void
355 {
356  answer = {D_w, R_u, R_v};
357 }
358 
359 
360 void
362 // returns normal vector to midPlane in GaussPoinr gp of receiver
363 {
364  FloatArray u, v;
365  u.beDifferenceOf( * this->giveNode(2)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
366  v.beDifferenceOf( * this->giveNode(3)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
367 
368  answer.beVectorProductOf(u, v);
369  answer.normalize();
370 }
371 
372 
373 double
375 //
376 // returns receiver's characteristic length for crack band models
377 // for a crack formed in the plane with normal normalToCrackPlane.
378 //
379 {
380  return this->giveCharacteristicLengthForPlaneElements(normalToCrackPlane);
381 }
382 
383 
384 double
386 // Returns the portion of the receiver which is attached to gp.
387 {
388  double detJ, weight;
389 
390  weight = gp->giveWeight();
392  return detJ * weight;
393 }
394 
395 
396 Interface *
398 {
399  if ( interface == LayeredCrossSectionInterfaceType ) {
400  return static_cast< LayeredCrossSectionInterface * >(this);
401  } else if ( interface == ZZNodalRecoveryModelInterfaceType ) {
402  return static_cast< ZZNodalRecoveryModelInterface * >(this);
403  } else if ( interface == SPRNodalRecoveryModelInterfaceType ) {
404  return static_cast< SPRNodalRecoveryModelInterface * >(this);
405  } else if ( interface == ZZErrorEstimatorInterfaceType ) {
406  return static_cast< ZZErrorEstimatorInterface * >(this);
407  }
408 
409 
410  return NULL;
411 }
412 
413 void
415 // Computes numerically the load vector of the receiver due to the body
416 // loads, at tStep.
417 // load is assumed to be in global cs.
418 // load vector is then transformed to coordinate system in each node.
419 // (should be global coordinate system, but there may be defined
420 // different coordinate system in each node)
421 {
422  double dens, dV;
423  FloatArray force, ntf;
424  FloatMatrix n, T;
425 
426  if ( ( forLoad->giveBCGeoType() != BodyLoadBGT ) || ( forLoad->giveBCValType() != ForceLoadBVT ) ) {
427  OOFEM_ERROR("unknown load type");
428  }
429 
430  // note: force is assumed to be in global coordinate system.
431  forLoad->computeComponentArrayAt(force, tStep, mode);
432  // transform from global to element local c.s
433  if ( this->computeLoadGToLRotationMtrx(T) ) {
434  force.rotatedWith(T, 'n');
435  }
436 
437  answer.clear();
438 
439  if ( force.giveSize() ) {
440  for ( GaussPoint *gp: *this->giveDefaultIntegrationRulePtr() ) {
441  this->computeNmatrixAt(gp->giveSubPatchCoordinates(), n);
442  dV = this->computeVolumeAround(gp) * this->giveCrossSection()->give(CS_Thickness, gp);
443  dens = this->giveCrossSection()->give('d', gp);
444  ntf.beTProductOf(n, force);
445  answer.add(dV * dens, ntf);
446  }
447  } else {
448  return;
449  }
450 }
451 
452 
453 #define POINT_TOL 1.e-3
454 
455 bool
457 //converts global coordinates to local planar area coordinates,
458 //does not return a coordinate in the thickness direction, but
459 //does check that the point is in the element thickness
460 {
461  // get node coordinates
462  double x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4;
463  this->giveNodeCoordinates(x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4);
464 
465  // Fetch local coordinates.
466  bool ok = this->interp_lin.global2local( answer, coords, FEIElementGeometryWrapper(this) ) > 0;
467 
468  //check that the point is in the element and set flag
469  for ( int i = 1; i <= 4; i++ ) {
470  if ( answer.at(i) < ( 0. - POINT_TOL ) ) {
471  return false;
472  }
473 
474  if ( answer.at(i) > ( 1. + POINT_TOL ) ) {
475  return false;
476  }
477  }
478 
479  //get midplane location at this point
480  double midplZ;
481  midplZ = z1 * answer.at(1) + z2 * answer.at(2) + z3 * answer.at(3) + z4 * answer.at(4);
482 
483  //check that the z is within the element
485  double elthick = cs->give(CS_Thickness, answer, this);
486 
487  if ( elthick / 2.0 + midplZ - fabs( coords.at(3) ) < -POINT_TOL ) {
488  answer.zero();
489  return false;
490  }
491 
492 
493  return ok;
494 }
495 
496 
497 int
499 {
500  FloatArray help;
501  answer.resize(6);
502  if ( type == IST_ShellForceTensor || type == IST_ShellStrainTensor ) {
503  if ( type == IST_ShellForceTensor ) {
504  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
505  } else {
506  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
507  }
508  answer.at(1) = 0.0; // nx
509  answer.at(2) = 0.0; // ny
510  answer.at(3) = 0.0; // nz
511  answer.at(4) = help.at(5); // vyz
512  answer.at(5) = help.at(4); // vxz
513  answer.at(6) = 0.0; // vxy
514  return 1;
515  } else if ( type == IST_ShellMomentTensor || type == IST_CurvatureTensor ) {
516  if ( type == IST_ShellMomentTensor ) {
517  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
518  } else {
519  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
520  }
521  answer.at(1) = help.at(1); // mx
522  answer.at(2) = help.at(2); // my
523  answer.at(3) = 0.0; // mz
524  answer.at(4) = 0.0; // myz
525  answer.at(5) = 0.0; // mxz
526  answer.at(6) = help.at(3); // mxy
527  return 1;
528  } else {
529  return NLStructuralElement :: giveIPValue(answer, gp, type, tStep);
530  }
531 }
532 
533 
534 void
536 {
537  pap.resize(4);
538  pap.at(1) = this->giveNode(1)->giveNumber();
539  pap.at(2) = this->giveNode(2)->giveNumber();
540  pap.at(3) = this->giveNode(3)->giveNumber();
541  pap.at(4) = this->giveNode(4)->giveNumber();
542 }
543 
544 
545 void
547 {
548  answer.resize(1);
549  if ( ( pap == this->giveNode(1)->giveNumber() ) ||
550  ( pap == this->giveNode(2)->giveNumber() ) ||
551  ( pap == this->giveNode(3)->giveNumber() ) ||
552  ( pap == this->giveNode(4)->giveNumber() ) ) {
553  answer.at(1) = pap;
554  } else {
555  OOFEM_ERROR("node unknown");
556  }
557 }
558 
559 
562 {
563  return SPRPatchType_2dxy;
564 }
565 
566 
567 //
568 // layered cross section support functions
569 //
570 void
572  GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep)
573 // returns full 3d strain vector of given layer (whose z-coordinate from center-line is
574 // stored in slaveGp) for given tStep
575 {
576  double layerZeta, layerZCoord, top, bottom;
577 
578  top = this->giveCrossSection()->give(CS_TopZCoord, masterGp);
579  bottom = this->giveCrossSection()->give(CS_BottomZCoord, masterGp);
580  layerZeta = slaveGp->giveNaturalCoordinate(3);
581  layerZCoord = 0.5 * ( ( 1. - layerZeta ) * bottom + ( 1. + layerZeta ) * top );
582 
583  answer.resize(5); // {Exx,Eyy,GMyz,GMzx,GMxy}
584 
585  answer.at(1) = masterGpStrain.at(1) * layerZCoord;
586  answer.at(2) = masterGpStrain.at(2) * layerZCoord;
587  answer.at(5) = masterGpStrain.at(3) * layerZCoord;
588  answer.at(3) = masterGpStrain.at(5);
589  answer.at(4) = masterGpStrain.at(4);
590 }
591 
592 void
593 QDKTPlate :: giveEdgeDofMapping(IntArray &answer, int iEdge) const
594 {
595  if ( iEdge == 1 ) { // edge between nodes 1 2
596  answer = {1, 2, 3, 4, 5, 6};
597  } else if ( iEdge == 2 ) { // edge between nodes 2 3
598  answer = {4, 5, 6, 7, 8, 9};
599  } else if ( iEdge == 3 ) { // edge between nodes 3 4
600  answer = {7, 8, 9, 10, 11, 12};
601  } else if ( iEdge == 4 ) { // edge between nodes 4 1
602  answer = {10, 11, 12, 1, 2, 3};
603  } else {
604  OOFEM_ERROR("wrong edge number");
605  }
606 }
607 
608 double
610 {
612  return detJ *gp->giveWeight();
613 }
614 
615 
616 int
618 {
619  // returns transformation matrix from
620  // edge local coordinate system
621  // to element local c.s
622  // (same as global c.s in this case)
623  //
624  // i.e. f(element local) = T * f(edge local)
625  //
626  double dx, dy, length;
627  IntArray edgeNodes;
628  Node *nodeA, *nodeB;
629 
630  answer.resize(3, 3);
631  answer.zero();
632 
633  this->interp_lin.computeLocalEdgeMapping(edgeNodes, iEdge);
634 
635  nodeA = this->giveNode( edgeNodes.at(1) );
636  nodeB = this->giveNode( edgeNodes.at(2) );
637 
638  dx = nodeB->giveCoordinate(1) - nodeA->giveCoordinate(1);
639  dy = nodeB->giveCoordinate(2) - nodeA->giveCoordinate(2);
640  length = sqrt(dx * dx + dy * dy);
641 
642  answer.at(1, 1) = 1.0;
643  answer.at(2, 2) = dx / length;
644  answer.at(2, 3) = -dy / length;
645  answer.at(3, 2) = dy / length;
646  answer.at(3, 3) = dx / length;
647 
648  return 1;
649 }
650 
651 
652 void
654 {
655  this->computeNmatrixAt(sgp->giveNaturalCoordinates(), answer);
656 }
657 
658 void
659 QDKTPlate::computeSurfaceNMatrix (FloatMatrix &answer, int boundaryID, const FloatArray& lcoords)
660 {
661  FloatArray n_vec;
662  this->giveInterpolation()->boundarySurfaceEvalN(n_vec, boundaryID, lcoords, FEIElementGeometryWrapper(this) );
663  answer.beNMatrixOf(n_vec, 3);
664 }
665 
666 void
668 {
669  answer.resize(12);
670  answer.zero();
671  if ( iSurf == 1 ) {
672  for (int i = 1; i<=12; i++) {
673  answer.at(i) = i;
674  }
675  } else {
676  OOFEM_ERROR("wrong surface number");
677  }
678 }
679 
682 {
683  IntegrationRule *iRule = new GaussIntegrationRule(1, this, 1, 1);
684  int npoints = iRule->getRequiredNumberOfIntegrationPoints(_Square, approxOrder);
685  iRule->SetUpPointsOnSquare(npoints, _Unknown);
686  return iRule;
687 }
688 
689 double
691 {
692  return this->computeVolumeAround(gp);
693 }
694 
695 
696 int
698 {
699  return 0;
700 }
701 
702 
703 //
704 // io routines
705 //
706 #ifdef __OOFEG
707 void
709 {
710  WCRec p [ 4 ];
711  GraphicObj *go;
712 
713  if ( !gc.testElementGraphicActivity(this) ) {
714  return;
715  }
716 
717  if ( this->giveMaterial()->isActivated(tStep) ) {
718  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
719  EASValsSetColor( gc.getElementColor() );
720  EASValsSetEdgeColor( gc.getElementEdgeColor() );
721  EASValsSetEdgeFlag(true);
722  EASValsSetFillStyle(FILL_SOLID);
723  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
724  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
725  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveCoordinate(2);
726  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveCoordinate(3);
727  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
728  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveCoordinate(2);
729  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveCoordinate(3);
730  p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveCoordinate(1);
731  p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveCoordinate(2);
732  p [ 2 ].z = ( FPNum ) this->giveNode(3)->giveCoordinate(3);
733  p [ 3 ].x = ( FPNum ) this->giveNode(4)->giveCoordinate(1);
734  p [ 3 ].y = ( FPNum ) this->giveNode(4)->giveCoordinate(2);
735  p [ 3 ].z = 0.;
736 
737  go = CreateQuad3D(p);
738  EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
739  EGAttachObject(go, ( EObjectP ) this);
740  EMAddGraphicsToModel(ESIModel(), go);
741  }
742 }
743 
744 
745 void
747 {
748  WCRec p [ 4 ];
749  GraphicObj *go;
750  double defScale = gc.getDefScale();
751 
752  if ( !gc.testElementGraphicActivity(this) ) {
753  return;
754  }
755 
756  if ( this->giveMaterial()->isActivated(tStep) ) {
757  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
758  EASValsSetColor( gc.getDeformedElementColor() );
759  EASValsSetEdgeColor( gc.getElementEdgeColor() );
760  EASValsSetEdgeFlag(true);
761  EASValsSetFillStyle(FILL_SOLID);
762  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
763  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
764  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale);
765  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(3, tStep, defScale);
766  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
767  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale);
768  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(3, tStep, defScale);
769  p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(1, tStep, defScale);
770  p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(2, tStep, defScale);
771  p [ 2 ].z = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(3, tStep, defScale);
772  p [ 3 ].x = ( FPNum ) this->giveNode(4)->giveUpdatedCoordinate(1, tStep, defScale);
773  p [ 3 ].y = ( FPNum ) this->giveNode(4)->giveUpdatedCoordinate(2, tStep, defScale);
774  p [ 3 ].z = ( FPNum ) this->giveNode(4)->giveUpdatedCoordinate(3, tStep, defScale);
775 
776  go = CreateQuad3D(p);
777  EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
778  EMAddGraphicsToModel(ESIModel(), go);
779  }
780 }
781 
782 
783 void
785 {
786  int i, indx, result = 0;
787  WCRec p [ 4 ];
788  GraphicObj *tr;
789  FloatArray v [ 4 ];
790  double s [ 4 ], defScale;
791 
792  if ( !gc.testElementGraphicActivity(this) ) {
793  return;
794  }
795 
796  EASValsSetLayer(OOFEG_VARPLOT_PATTERN_LAYER);
797  if ( gc.giveIntVarMode() == ISM_recovered ) {
798  for ( i = 1; i <= 4; i++ ) {
799  result += this->giveInternalStateAtNode(v [ i - 1 ], gc.giveIntVarType(), gc.giveIntVarMode(), i, tStep);
800  }
801 
802  if ( result != 4 ) {
803  return;
804  }
805 
806  indx = gc.giveIntVarIndx();
807 
808  for ( i = 1; i <= 4; i++ ) {
809  s [ i - 1 ] = v [ i - 1 ].at(indx);
810  }
811 
812  if ( gc.getScalarAlgo() == SA_ISO_SURF ) {
813  for ( i = 0; i < 4; i++ ) {
814  if ( gc.getInternalVarsDefGeoFlag() ) {
815  // use deformed geometry
816  defScale = gc.getDefScale();
817  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
818  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
819  p [ i ].z = 0.;
820  } else {
821  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
822  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
823  p [ i ].z = 0.;
824  }
825  }
826 
827  //EASValsSetColor(gc.getYieldPlotColor(ratio));
828  gc.updateFringeTableMinMax(s, 4);
829  tr = CreateQuadWD3D(p, s [ 0 ], s [ 1 ], s [ 2 ], s [ 3 ]);
830  EGWithMaskChangeAttributes(LAYER_MASK, tr);
831  EMAddGraphicsToModel(ESIModel(), tr);
832 
833  } else if ( ( gc.getScalarAlgo() == SA_ZPROFILE ) || ( gc.getScalarAlgo() == SA_COLORZPROFILE ) ) {
834  double landScale = gc.getLandScale();
835 
836  for ( i = 0; i < 4; i++ ) {
837  if ( gc.getInternalVarsDefGeoFlag() ) {
838  // use deformed geometry
839  defScale = gc.getDefScale();
840  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
841  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
842  p [ i ].z = s [ i ] * landScale;
843  } else {
844  p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
845  p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
846  p [ i ].z = s [ i ] * landScale;
847  }
848 
849  // this fixes a bug in ELIXIR
850  if ( fabs(s [ i ]) < 1.0e-6 ) {
851  s [ i ] = 1.0e-6;
852  }
853  }
854 
855  if ( gc.getScalarAlgo() == SA_ZPROFILE ) {
856  EASValsSetColor( gc.getDeformedElementColor() );
857  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
858  tr = CreateQuad3D(p);
859  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, tr);
860  } else {
861  gc.updateFringeTableMinMax(s, 4);
862  tr = CreateQuadWD3D(p, s [ 0 ], s [ 1 ], s [ 2 ], s [ 3 ]);
863  EGWithMaskChangeAttributes(LAYER_MASK, tr);
864  }
865 
866  EMAddGraphicsToModel(ESIModel(), tr);
867  }
868  } else if ( gc.giveIntVarMode() == ISM_local ) {
869  if ( integrationRulesArray [ 0 ]->giveNumberOfIntegrationPoints() != 4 ) {
870  return;
871  }
872 
873  IntArray ind(4);
874  WCRec pp [ 9 ];
875 
876  for ( i = 0; i < 4; i++ ) {
877  if ( gc.getInternalVarsDefGeoFlag() ) {
878  // use deformed geometry
879  defScale = gc.getDefScale();
880  pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
881  pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
882  pp [ i ].z = 0.;
883  } else {
884  pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
885  pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
886  pp [ i ].z = 0.;
887  }
888  }
889 
890  for ( i = 0; i < 3; i++ ) {
891  pp [ i + 4 ].x = 0.5 * ( pp [ i ].x + pp [ i + 1 ].x );
892  pp [ i + 4 ].y = 0.5 * ( pp [ i ].y + pp [ i + 1 ].y );
893  pp [ i + 4 ].z = 0.5 * ( pp [ i ].z + pp [ i + 1 ].z );
894  }
895 
896  pp [ 7 ].x = 0.5 * ( pp [ 3 ].x + pp [ 0 ].x );
897  pp [ 7 ].y = 0.5 * ( pp [ 3 ].y + pp [ 0 ].y );
898  pp [ 7 ].z = 0.5 * ( pp [ 3 ].z + pp [ 0 ].z );
899 
900  pp [ 8 ].x = 0.25 * ( pp [ 0 ].x + pp [ 1 ].x + pp [ 2 ].x + pp [ 3 ].x );
901  pp [ 8 ].y = 0.25 * ( pp [ 0 ].y + pp [ 1 ].y + pp [ 2 ].y + pp [ 3 ].y );
902  pp [ 8 ].z = 0.25 * ( pp [ 0 ].z + pp [ 1 ].z + pp [ 2 ].z + pp [ 3 ].z );
903 
904  for ( GaussPoint *gp: *this->giveDefaultIntegrationRulePtr() ) {
905  const FloatArray &gpCoords = gp->giveNaturalCoordinates();
906  if ( ( gpCoords.at(1) > 0. ) && ( gpCoords.at(2) > 0. ) ) {
907  ind.at(1) = 0;
908  ind.at(2) = 4;
909  ind.at(3) = 8;
910  ind.at(4) = 7;
911  } else if ( ( gpCoords.at(1) < 0. ) && ( gpCoords.at(2) > 0. ) ) {
912  ind.at(1) = 4;
913  ind.at(2) = 1;
914  ind.at(3) = 5;
915  ind.at(4) = 8;
916  } else if ( ( gpCoords.at(1) < 0. ) && ( gpCoords.at(2) < 0. ) ) {
917  ind.at(1) = 5;
918  ind.at(2) = 2;
919  ind.at(3) = 6;
920  ind.at(4) = 8;
921  } else {
922  ind.at(1) = 6;
923  ind.at(2) = 3;
924  ind.at(3) = 7;
925  ind.at(4) = 8;
926  }
927 
928  if ( giveIPValue(v [ 0 ], gp, gc.giveIntVarType(), tStep) == 0 ) {
929  return;
930  }
931 
932  indx = gc.giveIntVarIndx();
933 
934  for ( i = 1; i <= 4; i++ ) {
935  s [ i - 1 ] = v [ 0 ].at(indx);
936  }
937 
938  for ( i = 0; i < 4; i++ ) {
939  p [ i ].x = pp [ ind.at(i + 1) ].x;
940  p [ i ].y = pp [ ind.at(i + 1) ].y;
941  p [ i ].z = pp [ ind.at(i + 1) ].z;
942  }
943 
944  gc.updateFringeTableMinMax(s, 4);
945  tr = CreateQuadWD3D(p, s [ 0 ], s [ 1 ], s [ 2 ], s [ 3 ]);
946  EGWithMaskChangeAttributes(LAYER_MASK, tr);
947  EMAddGraphicsToModel(ESIModel(), tr);
948  }
949  }
950 }
951 
952 #endif
953 } // end namespace oofem
CrossSection * giveCrossSection()
Definition: element.C:495
int testElementGraphicActivity(Element *)
Test if particular element passed fulfills various filtering criteria for its graphics output...
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
The element interface required by ZZNodalRecoveryModel.
virtual int global2local(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Default implementation using Newton&#39;s method to find the local coordinates.
Definition: fei2dquadlin.C:120
virtual bool computeLocalCoordinates(FloatArray &answer, const FloatArray &gcoords)
Computes the element local coordinates from given global coordinates.
Definition: qdkt.C:456
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
virtual void evalN(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)=0
Evaluates the array of interpolation functions (shape functions) at given point.
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 void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: qdkt.C:83
virtual IntegrationRule * giveDefaultIntegrationRulePtr()
Access method for default integration rule.
Definition: element.h:822
ScalarAlgorithmType getScalarAlgo()
Bottom z coordinate.
Definition: crosssection.h:72
virtual double giveTransformationJacobian(const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation.
Definition: feinterpol.C:43
virtual FEInterpolation * giveInterpolation() const
Definition: qdkt.C:72
The element interface required by ZZNodalRecoveryModel.
Abstract base class for "structural" finite elements with geometrical nonlinearities.
virtual void SPRNodalRecoveryMI_giveDofMansDeterminedByPatch(IntArray &answer, int pap)
Definition: qdkt.C:546
void zero()
Sets all component to zero.
Definition: intarray.C:52
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...
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
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual void computeMidPlaneNormal(FloatArray &answer, const GaussPoint *gp)
Computes mid-plane normal of receiver at integration point.
Definition: qdkt.C:361
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: qdkt.C:286
static FEI2dQuadLin interp_lin
Element geometry approximation.
Definition: qdkt.h:71
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: qdkt.C:347
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
virtual void boundarySurfaceEvalN(FloatArray &answer, int isurf, const FloatArray &lcoords, const FEICellGeometry &cellgeo)=0
Evaluates the array of edge interpolation functions (shape functions) at given point.
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 SetUpPointsOnSquare(int, MaterialMode mode)
Sets up receiver&#39;s integration points on unit square integration domain.
virtual double giveCoordinate(int i)
Definition: node.C:82
virtual double giveCharacteristicLength(const FloatArray &normalToCrackPlane)
Returns the size of element in the given direction, in some cases adjusted (e.g.
Definition: qdkt.C:374
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.
#define OOFEG_DEFORMED_GEOMETRY_LAYER
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: qdkt.C:354
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: qdkt.C:708
void rotatedWith(FloatMatrix &r, char mode)
Returns the receiver a rotated according the change-of-base matrix r.
Definition: floatarray.C:799
#define POINT_TOL
Definition: qdkt.C:453
Abstract base class representing integration rule.
Distributed body load.
Definition: bcgeomtype.h:43
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: qdkt.C:617
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
Definition: qdkt.C:94
void beDifferenceOf(const FloatArray &a, const FloatArray &b)
Sets receiver to be a - b.
Definition: floatarray.C:341
QDKTPlate(int n, Domain *d)
Definition: qdkt.C:60
virtual SPRPatchType SPRNodalRecoveryMI_givePatchType()
Definition: qdkt.C:561
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
InternalStateType giveIntVarType()
double giveNaturalCoordinate(int i) const
Returns i-th natural element coordinate of receiver.
Definition: gausspoint.h:136
virtual double computeSurfaceVolumeAround(GaussPoint *gp, int iSurf)
Computes volume related to integration point on local surface.
Definition: qdkt.C:690
The element interface corresponding to ZZErrorEstimator.
virtual void computeLocalEdgeMapping(IntArray &edgeNodes, int iedge)
Definition: fei2dquadlin.C:295
#define OOFEM_ERROR(...)
Definition: error.h:61
REGISTER_Element(LSpace)
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
virtual int computeLoadGToLRotationMtrx(FloatMatrix &answer)
Returns transformation matrix from global coordinate system to local element coordinate system for el...
#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.
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
Definition: qdkt.C:784
virtual void giveGeneralizedStress_Plate(FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)=0
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
#define N(p, q)
Definition: mdm.C:367
void beTProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix^T and anArray.
Definition: floatarray.C:708
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 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: qdkt.C:304
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: qdkt.C:311
virtual void giveEdgeDofMapping(IntArray &answer, int iEdge) const
Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element ...
Definition: qdkt.C:593
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()
void beNMatrixOf(const FloatArray &n, int nsd)
Assigns the receiver to be a repeated diagonal matrix.
Definition: floatmatrix.C:505
virtual void giveSurfaceDofMapping(IntArray &answer, int iSurf) const
Assembles surface dof mapping mask, which provides mapping between surface local DOFs and "global" el...
Definition: qdkt.C:667
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: qdkt.C:414
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
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: qdkt.C:697
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
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: qdkt.C:498
virtual IntegrationRule * GetSurfaceIntegrationRule(int iSurf)
Definition: qdkt.C:681
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 int getRequiredNumberOfIntegrationPoints(integrationDomain dType, int approxOrder)
Abstract service.
Class Interface.
Definition: interface.h:82
virtual void computeSurfaceNMatrix(FloatMatrix &answer, int boundaryID, const FloatArray &lcoords)
Computes surface interpolation matrix.
Definition: qdkt.C:659
virtual void giveNodeCoordinates(double &x1, double &x2, double &x3, double &x4, double &y1, double &y2, double &y3, double &y4, double &z1, double &z2, double &z3, double &z4)
Definition: qdkt.C:318
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
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 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
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: qdkt.C:385
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
void updateFringeTableMinMax(double *s, int size)
virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge)
Computes volume related to integration point on local edge.
Definition: qdkt.C:609
Load is base abstract class for all loads.
Definition: load.h:61
The element interface required by LayeredCrossSection.
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
Abstract base class for all structural cross section models.
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.
Class implementing node in finite element mesh.
Definition: node.h:87
int giveNumber() const
Definition: femcmpnn.h:107
double normalize()
Normalizes receiver.
Definition: floatarray.C:828
virtual Interface * giveInterface(InterfaceType it)
Interface requesting service.
Definition: qdkt.C:397
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
virtual void SPRNodalRecoveryMI_giveSPRAssemblyPoints(IntArray &pap)
Definition: qdkt.C:535
#define OOFEG_VARPLOT_PATTERN_LAYER
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType type)
Definition: qdkt.C:746
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
void add(const FloatArray &src)
Adds array src to receiver.
Definition: floatarray.C:156
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
virtual Material * giveMaterial()
Definition: element.C:484
Class representing Gaussian-quadrature integration rule.
virtual void computeSurfaceNMatrixAt(FloatMatrix &answer, int iSurf, GaussPoint *gp)
Definition: qdkt.C:653
virtual void computeStrainVectorInLayer(FloatArray &answer, const FloatArray &masterGpStrain, GaussPoint *masterGp, GaussPoint *slaveGp, TimeStep *tStep)
Computes full 3D strain vector in element layer.
Definition: qdkt.C:571
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