OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
libeam3dnl.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/Beams/libeam3dnl.h"
36 #include "../sm/CrossSections/structuralcrosssection.h"
37 #include "../sm/Materials/structuralms.h"
38 #include "node.h"
39 #include "material.h"
40 #include "gausspoint.h"
41 #include "gaussintegrationrule.h"
42 #include "floatmatrix.h"
43 #include "intarray.h"
44 #include "floatarray.h"
45 #include "mathfem.h"
46 #include "timestep.h"
47 #include "classfactory.h"
48 
49 #ifdef __OOFEG
50  #include "oofeggraphiccontext.h"
51 #endif
52 
53 namespace oofem {
54 REGISTER_Element(LIBeam3dNL);
55 
56 LIBeam3dNL :: LIBeam3dNL(int n, Domain *aDomain) : NLStructuralElement(n, aDomain), tc(3, 3), tempTc(3, 3) //, kappa (3)
57 {
58  numberOfDofMans = 2;
59  l0 = 0.;
60  tempTcCounter = 0;
61  referenceNode = 0;
62  // init kappa vector at centre
63  // kappa.zero();
64 }
65 
66 
67 void
69 {
70  if ( vec.giveSize() != 3 ) {
71  OOFEM_ERROR("vec param size mismatch");
72  }
73 
74  answer.resize(3, 3);
75 
76  answer.at(1, 1) = answer.at(2, 2) = answer.at(3, 3) = 0.;
77  answer.at(1, 2) = -vec.at(3);
78  answer.at(1, 3) = vec.at(2);
79  answer.at(2, 1) = vec.at(3);
80  answer.at(2, 3) = -vec.at(1);
81  answer.at(3, 1) = -vec.at(2);
82  answer.at(3, 2) = vec.at(1);
83 }
84 
85 
86 void
88 {
89  FloatMatrix S(3, 3), SS(3, 3);
90  double psiSize;
91 
92  if ( psi.giveSize() != 3 ) {
93  OOFEM_ERROR("psi param size mismatch");
94  }
95 
96  answer.resize(3, 3);
97  answer.zero();
98 
99  psiSize = psi.computeNorm();
100  answer.at(1, 1) = answer.at(2, 2) = answer.at(3, 3) = 1.;
101 
102  if ( psiSize <= 1.e-40 ) {
103  return;
104  }
105 
106  this->computeSMtrx(S, psi);
107  SS.beProductOf(S, S);
108  S.times(sin(psiSize) / psiSize);
109  SS.times( ( 1. - cos(psiSize) ) / ( psiSize * psiSize ) );
110 
111  answer.add(S);
112  answer.add(SS);
113 }
114 
115 
116 void
118 {
119  // test if not previously done
120  if ( tStep->giveSolutionStateCounter() == tempTcCounter ) {
121  return;
122  }
123 
124  FloatArray u, centreSpin(3);
125  FloatMatrix dR(3, 3);
126 
127  // ask element's displacement increments
128  this->computeVectorOf(VM_Incremental, tStep, u);
129 
130  // interpolate spin at the centre
131  centreSpin.at(1) = 0.5 * ( u.at(4) + u.at(10) );
132  centreSpin.at(2) = 0.5 * ( u.at(5) + u.at(11) );
133  centreSpin.at(3) = 0.5 * ( u.at(6) + u.at(12) );
134 
135  // compute rotation matrix from centreSpin pseudovector
136  this->computeRotMtrx(dR, centreSpin);
137 
138  // update triad
139  tempTc.beProductOf(dR, tc);
140  // remember timestamp
142 }
143 
144 
145 void
147 {
148  FloatArray xd(3), eps(3), curv(3);
149 
150  // update temp triad
151  this->updateTempTriad(tStep);
152 
153  this->computeXdVector(xd, tStep);
154 
155  // compute eps_xl, gamma_l2, gamma_l3
156  eps.beTProductOf(tempTc, xd);
157  eps.times(1. / this->l0);
158  eps.at(1) -= 1.0;
159 
160  // update curvature at midpoint
161  this->computeTempCurv(curv, tStep);
162 
163  answer.resize(6);
164  answer.at(1) = eps.at(1); // eps_xl
165  answer.at(2) = eps.at(2); // gamma_l2
166  answer.at(3) = eps.at(3); // gamma_l3
167  answer.at(4) = curv.at(1); // kappa_1
168  answer.at(5) = curv.at(2); // kappa_2
169  answer.at(6) = curv.at(3); // kappa_3
170 }
171 
172 
173 void
175 {
176  FloatArray xd(3);
177  FloatMatrix s(3, 3);
178 
179  this->computeXdVector(xd, tStep);
180  this->computeSMtrx(s, xd);
181 
182  answer.resize(12, 6);
183  answer.zero();
184 
185  for ( int i = 1; i < 4; i++ ) {
186  answer.at(i, i) = -1.0;
187  answer.at(i + 6, i) = 1.0;
188  answer.at(i + 3, i + 3) = -1.0;
189  answer.at(i + 9, i + 3) = 1.0;
190 
191  for ( int j = 1; j < 4; j++ ) {
192  answer.at(i + 3, j) = answer.at(i + 9, j) = 0.5 * s.at(j, i);
193  }
194  }
195 }
196 
197 
198 void
199 LIBeam3dNL :: giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord)
200 {
202  FloatArray nm(6), stress, strain;
203  FloatMatrix x;
204  double s1, s2;
205 
206  // update temp triad
207  this->updateTempTriad(tStep);
208 
209  if ( useUpdatedGpRecord == 1 ) {
210  stress = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
211  } else {
212  this->computeStrainVector(strain, gp, tStep);
213  this->computeStressVector(stress, strain, gp, tStep);
214  }
215 
216  for ( int i = 1; i <= 3; i++ ) {
217  s1 = s2 = 0.0;
218  for ( int j = 1; j <= 3; j++ ) {
219  s1 += tempTc.at(i, j) * stress.at(j);
220  s2 += tempTc.at(i, j) * stress.at(j + 3);
221  }
222 
223  nm.at(i) = s1;
224  nm.at(i + 3) = s2;
225  }
226 
227  this->computeXMtrx(x, tStep);
228  answer.beProductOf(x, nm);
229 }
230 
231 
232 void
234 {
235  FloatArray u(3);
236 
237  answer.resize(3);
238  // ask element's displacements
239  this->computeVectorOf(VM_Total, tStep, u);
240 
241  answer.at(1) = ( this->giveNode(2)->giveCoordinate(1) + u.at(7) ) -
242  ( this->giveNode(1)->giveCoordinate(1) + u.at(1) );
243  answer.at(2) = ( this->giveNode(2)->giveCoordinate(2) + u.at(8) ) -
244  ( this->giveNode(1)->giveCoordinate(2) + u.at(2) );
245  answer.at(3) = ( this->giveNode(2)->giveCoordinate(3) + u.at(9) ) -
246  ( this->giveNode(1)->giveCoordinate(3) + u.at(3) );
247 }
248 
249 
250 void
252 {
253  double s1, s2;
254  FloatMatrix d, x, xt(12, 6), dxt, sn, sm, sxd, y;
255  FloatArray n(3), m(3), xd(3), stress, strain;
257 
258  answer.clear();
259 
260  // linear part
261 
262  this->updateTempTriad(tStep); // update temp triad
263  this->computeXMtrx(x, tStep);
264  xt.zero();
265  for ( int i = 1; i <= 12; i++ ) {
266  for ( int j = 1; j <= 3; j++ ) {
267  for ( int k = 1; k <= 3; k++ ) {
268  // compute x*Tbar, taking into account sparsity of Tbar
269  xt.at(i, j) += x.at(i, k) * tempTc.at(k, j);
270  xt.at(i, j + 3) += x.at(i, k + 3) * tempTc.at(k, j);
271  }
272  }
273  }
274 
275  this->computeConstitutiveMatrixAt(d, rMode, gp, tStep);
276  dxt.beProductTOf(d, xt);
277  answer.beProductOf(xt, dxt);
278  answer.times(1. / this->l0);
279 
280  // geometric stiffness ks = ks1+ks2
281  // ks1
282  this->computeStrainVector(strain, gp, tStep);
283  this->computeStressVector(stress, strain, gp, tStep);
284 
285  for ( int i = 1; i <= 3; i++ ) {
286  s1 = s2 = 0.0;
287  for ( int j = 1; j <= 3; j++ ) {
288  s1 += tempTc.at(i, j) * stress.at(j);
289  s2 += tempTc.at(i, j) * stress.at(j + 3);
290  }
291 
292  n.at(i) = s1;
293  m.at(i) = s2;
294  }
295 
296  this->computeSMtrx(sn, n);
297  this->computeSMtrx(sm, m);
298 
299  for ( int i = 1; i <= 3; i++ ) {
300  for ( int j = 1; j <= 3; j++ ) {
301  answer.at(i, j + 3) += sn.at(i, j);
302  answer.at(i, j + 9) += sn.at(i, j);
303  answer.at(i + 3, j + 3) += sm.at(i, j);
304  answer.at(i + 3, j + 9) += sm.at(i, j);
305 
306  answer.at(i + 6, j + 3) -= sn.at(i, j);
307  answer.at(i + 6, j + 9) -= sn.at(i, j);
308  answer.at(i + 9, j + 3) -= sm.at(i, j);
309  answer.at(i + 9, j + 9) -= sm.at(i, j);
310  }
311  }
312 
313  // ks2
314  this->computeXdVector(xd, tStep);
315  this->computeSMtrx(sxd, xd);
316 
317  y.beProductOf(sxd, sn);
318  y.times(0.5);
319 
320  for ( int i = 1; i <= 3; i++ ) {
321  for ( int j = 1; j <= 3; j++ ) {
322  answer.at(i + 3, j) -= sn.at(i, j);
323  answer.at(i + 3, j + 3) += y.at(i, j);
324  answer.at(i + 3, j + 6) += sn.at(i, j);
325  answer.at(i + 3, j + 9) += y.at(i, j);
326 
327  answer.at(i + 9, j) -= sn.at(i, j);
328  answer.at(i + 9, j + 3) += y.at(i, j);
329  answer.at(i + 9, j + 6) += sn.at(i, j);
330  answer.at(i + 9, j + 9) += y.at(i, j);
331  }
332  }
333 }
334 
335 
336 void
338 // Sets up the array of Gauss Points of the receiver.
339 {
340  if ( integrationRulesArray.size() == 0 ) {
341  integrationRulesArray.resize( 1 );
342  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 2) );
344  }
345 }
346 
347 
350 {
351  IRResultType result; // Required by IR_GIVE_FIELD macro
352 
353  // first call parent
355  if ( result != IRRT_OK ) {
356  return result;
357  }
358 
360  if ( referenceNode == 0 ) {
361  OOFEM_ERROR("wrong reference node specified");
362  }
363 
364  /*
365  * if (this->hasString (initString, "dofstocondense")) {
366  * dofsToCondense = this->ReadIntArray (initString, "dofstocondense");
367  * if (dofsToCondense->giveSize() >= 12)
368  * OOFEM_ERROR("wrong input data for condensed dofs");
369  * } else {
370  * dofsToCondense = NULL;
371  * }
372  */
373 
375  // compute initial triad at centre - requires nodal coordinates
376  FloatMatrix lcs;
377  this->giveLocalCoordinateSystem(lcs);
378  this->tc.beTranspositionOf(lcs);
379 
380  return IRRT_OK;
381 }
382 
383 
384 double
386 // Returns the original length (l0) of the receiver.
387 {
388  double dx, dy, dz;
389  Node *nodeA, *nodeB;
390 
391  if ( l0 == 0. ) {
392  nodeA = this->giveNode(1);
393  nodeB = this->giveNode(2);
394  dx = nodeB->giveCoordinate(1) - nodeA->giveCoordinate(1);
395  dy = nodeB->giveCoordinate(2) - nodeA->giveCoordinate(2);
396  dz = nodeB->giveCoordinate(3) - nodeA->giveCoordinate(3);
397  l0 = sqrt(dx * dx + dy * dy + dz * dz);
398  }
399 
400  return l0;
401 }
402 
403 
404 void
406 // Returns the lumped mass matrix of the receiver. This expression is
407 // valid in both local and global axes.
408 {
409  GaussPoint *gp = integrationRulesArray [ 0 ]->getIntegrationPoint(0);
410  double density = this->giveStructuralCrossSection()->give('d', gp);
411  double halfMass = density * this->giveCrossSection()->give(CS_Area, gp) * this->computeLength() / 2.;
412  answer.resize(12, 12);
413  answer.zero();
414  answer.at(1, 1) = answer.at(2, 2) = answer.at(3, 3) = halfMass;
415  answer.at(7, 7) = answer.at(8, 8) = answer.at(9, 9) = halfMass;
416 
417  double Ik, Iy, Iz;
418  Ik = this->giveCrossSection()->give(CS_TorsionMomentX, gp);
419  Iy = this->giveCrossSection()->give(CS_InertiaMomentY, gp);
420  Iz = this->giveCrossSection()->give(CS_InertiaMomentZ, gp);
421  halfMass = density * this->computeLength() / 2.;
422  answer.at(4, 4) = answer.at(10, 10) = Ik * halfMass;
423  answer.at(5, 5) = answer.at(11, 11) = Iy * halfMass;
424  answer.at(6, 6) = answer.at(12, 12) = Iz * halfMass;
425 }
426 
427 
428 void
430 // Returns the displacement interpolation matrix {N} of the receiver, eva-
431 // luated at gp.
432 {
433  double ksi, n1, n2;
434 
435  ksi = iLocCoord.at(1);
436  n1 = ( 1. - ksi ) * 0.5;
437  n2 = ( 1. + ksi ) * 0.5;
438 
439  answer.resize(6, 12);
440  answer.zero();
441 
442  // u
443  answer.at(1, 1) = n1;
444  answer.at(1, 7) = n2;
445  // v
446  answer.at(2, 2) = n1;
447  answer.at(2, 8) = n2;
448  // w
449  answer.at(3, 3) = n1;
450  answer.at(3, 9) = n2;
451  // fi_x
452  answer.at(4, 4) = n1;
453  answer.at(4, 10) = n2;
454  // fi_y
455  answer.at(5, 5) = n1;
456  answer.at(5, 11) = n2;
457  // fi_z
458  answer.at(6, 6) = n1;
459  answer.at(6, 12) = n2;
460 }
461 
462 
463 double
465 // Returns the length of the receiver. This method is valid only if 1
466 // Gauss point is used.
467 {
468  double weight = gp->giveWeight();
469  return weight * 0.5 * this->computeLength();
470 }
471 
472 
473 int
475 {
478  if ( type == IST_BeamForceMomentTensor ) {
479  answer = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
480  return 1;
481  } else if ( type == IST_BeamStrainCurvatureTensor ) {
482  answer = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStrainVector();
483  return 1;
484  } else {
485  return StructuralElement :: giveIPValue(answer, gp, type, tStep);
486  }
487 }
488 
489 
490 void
492 {
493  answer = {D_u, D_v, D_w, R_u, R_v, R_w};
494 }
495 
496 
497 int
499 {
500  double ksi, n1, n2;
501 
502  ksi = lcoords.at(1);
503  n1 = ( 1. - ksi ) * 0.5;
504  n2 = ( 1. + ksi ) * 0.5;
505 
506  answer.resize(3);
507  answer.at(1) = n1 * this->giveNode(1)->giveCoordinate(1) + n2 *this->giveNode(2)->giveCoordinate(1);
508  answer.at(2) = n1 * this->giveNode(1)->giveCoordinate(2) + n2 *this->giveNode(2)->giveCoordinate(2);
509  answer.at(3) = n1 * this->giveNode(1)->giveCoordinate(3) + n2 *this->giveNode(2)->giveCoordinate(3);
510 
511  return 1;
512 }
513 
514 
515 void
517 {
518  this->giveStructuralCrossSection()->give3dBeamStiffMtrx(answer, rMode, gp, tStep);
519 }
520 
521 
522 void
524 {
525  this->giveStructuralCrossSection()->giveGeneralizedStress_Beam3d(answer, gp, strain, tStep);
526 }
527 
528 
529 void
531 {
532  /*
533  * provides dof mapping of local edge dofs (only nonzero are taken into account)
534  * to global element dofs
535  */
536  if ( iEdge != 1 ) {
537  OOFEM_ERROR("wrong edge number");
538  }
539 
540  answer.resize(12);
541  for ( int i = 1; i <= 12; i++ ) {
542  answer.at(i) = i;
543  }
544 }
545 
546 
547 double
549 {
550  if ( iEdge != 1 ) { // edge between nodes 1 2
551  OOFEM_ERROR("wrong egde number");
552  }
553 
554  double weight = gp->giveWeight();
555  return 0.5 * this->computeLength() * weight;
556 }
557 
558 
559 int
561 //
562 // returns a unit vectors of local coordinate system at element
563 // stored rowwise (mainly used by some materials with ortho and anisotrophy)
564 //
565 {
566  FloatArray lx(3), ly(3), lz(3), help(3);
567  double length = this->computeLength();
568  Node *nodeA, *nodeB, *refNode;
569 
570  answer.resize(3, 3);
571  answer.zero();
572  nodeA = this->giveNode(1);
573  nodeB = this->giveNode(2);
574  refNode = this->giveDomain()->giveNode(this->referenceNode);
575 
576  for ( int i = 1; i <= 3; i++ ) {
577  lx.at(i) = ( nodeB->giveCoordinate(i) - nodeA->giveCoordinate(i) ) / length;
578  help.at(i) = ( refNode->giveCoordinate(i) - nodeA->giveCoordinate(i) );
579  }
580 
581  lz.beVectorProductOf(lx, help);
582  lz.normalize();
583  ly.beVectorProductOf(lz, lx);
584  ly.normalize();
585 
586  for ( int i = 1; i <= 3; i++ ) {
587  answer.at(1, i) = lx.at(i);
588  answer.at(2, i) = ly.at(i);
589  answer.at(3, i) = lz.at(i);
590  }
591 
592  return 1;
593 }
594 
595 
596 /*
597  * int
598  * LIBeam3dNL :: computeGtoLRotationMatrix (FloatMatrix& answer)
599  * // Returns the rotation matrix of the receiver.
600  * // rotation matrix is computed for original configuration only
601  * {
602  * FloatMatrix lcs;
603  * int i,j;
604  *
605  * answer.resize(12,12);
606  * answer.zero();
607  *
608  * this->giveLocalCoordinateSystem (lcs);
609  *
610  * for (i=1; i <= 3; i++)
611  * for (j=1; j <= 3; j++) {
612  * answer.at(i,j) = lcs.at(i,j);
613  * answer.at(i+3, j+3) = lcs.at(i,j);
614  * answer.at(i+6, j+6) = lcs.at(i,j);
615  * answer.at(i+9, j+9) = lcs.at(i,j);
616  * }
617  *
618  * return 1 ;
619  * }
620  */
621 
622 int
624 {
625  /*
626  * Returns transformation matrix from global coordinate system to local
627  * element coordinate system for element load vector components.
628  * If no transformation is necessary, answer is empty matrix (default);
629  *
630  * Does not support follower load
631  */
632 
633  FloatMatrix lcs;
634 
635  answer.resize(6, 6);
636  answer.zero();
637 
638  this->giveLocalCoordinateSystem(lcs);
639 
640  for ( int i = 1; i <= 3; i++ ) {
641  for ( int j = 1; j <= 3; j++ ) {
642  answer.at(i, j) = lcs.at(i, j);
643  answer.at(3 + i, 3 + j) = lcs.at(i, j);
644  }
645  }
646 
647  return 1;
648 }
649 
650 
651 int
653 {
654  // returns transformation matrix from
655  // edge local coordinate system
656  // to element local c.s
657  // (same as global c.s in this case)
658  //
659  // i.e. f(element local) = T * f(edge local)
660  //
661  answer.clear();
662  return 0;
663 }
664 
665 
666 void
668 {
669  FloatArray lc(1);
670  NLStructuralElement :: computeBodyLoadVectorAt(answer, load, tStep, mode);
671  answer.times( this->giveCrossSection()->give(CS_Area, lc, this) );
672 }
673 
674 
675 
677 // Updates the receiver at end of step.
678 {
680 
681  // update triad
682  this->updateTempTriad(tStep);
683  tc = tempTc;
684 
685  // update curvature
686  //FloatArray curv;
687  //this->computeTempCurv (curv, tStep);
688  //kappa = curv;
689 }
690 
691 void
693 // initializes receiver to new time step or can be used
694 // if current time step must be restarted
695 {
697  tempTc = tc;
698 }
699 
700 
701 void
703 {
705  GaussPoint *gp = iRule->getIntegrationPoint(0);
706 
707  FloatArray ui(3), ac(3), PrevEpsilon;
708  FloatMatrix sc(3, 3), tmid(3, 3);
709 
710  // update curvature at midpoint
711  // first, compute Tmid
712  // ask increments
713  this->computeVectorOf(VM_Incremental, tStep, ui);
714 
715  ac.at(1) = 0.5 * ( ui.at(10) - ui.at(4) );
716  ac.at(2) = 0.5 * ( ui.at(11) - ui.at(5) );
717  ac.at(3) = 0.5 * ( ui.at(12) - ui.at(6) );
718  this->computeSMtrx(sc, ac);
719  sc.times(1. / 2.);
720  // compute I+sc
721  sc.at(1, 1) += 1.0;
722  sc.at(2, 2) += 1.0;
723  sc.at(3, 3) += 1.0;
724  tmid.beProductOf(sc, this->tc);
725 
726  // update curvature at centre
727  ac.at(1) = ( ui.at(10) - ui.at(4) );
728  ac.at(2) = ( ui.at(11) - ui.at(5) );
729  ac.at(3) = ( ui.at(12) - ui.at(6) );
730 
731  answer.beTProductOf(tmid, ac);
732  answer.times(1 / this->l0);
733 
734  // ask for previous kappa
735  StructuralMaterialStatus *matStat = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() );
736  if ( matStat ) {
737  PrevEpsilon = matStat->giveStrainVector();
738  }
739  if ( PrevEpsilon.giveSize() ) {
740  answer.at(1) += PrevEpsilon.at(4);
741  answer.at(2) += PrevEpsilon.at(5);
742  answer.at(3) += PrevEpsilon.at(6);
743  }
744 }
745 
746 
747 #ifdef __OOFEG
749 {
750  GraphicObj *go;
751 
752  if ( !gc.testElementGraphicActivity(this) ) {
753  return;
754  }
755 
756  // if (!go) { // create new one
757  WCRec p [ 2 ]; /* poin */
758  EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
759  EASValsSetColor( gc.getElementColor() );
760  EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
761  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
762  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveCoordinate(2);
763  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveCoordinate(3);
764  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
765  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveCoordinate(2);
766  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveCoordinate(3);
767  go = CreateLine3D(p);
768  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
769  EGAttachObject(go, ( EObjectP ) this);
770  EMAddGraphicsToModel(ESIModel(), go);
771 }
772 
773 
775 {
776  GraphicObj *go;
777 
778  if ( !gc.testElementGraphicActivity(this) ) {
779  return;
780  }
781 
782  double defScale = gc.getDefScale();
783  // if (!go) { // create new one
784  WCRec p [ 2 ]; /* poin */
785  EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
786  EASValsSetColor( gc.getDeformedElementColor() );
787  EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
788  p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
789  p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale);
790  p [ 0 ].z = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(3, tStep, defScale);
791 
792  p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
793  p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale);
794  p [ 1 ].z = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(3, tStep, defScale);
795  go = CreateLine3D(p);
796  EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, go);
797  EMAddGraphicsToModel(ESIModel(), go);
798 }
799 #endif
800 } // 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 int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
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: libeam3dnl.C:652
Class and object Domain.
Definition: domain.h:115
void computeVectorOf(ValueModeType u, TimeStep *tStep, FloatArray &answer)
Returns local vector of unknowns.
Definition: element.C:86
virtual IntegrationRule * giveDefaultIntegrationRulePtr()
Access method for default integration rule.
Definition: element.h:822
virtual void giveInternalForcesVector(FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord=0)
Evaluates nodal representation of real internal forces.
Definition: libeam3dnl.C:199
FloatMatrix tc
Last equilibrium triad at the centre.
Definition: libeam3dnl.h:60
virtual int computeLoadGToLRotationMtrx(FloatMatrix &answer)
Returns transformation matrix from global coordinate system to local element coordinate system for el...
Definition: libeam3dnl.C:623
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Definition: libeam3dnl.C:337
Abstract base class for "structural" finite elements with geometrical nonlinearities.
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
Definition: libeam3dnl.C:405
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
double l0
Initial length.
Definition: libeam3dnl.h:58
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
virtual void initForNewStep()
Initializes receivers state to new time step.
Definition: libeam3dnl.C:692
#define OOFEG_RAW_GEOMETRY_LAYER
This class implements a structural material status information.
Definition: structuralms.h:65
void computeRotMtrx(FloatMatrix &answer, FloatArray &psi)
Evaluates the rotation matrix for large rotations according to Rodrigues formula for given pseudovect...
Definition: libeam3dnl.C:87
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
Definition: libeam3dnl.C:516
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual int giveLocalCoordinateSystem(FloatMatrix &answer)
Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy...
Definition: libeam3dnl.C:560
FloatMatrix tempTc
Temporary triad at the centre.
Definition: libeam3dnl.h:64
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
Definition: libeam3dnl.C:429
virtual void initForNewStep()
Initializes receivers state to new time step.
Definition: element.C:846
virtual void giveGeneralizedStress_Beam3d(FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)=0
Moment of inertia around z-axis.
Definition: crosssection.h:64
#define S(p)
Definition: mdm.C:481
virtual double giveCoordinate(int i)
Definition: node.C:82
virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
Definition: libeam3dnl.C:774
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 void updateYourself(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
Definition: libeam3dnl.C:676
#define OOFEG_DEFORMED_GEOMETRY_LAYER
virtual void giveEdgeDofMapping(IntArray &answer, int iEdge) const
Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element ...
Definition: libeam3dnl.C:530
Abstract base class representing integration rule.
virtual int computeGlobalCoordinates(FloatArray &answer, const FloatArray &lcoords)
Computes the global coordinates from given element&#39;s local coordinates.
Definition: libeam3dnl.C:498
virtual void computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
Computes the stiffness matrix of receiver.
Definition: libeam3dnl.C:251
#define _IFT_LIBeam3dNL_refnode
Definition: libeam3dnl.h:43
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. ...
void updateTempTriad(TimeStep *tStep)
Updates the temporary triad at the centre to the state identified by given solution step...
Definition: libeam3dnl.C:117
StateCounterType giveSolutionStateCounter()
Returns current solution state counter.
Definition: timestep.h:188
#define OOFEM_ERROR(...)
Definition: error.h:61
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
Definition: libeam3dnl.C:474
virtual double computeLength()
Computes the length (zero for all but 1D geometries)
Definition: libeam3dnl.C:385
REGISTER_Element(LSpace)
#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
Moment of inertia around x-axis.
Definition: crosssection.h:65
void computeTempCurv(FloatArray &answer, TimeStep *tStep)
Compute the temporary curvature at the centre to the state identified by given solution step...
Definition: libeam3dnl.C:702
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
virtual double giveUpdatedCoordinate(int ic, TimeStep *tStep, double scale=1.)
Returns updated ic-th coordinate of receiver.
Definition: node.C:245
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 give3dBeamStiffMtrx(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)=0
Computes the stiffness matrix for 2d beams.
void computeSMtrx(FloatMatrix &answer, FloatArray &vec)
Evaluates the S matrix from given vector vec.
Definition: libeam3dnl.C:68
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: libeam3dnl.C:349
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
Definition: libeam3dnl.C:748
virtual int setupIntegrationPoints(IntegrationRule &irule, int npoints, Element *element)
Sets up integration rule for the given element.
Definition: crosssection.C:54
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 void updateYourself(TimeStep *tStep)
Updates element state after equilibrium in time step has been reached.
IntegrationPointStatus * giveMaterialStatus()
Returns reference to associated material status (NULL if not defined).
Definition: gausspoint.h:205
LIBeam3dNL(int n, Domain *d)
Definition: libeam3dnl.C:56
double computeNorm() const
Computes the norm (or length) of the vector.
Definition: floatarray.C:840
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 add(const FloatMatrix &a)
Adds matrix to the receiver.
Definition: floatmatrix.C:1023
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
Definition: libeam3dnl.C:464
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: libeam3dnl.C:523
virtual void computeStrainVector(FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
Compute strain vector of receiver evaluated at given integration point at given time step from elemen...
Definition: libeam3dnl.C:146
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
Definition: libeam3dnl.C:491
void times(double s)
Multiplies receiver with scalar.
Definition: floatarray.C:818
void beTranspositionOf(const FloatMatrix &src)
Assigns to the receiver the transposition of parameter.
Definition: floatmatrix.C:323
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: libeam3dnl.C:667
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.
void computeXdVector(FloatArray &answer, TimeStep *tStep)
Computes x_21&#39; vector for given solution state.
Definition: libeam3dnl.C:233
void beProductTOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:397
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
Domain * giveDomain() const
Definition: femcmpnn.h:100
Load is base abstract class for all loads.
Definition: load.h:61
int referenceNode
Reference node.
Definition: libeam3dnl.h:68
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:337
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
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.
virtual double computeEdgeVolumeAround(GaussPoint *gp, int iEdge)
Computes volume related to integration point on local edge.
Definition: libeam3dnl.C:548
void clear()
Sets size of receiver to be an empty matrix. It will have zero rows and zero columns size...
Definition: floatmatrix.h:516
Class implementing node in finite element mesh.
Definition: node.h:87
#define IR_GIVE_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:69
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
const FloatArray & giveStrainVector() const
Returns the const pointer to receiver&#39;s strain vector.
Definition: structuralms.h:105
Moment of inertia around y-axis.
Definition: crosssection.h:63
Class representing integration point in finite element program.
Definition: gausspoint.h:93
StateCounterType tempTcCounter
Time stamp of temporary centre triad.
Definition: libeam3dnl.h:66
Class representing solution step.
Definition: timestep.h:80
int numberOfDofMans
Number of dofmanagers.
Definition: element.h:149
Class representing Gaussian-quadrature integration rule.
void computeXMtrx(FloatMatrix &answer, TimeStep *tStep)
Computes X mtrx at given solution state.
Definition: libeam3dnl.C:174
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