OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
quad1platesubsoil.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 - 2014 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/quad1platesubsoil.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 namespace oofem {
52 REGISTER_Element(Quad1PlateSubSoil);
53 
54 FEI2dQuadLin Quad1PlateSubSoil :: interp_lin(1, 2);
55 
59 {
61  numberOfDofMans = 4;
62 }
63 
64 
67 
68 
71 {
72  return & interp_lin;
73 }
74 
75 
76 void
78 // Sets up the array containing the four Gauss points of the receiver.
79 {
80  if ( integrationRulesArray.size() == 0 ) {
81  integrationRulesArray.resize( 1 );
82  integrationRulesArray [ 0 ].reset( new GaussIntegrationRule(1, this, 1, 5) );
84  }
85 }
86 
87 
88 void
90 {
91  OOFEM_ERROR("Body load not supported, use surface load instead");
92 }
93 
94 
95 void
97 // Returns the [3x4] strain-displacement matrix {B} of the receiver,
98 // evaluated at gp.
99 {
100  FloatArray n;
101  FloatMatrix dn;
102 
105 
106  answer.resize(3, 4);
107  answer.zero();
108 
110  for ( int i = 0; i < 4; ++i ) {
111  answer(0, i) = n(i); // eps_z
112  answer(1, i) = dn(i, 0); // gamma_xz
113  answer(2, i) = dn(i, 1); // gamma_yz
114  }
115 }
116 
117 
118 void
120 {
121  this->giveStructuralCrossSection()->giveGeneralizedStress_PlateSubSoil(answer, gp, strain, tStep);
122 }
123 
124 
125 void
127 {
128  this->giveStructuralCrossSection()->give2dPlateSubSoilStiffMtrx(answer, rMode, gp, tStep);
129 }
130 
131 
134 {
135  this->numberOfGaussPoints = 4;
137 }
138 
139 
140 void
142 {
143  answer = {D_w};
144 }
145 
146 
147 void
149 {
150  FloatArray u, v;
151  u.beDifferenceOf( * this->giveNode(2)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
152  v.beDifferenceOf( * this->giveNode(3)->giveCoordinates(), * this->giveNode(1)->giveCoordinates() );
153 
154  answer.beVectorProductOf(u, v);
155  answer.normalize();
156 }
157 
158 
159 double
161 //
162 // returns receiver's characteristic length for crack band models
163 // for a crack formed in the plane with normal normalToCrackPlane.
164 //
165 {
166  return this->giveCharacteristicLengthForPlaneElements(normalToCrackPlane);
167 }
168 
169 
170 double
172 {
173  double detJ, weight;
174 
175  weight = gp->giveWeight();
177  return detJ * weight;
178 }
179 
180 
181 void
183 // Returns the lumped mass matrix of the receiver.
184 {
185  OOFEM_ERROR("Mass matrix not provided");
186 }
187 
188 
189 int
191 {
192  if ( type == IST_ShellForceTensor ){
193  FloatArray help;
194  answer.resize(6);
195  help = static_cast< StructuralMaterialStatus * >( gp->giveMaterialStatus() )->giveStressVector();
196  answer.at(1) = 0.0; // nx
197  answer.at(2) = 0.0; // ny
198  answer.at(3) = help.at(1); // nz
199  answer.at(4) = help.at(3); // vyz
200  answer.at(5) = help.at(2); // vxz
201  answer.at(6) = 0.0; // vxy
202  return 1;
203  }
204  return StructuralElement :: giveIPValue(answer, gp, type, tStep);
205 }
206 
207 Interface *
209 {
210  if ( interface == ZZNodalRecoveryModelInterfaceType ) {
211  return static_cast< ZZNodalRecoveryModelInterface * >(this);
212  } else if ( interface == SPRNodalRecoveryModelInterfaceType ) {
213  return static_cast< SPRNodalRecoveryModelInterface * >(this);
214  }
215 
216  return NULL;
217 }
218 
219 void
221 {
222  pap.resize(4);
223  for ( int i = 1; i < 5; i++ ) {
224  pap.at(i) = this->giveNode(i)->giveNumber();
225  }
226 }
227 
228 void
230 {
231  int found = 0;
232  answer.resize(1);
233 
234  for ( int i = 1; i < 5; i++ ) {
235  if ( pap == this->giveNode(i)->giveNumber() ) {
236  found = 1;
237  }
238  }
239 
240  if ( found ) {
241  answer.at(1) = pap;
242  } else {
243  OOFEM_ERROR("node unknown");
244  }
245 }
246 
247 
248 void
250 {
251  this->computeNmatrixAt(sgp->giveNaturalCoordinates(), answer);
252 }
253 
254 void
256 {
257  answer.resize(4);
258  answer.zero();
259  if ( iSurf == 1 ) {
260  for (int i = 1; i<=4; i++) {
261  answer.at(i) = i;
262  }
263  } else {
264  OOFEM_ERROR("wrong surface number");
265  }
266 }
267 
270 {
271  IntegrationRule *iRule = new GaussIntegrationRule(1, this, 1, 1);
272  int npoints = iRule->getRequiredNumberOfIntegrationPoints(_Square, approxOrder);
273  iRule->SetUpPointsOnSquare(npoints, _Unknown);
274  return iRule;
275 }
276 
277 double
279 {
280  return this->computeVolumeAround(gp);
281 }
282 
283 
284 int
286 {
287  return 0;
288 }
289 
290 
291 } // end namespace oofem
CrossSection * giveCrossSection()
Definition: element.C:495
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
The element interface required by ZZNodalRecoveryModel.
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
virtual void computeStressVector(FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)
Computes the stress vector of receiver at given integration point, at time step tStep.
virtual void giveSurfaceDofMapping(IntArray &answer, int iSurf) const
Assembles surface dof mapping mask, which provides mapping between surface local DOFs and "global" el...
virtual void computeSurfaceNMatrixAt(FloatMatrix &answer, int iSurf, GaussPoint *gp)
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 double computeSurfaceVolumeAround(GaussPoint *gp, int iSurf)
Computes volume related to integration point on local surface.
virtual void giveGeneralizedStress_PlateSubSoil(FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)=0
virtual double giveTransformationJacobian(const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the determinant of the transformation.
Definition: feinterpol.C:43
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...
The element interface required by ZZNodalRecoveryModel.
void zero()
Sets all component to zero.
Definition: intarray.C:52
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
virtual IntegrationRule * GetSurfaceIntegrationRule(int iSurf)
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in full form.
This class implements a structural material status information.
Definition: structuralms.h:65
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. ...
virtual void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
Computes lumped mass matrix of receiver.
virtual void giveDofManDofIDMask(int inode, IntArray &) const
Returns dofmanager dof mask for node.
virtual int SetUpPointsOnSquare(int, MaterialMode mode)
Sets up receiver&#39;s integration points on unit square integration domain.
virtual void computeMidPlaneNormal(FloatArray &answer, const GaussPoint *gp)
Computes mid-plane normal of receiver at integration point.
virtual void computeNmatrixAt(const FloatArray &iLocCoord, FloatMatrix &answer)
Computes interpolation matrix for element unknowns.
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 Interface * giveInterface(InterfaceType it)
Interface requesting service.
Abstract base class representing integration rule.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
void beDifferenceOf(const FloatArray &a, const FloatArray &b)
Sets receiver to be a - b.
Definition: floatarray.C:341
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
Abstract base class for all "structural" finite elements.
virtual double giveCharacteristicLength(const FloatArray &normalToCrackPlane)
Returns the size of element in the given direction, in some cases adjusted (e.g.
#define OOFEM_ERROR(...)
Definition: error.h:61
Quad1PlateSubSoil(int n, Domain *d)
REGISTER_Element(LSpace)
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
virtual void give2dPlateSubSoilStiffMtrx(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)=0
Method for computing subsoil stiffness matrix for plates.
virtual double giveWeight()
Returns integration weight of receiver.
Definition: gausspoint.h:181
StructuralCrossSection * giveStructuralCrossSection()
Helper function which returns the structural cross-section for the element.
Wrapper around element definition to provide FEICellGeometry interface.
Definition: feinterpol.h:95
void resize(int n)
Checks size of receiver towards requested bounds.
Definition: intarray.C:124
virtual double computeVolumeAround(GaussPoint *gp)
Returns volume related to given integration point.
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
virtual void SPRNodalRecoveryMI_giveSPRAssemblyPoints(IntArray &pap)
virtual void SPRNodalRecoveryMI_giveDofMansDeterminedByPatch(IntArray &answer, int pap)
Class representing vector of real numbers.
Definition: floatarray.h:82
virtual void computeBmatrixAt(GaussPoint *gp, FloatMatrix &answer, int=1, int=ALL_STRAINS)
Computes the geometrical matrix of receiver in given integration point.
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
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
static FEI2dQuadLin interp_lin
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 FEInterpolation * giveInterpolation() const
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
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
InterfaceType
Enumerative type, used to identify interface type.
Definition: interfacetype.h:43
Load is base abstract class for all loads.
Definition: load.h:61
the oofem namespace is to define a context or scope in which all oofem names are defined.
int giveNumber() const
Definition: femcmpnn.h:107
double normalize()
Normalizes receiver.
Definition: floatarray.C:828
Node * giveNode(int i) const
Returns reference to the i-th node of element.
Definition: element.h:610
virtual double evaldNdx(FloatMatrix &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the matrix of derivatives of interpolation functions (shape functions) at given point...
Definition: fei2dquadlin.C:75
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Class representing integration point in finite element program.
Definition: gausspoint.h:93
Class representing solution step.
Definition: timestep.h:80
virtual void computeConstitutiveMatrixAt(FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
Computes constitutive matrix of receiver.
int numberOfDofMans
Number of dofmanagers.
Definition: element.h:149
virtual void evalN(FloatArray &answer, const FloatArray &lcoords, const FEICellGeometry &cellgeo)
Evaluates the array of interpolation functions (shape functions) at given point.
Definition: fei2dquadlin.C:59
const FloatArray & giveNaturalCoordinates()
Returns coordinate array of receiver.
Definition: gausspoint.h:138
Class representing Gaussian-quadrature integration rule.
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:31 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011