OOFEM 3.0
Loading...
Searching...
No Matches
qplanstrss.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 - 2025 Borek Patzak
14 *
15 *
16 *
17 * Czech Technical University, Faculty of Civil Engineering,
18 * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19 *
20 * This library is free software; you can redistribute it and/or
21 * modify it under the terms of the GNU Lesser General Public
22 * License as published by the Free Software Foundation; either
23 * version 2.1 of the License, or (at your option) any later version.
24 *
25 * This program is distributed in the hope that it will be useful,
26 * but WITHOUT ANY WARRANTY; without even the implied warranty of
27 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28 * Lesser General Public License for more details.
29 *
30 * You should have received a copy of the GNU Lesser General Public
31 * License along with this library; if not, write to the Free Software
32 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33 */
34
35#include "sm/Elements/PlaneStress/qplanstrss.h"
36#include "fei2dquadquad.h"
37#include "crosssection.h"
38#include "gausspoint.h"
39#include "gaussintegrationrule.h"
40#include "floatmatrix.h"
41#include "floatarray.h"
42#include "intarray.h"
43#include "mathfem.h"
44#include "classfactory.h"
45
46#ifdef __OOFEG
47 #include "oofeggraphiccontext.h"
48#endif
49
50namespace oofem {
51REGISTER_Element(QPlaneStress2d);
52
53FEI2dQuadQuad QPlaneStress2d :: interpolation(1, 2);
54
55QPlaneStress2d :: QPlaneStress2d(int n, Domain *aDomain) :
56 PlaneStressElement(n, aDomain), ZZNodalRecoveryModelInterface(this), NodalAveragingRecoveryModelInterface(), SpatialLocalizerInterface(this)
57 // Constructor.
58{
59 numberOfDofMans = 8;
60 numberOfGaussPoints = 4;
61}
62
63Interface *
64QPlaneStress2d :: giveInterface(InterfaceType interface)
65{
66 if ( interface == ZZNodalRecoveryModelInterfaceType ) {
67 return static_cast< ZZNodalRecoveryModelInterface * >(this);
68 } else if ( interface == NodalAveragingRecoveryModelInterfaceType ) {
69 return static_cast< NodalAveragingRecoveryModelInterface * >(this);
70 } else if ( interface == SpatialLocalizerInterfaceType ) {
71 return static_cast< SpatialLocalizerInterface * >(this);
72 }
73
74 return NULL;
75}
76
77FEInterpolation *QPlaneStress2d :: giveInterpolation() const { return & interpolation; }
78
79
80
81#ifdef __OOFEG
82void QPlaneStress2d :: drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
83{
84 WCRec p [ 4 ];
85 GraphicObj *go;
86
87 if ( !gc.testElementGraphicActivity(this) ) {
88 return;
89 }
90
91 EASValsSetLineWidth(OOFEG_RAW_GEOMETRY_WIDTH);
92 EASValsSetColor( gc.getElementColor() );
93 EASValsSetEdgeColor( gc.getElementEdgeColor() );
94 EASValsSetEdgeFlag(true);
95 EASValsSetLayer(OOFEG_RAW_GEOMETRY_LAYER);
96 EASValsSetFillStyle(FILL_HOLLOW);
97 p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveCoordinate(1);
98 p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveCoordinate(2);
99 p [ 0 ].z = 0.;
100 p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveCoordinate(1);
101 p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveCoordinate(2);
102 p [ 1 ].z = 0.;
103 p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveCoordinate(1);
104 p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveCoordinate(2);
105 p [ 2 ].z = 0.;
106 p [ 3 ].x = ( FPNum ) this->giveNode(4)->giveCoordinate(1);
107 p [ 3 ].y = ( FPNum ) this->giveNode(4)->giveCoordinate(2);
108 p [ 3 ].z = 0.;
109
110 go = CreateQuad3D(p);
111 EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
112 EGAttachObject(go, ( EObjectP ) this);
113 EMAddGraphicsToModel(ESIModel(), go);
114}
115
116
117void QPlaneStress2d :: drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType type)
118{
119 WCRec p [ 4 ];
120 GraphicObj *go;
121 double defScale = gc.getDefScale();
122
123 if ( !gc.testElementGraphicActivity(this) ) {
124 return;
125 }
126
127 EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
128 EASValsSetColor( gc.getDeformedElementColor() );
129 EASValsSetEdgeColor( gc.getElementEdgeColor() );
130 EASValsSetEdgeFlag(true);
131 EASValsSetLayer(OOFEG_DEFORMED_GEOMETRY_LAYER);
132 EASValsSetFillStyle(FILL_HOLLOW);
133 p [ 0 ].x = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(1, tStep, defScale);
134 p [ 0 ].y = ( FPNum ) this->giveNode(1)->giveUpdatedCoordinate(2, tStep, defScale);
135 p [ 0 ].z = 0.;
136 p [ 1 ].x = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(1, tStep, defScale);
137 p [ 1 ].y = ( FPNum ) this->giveNode(2)->giveUpdatedCoordinate(2, tStep, defScale);
138 p [ 1 ].z = 0.;
139 p [ 2 ].x = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(1, tStep, defScale);
140 p [ 2 ].y = ( FPNum ) this->giveNode(3)->giveUpdatedCoordinate(2, tStep, defScale);
141 p [ 2 ].z = 0.;
142 p [ 3 ].x = ( FPNum ) this->giveNode(4)->giveUpdatedCoordinate(1, tStep, defScale);
143 p [ 3 ].y = ( FPNum ) this->giveNode(4)->giveUpdatedCoordinate(2, tStep, defScale);
144 p [ 3 ].z = 0.;
145
146 go = CreateQuad3D(p);
147 EGWithMaskChangeAttributes(WIDTH_MASK | FILL_MASK | COLOR_MASK | EDGE_COLOR_MASK | EDGE_FLAG_MASK | LAYER_MASK, go);
148 EMAddGraphicsToModel(ESIModel(), go);
149}
150
151
152void QPlaneStress2d :: drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
153{
154 int i, indx, n [ 4 ], result = 0;
155 WCRec p [ 4 ], pp [ 9 ];
156 GraphicObj *tr;
157 FloatArray v [ 8 ];
158 double s [ 9 ], ss [ 4 ], defScale;
159
160 if ( !gc.testElementGraphicActivity(this) ) {
161 return;
162 }
163
164 EASValsSetLayer(OOFEG_VARPLOT_PATTERN_LAYER);
165 if ( gc.giveIntVarMode() == ISM_recovered ) {
166 // ============ plot the recovered values (smoothed data) ===============
167 for ( i = 1; i <= 8; i++ ) {
168 result += this->giveInternalStateAtNode(v [ i - 1 ], gc.giveIntVarType(), gc.giveIntVarMode(), i, tStep);
169 }
170
171 if ( result != 8 ) {
172 return;
173 }
174
175 indx = gc.giveIntVarIndx();
176
177 for ( i = 1; i <= 8; i++ ) {
178 s [ i - 1 ] = v [ i - 1 ].at(indx);
179 }
180
181 // auxiliary value at an added central node
182 // computed as average of the values at all Gauss points
183
184 s [ 8 ] = 0.;
185 for ( GaussPoint *gp: *this->giveDefaultIntegrationRulePtr() ) {
186 if ( giveIPValue(v [ 0 ], gp, gc.giveIntVarType(), tStep) == 0 ) {
187 return;
188 }
189
190 s [ 8 ] += v [ 0 ].at(indx);
191 }
192
193 s [ 8 ] /= integrationRulesArray [ 0 ]->giveNumberOfIntegrationPoints();
194 //s[8] = (s[4]+s[5]+s[6]+s[7])/4.;
195
196 for ( i = 0; i < 8; i++ ) {
197 if ( gc.getInternalVarsDefGeoFlag() ) {
198 // use deformed geometry
199 defScale = gc.getDefScale();
200 pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
201 pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
202 pp [ i ].z = 0.;
203 } else {
204 // use initial geometry
205 pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
206 pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
207 pp [ i ].z = 0.;
208 }
209 }
210
211 pp [ 8 ].x = ( pp [ 4 ].x + pp [ 5 ].x + pp [ 6 ].x + pp [ 7 ].x ) / 4.;
212 pp [ 8 ].y = ( pp [ 4 ].y + pp [ 5 ].y + pp [ 6 ].y + pp [ 7 ].y ) / 4.;
213 pp [ 8 ].z = 0.;
214
215
216 for ( int t = 1; t <= 4; t++ ) {
217 if ( t == 1 ) {
218 n [ 0 ] = 0;
219 n [ 1 ] = 4;
220 n [ 2 ] = 8;
221 n [ 3 ] = 7;
222 } else if ( t == 2 ) {
223 n [ 0 ] = 4;
224 n [ 1 ] = 1;
225 n [ 2 ] = 5;
226 n [ 3 ] = 8;
227 } else if ( t == 3 ) {
228 n [ 0 ] = 5;
229 n [ 1 ] = 2;
230 n [ 2 ] = 6;
231 n [ 3 ] = 8;
232 } else {
233 n [ 0 ] = 6;
234 n [ 1 ] = 3;
235 n [ 2 ] = 7;
236 n [ 3 ] = 8;
237 }
238
239 ss [ 0 ] = s [ n [ 0 ] ];
240 ss [ 1 ] = s [ n [ 1 ] ];
241 ss [ 2 ] = s [ n [ 2 ] ];
242 ss [ 3 ] = s [ n [ 3 ] ];
243
244
245 for ( i = 0; i < 4; i++ ) {
246 p [ i ].x = pp [ n [ i ] ].x;
247 p [ i ].y = pp [ n [ i ] ].y;
248 p [ i ].z = 0.;
249 }
250
251 if ( gc.getScalarAlgo() == SA_ISO_SURF ) {
252 /*
253 * for ( i = 0; i < 4; i++ ) {
254 * if ( gc.getInternalVarsDefGeoFlag() ) {
255 * // use deformed geometry
256 * defScale = gc.getDefScale();
257 * p [ i ].x = ( FPNum ) this->giveNode(n[i] + 1)->giveUpdatedCoordinate(1, tStep, defScale);
258 * p [ i ].y = ( FPNum ) this->giveNode(n[i] + 1)->giveUpdatedCoordinate(2, tStep, defScale);
259 * p [ i ].z = 0.;
260 * } else {
261 * // use initial geometry
262 * p [ i ].x = ( FPNum ) this->giveNode(n[i] + 1)->giveCoordinate(1);
263 * p [ i ].y = ( FPNum ) this->giveNode(n[i] + 1)->giveCoordinate(2);
264 * p [ i ].z = 0.;
265 * }
266 * }
267 */
268 //EASValsSetColor(gc.getYieldPlotColor(ratio));
269 gc.updateFringeTableMinMax(ss, 4);
270 tr = CreateQuadWD3D(p, ss [ 0 ], ss [ 1 ], ss [ 2 ], ss [ 3 ]);
271 EGWithMaskChangeAttributes(LAYER_MASK, tr);
272 EMAddGraphicsToModel(ESIModel(), tr);
273 } else if ( ( gc.getScalarAlgo() == SA_ZPROFILE ) || ( gc.getScalarAlgo() == SA_COLORZPROFILE ) ) {
274 //double landScale = gc.getLandScale();
275
276 for ( i = 0; i < 4; i++ ) {
277 /*
278 * if ( gc.getInternalVarsDefGeoFlag() ) {
279 * // use deformed geometry
280 * defScale = gc.getDefScale();
281 * p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
282 * p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
283 * p [ i ].z = ss [ i ] * landScale;
284 * } else {
285 * // use initial geometry
286 * p [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
287 * p [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
288 * p [ i ].z = ss [ i ] * landScale;
289 * }
290 */
291
292 // this fixes a bug in ELIXIR
293 if ( fabs(ss [ i ]) < 1.0e-6 ) {
294 ss [ i ] = 1.0e-6;
295 }
296 }
297
298 if ( gc.getScalarAlgo() == SA_ZPROFILE ) {
299 EASValsSetColor( gc.getDeformedElementColor() );
300 EASValsSetLineWidth(OOFEG_DEFORMED_GEOMETRY_WIDTH);
301 tr = CreateQuad3D(p);
302 EGWithMaskChangeAttributes(WIDTH_MASK | COLOR_MASK | LAYER_MASK, tr);
303 } else {
304 gc.updateFringeTableMinMax(s, 4);
305 tr = CreateQuadWD3D(p, ss [ 0 ], ss [ 1 ], ss [ 2 ], ss [ 3 ]);
306 EGWithMaskChangeAttributes(LAYER_MASK, tr);
307 }
308
309 EMAddGraphicsToModel(ESIModel(), tr);
310 }
311 }
312 } else if ( gc.giveIntVarMode() == ISM_local ) {
313 // ========== plot the local values (raw data) =====================
314 if ( numberOfGaussPoints != 4 ) {
315 return;
316 }
317
318 IntArray ind(4);
319 WCRec pp [ 9 ];
320
321 for ( i = 0; i < 8; i++ ) {
322 if ( gc.getInternalVarsDefGeoFlag() ) {
323 // use deformed geometry
324 defScale = gc.getDefScale();
325 pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(1, tStep, defScale);
326 pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveUpdatedCoordinate(2, tStep, defScale);
327 pp [ i ].z = 0.;
328 } else {
329 pp [ i ].x = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(1);
330 pp [ i ].y = ( FPNum ) this->giveNode(i + 1)->giveCoordinate(2);
331 pp [ i ].z = 0.;
332 }
333 }
334
335 pp [ 8 ].x = 0.25 * ( pp [ 0 ].x + pp [ 1 ].x + pp [ 2 ].x + pp [ 3 ].x );
336 pp [ 8 ].y = 0.25 * ( pp [ 0 ].y + pp [ 1 ].y + pp [ 2 ].y + pp [ 3 ].y );
337 pp [ 8 ].z = 0.;
338
339 for ( GaussPoint *gp: *this->giveDefaultIntegrationRulePtr() ) {
340 const FloatArray& gpCoords = gp->giveNaturalCoordinates();
341 if ( ( gpCoords.at(1) > 0. ) && ( gpCoords.at(2) > 0. ) ) {
342 ind.at(1) = 0;
343 ind.at(2) = 4;
344 ind.at(3) = 8;
345 ind.at(4) = 7;
346 } else if ( ( gpCoords.at(1) < 0. ) && ( gpCoords.at(2) > 0. ) ) {
347 ind.at(1) = 4;
348 ind.at(2) = 1;
349 ind.at(3) = 5;
350 ind.at(4) = 8;
351 } else if ( ( gpCoords.at(1) < 0. ) && ( gpCoords.at(2) < 0. ) ) {
352 ind.at(1) = 5;
353 ind.at(2) = 2;
354 ind.at(3) = 6;
355 ind.at(4) = 8;
356 } else {
357 ind.at(1) = 6;
358 ind.at(2) = 3;
359 ind.at(3) = 7;
360 ind.at(4) = 8;
361 }
362
363 if ( giveIPValue(v [ 0 ], gp, gc.giveIntVarType(), tStep) == 0 ) {
364 return;
365 }
366
367 indx = gc.giveIntVarIndx();
368
369 for ( i = 1; i <= 4; i++ ) {
370 s [ i - 1 ] = v [ 0 ].at(indx);
371 }
372
373 for ( i = 0; i < 4; i++ ) {
374 p [ i ].x = pp [ ind.at(i + 1) ].x;
375 p [ i ].y = pp [ ind.at(i + 1) ].y;
376 p [ i ].z = pp [ ind.at(i + 1) ].z;
377 }
378
379 gc.updateFringeTableMinMax(s, 4);
380 tr = CreateQuadWD3D(p, s [ 0 ], s [ 1 ], s [ 2 ], s [ 3 ]);
381 EGWithMaskChangeAttributes(LAYER_MASK, tr);
382 EMAddGraphicsToModel(ESIModel(), tr);
383 }
384 }
385}
386#endif
387
388void
389QPlaneStress2d :: NodalAveragingRecoveryMI_computeNodalValue(FloatArray &answer, int node,
390 InternalStateType type, TimeStep *tStep)
391{
392 if ( numberOfGaussPoints != 4 ) {
393 return;
394 }
395
396 GaussPoint *gp;
397
398 if ( node < 5 ) {
399 int i = 0;
400 switch ( node ) {
401 case 1: i = 4;
402 break;
403 case 2: i = 2;
404 break;
405 case 3: i = 1;
406 break;
407 case 4: i = 3;
408 break;
409 }
410
411 gp = integrationRulesArray [ 0 ]->getIntegrationPoint(i - 1);
412 this->giveIPValue(answer, gp, type, tStep);
413 } else {
414 int i1 = 0, i2 = 0;
415 switch ( node ) {
416 case 5: i1 = 4;
417 i2 = 2;
418 break;
419 case 6: i1 = 2;
420 i2 = 1;
421 break;
422 case 7: i1 = 1;
423 i2 = 3;
424 break;
425 case 8: i1 = 3;
426 i2 = 4;
427 break;
428 }
429
430 FloatArray contrib;
431 gp = integrationRulesArray [ 0 ]->getIntegrationPoint(i1 - 1);
432 this->giveIPValue(contrib, gp, type, tStep);
433 gp = integrationRulesArray [ 0 ]->getIntegrationPoint(i2 - 1);
434 this->giveIPValue(answer, gp, type, tStep);
435 answer.add(contrib);
436 answer.times(0.5);
437 }
438}
439
440
441} // end namespace oofem
REGISTER_Element
#define REGISTER_Element(class)
Definition classfactory.h:131
oofem::Element::giveNode
Node * giveNode(int i) const
Definition element.h:629
oofem::Element::numberOfDofMans
int numberOfDofMans
Number of dofmanagers.
Definition element.h:136
oofem::Element::integrationRulesArray
std::vector< std ::unique_ptr< IntegrationRule > > integrationRulesArray
Definition element.h:157
oofem::Element::numberOfGaussPoints
int numberOfGaussPoints
Definition element.h:175
oofem::Element::giveDefaultIntegrationRulePtr
virtual IntegrationRule * giveDefaultIntegrationRulePtr()
Definition element.h:886
oofem::FloatArray::at
double & at(Index i)
Definition floatarray.h:202
oofem::FloatArray::add
void add(const FloatArray &src)
Definition floatarray.C:218
oofem::FloatArray::times
void times(double s)
Definition floatarray.C:834
oofem::IntArray::at
int & at(std::size_t i)
Definition intarray.h:104
oofem::PlaneStressElement::PlaneStressElement
PlaneStressElement(int n, Domain *d)
Definition structural2delement.C:235
oofem::QPlaneStress2d::interpolation
static FEI2dQuadQuad interpolation
Definition qplanstrss.h:56
oofem::SpatialLocalizerInterface::SpatialLocalizerInterface
SpatialLocalizerInterface(Element *element)
Definition spatiallocalizer.h:65
oofem::StructuralElement::giveIPValue
int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep) override
Definition structuralelement.C:1139
oofem::StructuralElement::giveInternalStateAtNode
int giveInternalStateAtNode(FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *tStep) override
Definition structuralelement.C:1256
oofem::ZZNodalRecoveryModelInterface::ZZNodalRecoveryModelInterface
ZZNodalRecoveryModelInterface(Element *element)
Constructor.
Definition zznodalrecoverymodel.h:106
oofem::ZZNodalRecoveryModelInterfaceType
@ ZZNodalRecoveryModelInterfaceType
Definition interfacetype.h:49
oofem::SpatialLocalizerInterfaceType
@ SpatialLocalizerInterfaceType
Definition interfacetype.h:57
oofem::NodalAveragingRecoveryModelInterfaceType
@ NodalAveragingRecoveryModelInterfaceType
Definition interfacetype.h:50
gc
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
OOFEG_VARPLOT_PATTERN_LAYER
#define OOFEG_VARPLOT_PATTERN_LAYER
Definition oofeggraphiccontext.h:105
OOFEG_DEFORMED_GEOMETRY_LAYER
#define OOFEG_DEFORMED_GEOMETRY_LAYER
Definition oofeggraphiccontext.h:102
OOFEG_DEFORMED_GEOMETRY_WIDTH
#define OOFEG_DEFORMED_GEOMETRY_WIDTH
Definition oofeggraphiccontext.h:95
OOFEG_RAW_GEOMETRY_WIDTH
#define OOFEG_RAW_GEOMETRY_WIDTH
Definition oofeggraphiccontext.h:94
OOFEG_RAW_GEOMETRY_LAYER
#define OOFEG_RAW_GEOMETRY_LAYER
Definition oofeggraphiccontext.h:101
This page is part of the OOFEM-3.0 documentation. Copyright Copyright (C) 1994-2025 Borek Patzak Bořek Patzák
Project e-mail: oofem@fsv.cvut.cz
Generated at for OOFEM by doxygen 1.15.0 written by Dimitri van Heesch, © 1997-2011