OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
prescribedgenstrainshell7.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 
36 #include "dofiditem.h"
37 #include "dofmanager.h"
38 #include "dof.h"
39 #include "valuemodetype.h"
40 #include "floatarray.h"
41 #include "floatmatrix.h"
42 #include "function.h"
43 #include "engngm.h"
44 #include "set.h"
45 #include "node.h"
46 #include "element.h"
47 #include "classfactory.h"
48 #include "dynamicinputrecord.h"
49 #include "feinterpol.h"
50 #include "sparsemtrx.h"
51 #include "sparselinsystemnm.h"
53 
54 namespace oofem {
55 REGISTER_BoundaryCondition(PrescribedGenStrainShell7);
56 
57 double PrescribedGenStrainShell7 :: give(Dof *dof, ValueModeType mode, double time)
58 {
59  DofIDItem id = dof->giveDofID();
60  FloatArray *coords = dof->giveDofManager()->giveCoordinates();
61 
62  if ( coords->giveSize() != this->centerCoord.giveSize() ) {
63  OOFEM_ERROR("PrescribedGenStrainShell7 :: give - Size of coordinate system different from center coordinate in b.c.");
64  }
65 
66  double factor = 0;
67  if ( mode == VM_Total ) {
68  factor = this->giveTimeFunction()->evaluateAtTime(time);
69  } else if ( mode == VM_Velocity ) {
70  factor = this->giveTimeFunction()->evaluateVelocityAtTime(time);
71  } else if ( mode == VM_Acceleration ) {
72  factor = this->giveTimeFunction()->evaluateAccelerationAtTime(time);
73  } else {
74  OOFEM_ERROR("Should not be called for value mode type then total, velocity, or acceleration.");
75  }
76 
77  // Reminder: u_i = F_ij . (x_j - xb_j) = H_ij . dx_j
78  FloatArray dx;
79  dx.beDifferenceOf(* coords, this->centerCoord);
80 
81  // Assuming the coordinate system to be local, dx(3) = z
82  this->setDeformationGradient( dx.at(3) );
83 
84  FloatArray u, u2;
85  u.beProductOf(gradient, dx);
86 
87  // Add second order contribution, note only higher order in the thickness direction
88  this->evaluateHigherOrderContribution(u2, dx.at(3), dx);
89  u.add(u2);
90 
91  u.times( factor );
92 
93  switch ( id ) {
94  case D_u:
95  case V_u:
96  return u.at(1);
97 
98  case D_v:
99  case V_v:
100  return u.at(2);
101 
102  case D_w:
103  case V_w:
104  return u.at(3);
105 
106  default:
107  return 0.0;
108  }
109 }
110 
111 
112 
113 void
115 {
116  // Evaluates the covariant base vectors in the current configuration
117  FloatArray g1; FloatArray g2; FloatArray g3;
118 
119  FloatArray dxdxi1, dxdxi2, m, dmdxi1, dmdxi2;
120  double dgamdxi1, dgamdxi2, gam;
121  Shell7Base :: giveGeneralizedStrainComponents(genEps, dxdxi1, dxdxi2, dmdxi1, dmdxi2, m, dgamdxi1, dgamdxi2, gam);
122  double fac1 = ( zeta + 0.5 * gam * zeta * zeta );
123  double fac2 = ( 0.5 * zeta * zeta );
124  double fac3 = ( 1.0 + zeta * gam );
125 
126  g1 = dxdxi1 + fac1*dmdxi1 + fac2*dgamdxi1*m;
127  g2 = dxdxi2 + fac1*dmdxi2 + fac2*dgamdxi2*m;
128  g3 = fac3*m;
129  gcov.resize(3,3);
130  gcov.setColumn(g1,1); gcov.setColumn(g2,2); gcov.setColumn(g3,3);
131 }
132 
133 
134 void
136 {
137  // Evaluates the initial base vectors given the array of generalized strain
138  FloatArray G1(3), G2(3), G3(3);
139 
140  G1.at(1) = genEps.at(1) + zeta * genEps.at(7);
141  G1.at(2) = genEps.at(2) + zeta * genEps.at(8);
142  G1.at(3) = genEps.at(3) + zeta * genEps.at(9);
143 
144  G2.at(1) = genEps.at(4) + zeta * genEps.at(10);
145  G2.at(2) = genEps.at(5) + zeta * genEps.at(11);
146  G2.at(3) = genEps.at(6) + zeta * genEps.at(12);
147 
148  G3.at(1) = genEps.at(13);
149  G3.at(2) = genEps.at(14);
150  G3.at(3) = genEps.at(15);
151 
152 
153  Gcov.resize(3,3);
154  Gcov.setColumn(G1,1); Gcov.setColumn(G2,2); Gcov.setColumn(G3,3);
155 }
156 
157 
158 void
160 {
161  // Computes the deformation gradient in matrix form as open product(g_i, G^i) = gcov*Gcon^T
162  FloatMatrix gcov, Gcon, Gcov;
163 
164  this->evalCovarBaseVectorsAt(gcov, this->genEps, zeta);
165  this->evalInitialCovarBaseVectorsAt(Gcov, this->initialGenEps, zeta);
166  Shell7Base :: giveDualBase(Gcov, Gcon);
167 
168  this->gradient.beProductTOf(gcov, Gcon);
169  this->gradient.at(1,1) -= 1.0;
170  this->gradient.at(2,2) -= 1.0;
171  this->gradient.at(3,3) -= 1.0;
172 }
173 
174 void
176 {
177  // Computes the higher order contribution from the second gradient F2_ijk = (g3,3)_i * (G^3)_j * (G^3)_k
178  // Simplified version with only contribtion in the xi-direction
179  FloatMatrix gcov, Gcon, Gcov;
180 
181  this->evalCovarBaseVectorsAt(gcov, this->genEps, zeta);
182  this->evalInitialCovarBaseVectorsAt(Gcov, this->initialGenEps, zeta);
183  Shell7Base :: giveDualBase(Gcov, Gcon);
184 
185  FloatArray G3(3), g3prime(3), m(3);
186  G3.at(1) = Gcon.at(1,3);
187  G3.at(2) = Gcon.at(2,3);
188  G3.at(3) = Gcon.at(3,3);
189 
190  double factor = G3.dotProduct(dx);
191  double gamma = this->genEps.at(18);
192  m.at(1) = this->genEps.at(13);
193  m.at(2) = this->genEps.at(14);
194  m.at(3) = this->genEps.at(15);
195  g3prime = gamma*m;
196 
197  answer = 0.5*factor*factor * g3prime;
198 
199 
200 }
201 
202 #if 0
204 // This is written in a very general way, supporting both fm and sm problems.
205 // v_prescribed = C.d = (x-xbar).d;
206 // C = [x 0 y]
207 // [0 y x]
208 // [ ... ] in 2D, voigt form [d_11, d_22, d_12]
209 // C = [x 0 0 y z 0]
210 // [0 y 0 x 0 z]
211 // [0 0 z 0 x y]
212 // [ ......... ] in 3D, voigt form [d_11, d_22, d_33, d_23, d_13, d_12]
213 {
214  Domain *domain = this->giveDomain();
215  int nNodes = domain->giveNumberOfDofManagers();
216 
217  int nsd = domain->giveNumberOfSpatialDimensions();
219  C.resize(npeq, nsd * ( nsd + 1 ) / 2);
220  C.zero();
221 
222  FloatArray &cCoords = this->giveCenterCoordinate();
223  double xbar = cCoords.at(1), ybar = cCoords.at(2), zbar = 0.0;
224  if ( nsd == 3 ) {
225  zbar = cCoords.at(3);
226  }
227 
228  for ( int i = 1; i <= nNodes; i++ ) {
229  Node *n = domain->giveNode(i);
230  FloatArray *coords = n->giveCoordinates();
231  Dof *d1 = n->giveDofWithID( this->dofs(0) );
232  Dof *d2 = n->giveDofWithID( this->dofs(1) );
233  int k1 = d1->__givePrescribedEquationNumber();
234  int k2 = d2->__givePrescribedEquationNumber();
235  if ( nsd == 2 ) {
236  if ( k1 ) {
237  C.at(k1, 1) = coords->at(1) - xbar;
238  C.at(k1, 3) = coords->at(2) - ybar;
239  }
240 
241  if ( k2 ) {
242  C.at(k2, 2) = coords->at(2) - ybar;
243  C.at(k2, 3) = coords->at(1) - xbar;
244  }
245  } else { // nsd == 3
246  OOFEM_ERROR("PrescribedGenStrainShell7 :: updateCoefficientMatrix - 3D Not tested yet!");
247  Dof *d3 = n->giveDofWithID( this->dofs(2) );
248  int k3 = d3->__givePrescribedEquationNumber();
249 
250  if ( k1 ) {
251  C.at(k1, 1) = coords->at(1) - xbar;
252  C.at(k1, 4) = coords->at(2) - ybar;
253  C.at(k1, 5) = coords->at(3) - zbar;
254  }
255  if ( k2 ) {
256  C.at(k2, 2) = coords->at(2) - ybar;
257  C.at(k2, 4) = coords->at(1) - xbar;
258  C.at(k2, 6) = coords->at(3) - zbar;
259  }
260  if ( k3 ) {
261  C.at(k3, 3) = coords->at(3) - zbar;
262  C.at(k3, 5) = coords->at(1) - xbar;
263  C.at(k3, 6) = coords->at(2) - ybar;
264  }
265  }
266  }
267 }
268 #endif
269 
271 {
272  int nsd = this->domain->giveNumberOfSpatialDimensions();
273  double domain_size = 0.0;
274  // This requires the boundary to be consistent and ordered correctly.
275  Set *set = this->giveDomain()->giveSet(this->set);
276  const IntArray &boundaries = set->giveBoundaryList();
277 
278  for ( int pos = 1; pos <= boundaries.giveSize() / 2; ++pos ) {
279  Element *e = this->giveDomain()->giveElement( boundaries.at(pos * 2 - 1) );
280  int boundary = boundaries.at(pos * 2);
282  domain_size += fei->evalNXIntegral( boundary, FEIElementGeometryWrapper(e) );
283  }
284  return domain_size / nsd;
285 }
286 
287 
288 
289 
291 {
292  IRResultType result; // Required by IR_GIVE_FIELD macro
293 
296 
298  this->centerCoord.zero();
300 
302 }
303 
304 
306 {
311 }
312 } // end namespace oofem
void setField(int item, InputFieldType id)
int giveNumberOfColumns() const
Returns number of columns of receiver.
Definition: floatmatrix.h:158
REGISTER_BoundaryCondition(BoundaryCondition)
virtual double evaluateAccelerationAtTime(double t)=0
Returns the second time derivative of the function at given time.
Class and object Domain.
Definition: domain.h:115
static void giveGeneralizedStrainComponents(FloatArray genEps, FloatArray &dphidxi1, FloatArray &dphidxi2, FloatArray &dmdxi1, FloatArray &dmdxi2, FloatArray &m, double &dgamdxi1, double &dgamdxi2, double &gam)
Definition: shell7base.C:1604
virtual int giveNumberOfDomainEquations(int di, const UnknownNumberingScheme &num)
Returns number of equations for given domain in active (current time step) time step.
Definition: engngm.C:391
FloatArray centerCoord
Center coordinate .
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
int giveNumberOfDofManagers() const
Returns number of dof managers in domain.
Definition: domain.h:432
FloatArray genEps
Current generalized strain eps = [xbar_{,alpha}, m_{,alpha}, m, gamma_{,alpha}, gamma ] = [ (3 3) (3 ...
virtual FloatArray & giveCenterCoordinate()
Returns the center coordinate.
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
virtual FloatArray * giveCoordinates()
Definition: dofmanager.h:382
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
EngngModel * giveEngngModel()
Returns engineering model to which receiver is associated.
Definition: domain.C:433
virtual void giveInputRecord(DynamicInputRecord &input)
Setups the input record string of receiver.
Abstract base class for all finite elements.
Definition: element.h:145
int giveNumber()
Returns domain number.
Definition: domain.h:266
virtual double evalNXIntegral(int boundary, const FEICellGeometry &cellgeo)
Computes the integral .
Definition: feinterpol.h:420
int giveNumberOfSpatialDimensions()
Returns number of spatial dimensions.
Definition: domain.C:1067
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
virtual int __givePrescribedEquationNumber()=0
Returns prescribed equation number of receiver.
virtual FEInterpolation * giveInterpolation() const
Definition: element.h:629
void beDifferenceOf(const FloatArray &a, const FloatArray &b)
Sets receiver to be a - b.
Definition: floatarray.C:341
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Class representing a general abstraction for finite element interpolation class.
Definition: feinterpol.h:132
Element * giveElement(int n)
Service for accessing particular domain fe element.
Definition: domain.C:160
void evalCovarBaseVectorsAt(FloatMatrix &gcov, FloatArray &genEps, double zeta)
Constructs a coefficient matrix for all prescribed unknowns.
#define OOFEM_ERROR(...)
Definition: error.h:61
DofIDItem
Type representing particular dof type.
Definition: dofiditem.h:86
Set of elements, boundaries, edges and/or nodes.
Definition: set.h:66
DofIDItem giveDofID() const
Returns DofID value of receiver, which determines type of of unknown connected to receiver (e...
Definition: dof.h:276
void evalInitialCovarBaseVectorsAt(FloatMatrix &Gcov, FloatArray &genEps, double zeta)
Wrapper around element definition to provide FEICellGeometry interface.
Definition: feinterpol.h:95
Set * giveSet(int n)
Service for accessing particular domain set.
Definition: domain.C:363
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
IntArray dofs
Dofs that b.c. is applied to (relevant for Dirichlet type b.c.s).
double at(int i, int j) const
Coefficient access function.
Definition: floatmatrix.h:176
void evaluateHigherOrderContribution(FloatArray &answer, double zeta, FloatArray &dx)
#define _IFT_PrescribedGenStrainShell7_generalizedstrain
virtual double evaluateVelocityAtTime(double t)=0
Returns the first time derivative of the function at given time.
Class representing vector of real numbers.
Definition: floatarray.h:82
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
void updateCoefficientMatrix(FloatMatrix &C)
Constructs a coefficient matrix for all prescribed unknowns.
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
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
The representation of EngngModel default prescribed unknown numbering.
Dof * giveDofWithID(int dofID) const
Returns DOF with given dofID; issues error if not present.
Definition: dofmanager.C:119
FloatMatrix gradient
Prescribed gradient .
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
void setColumn(const FloatArray &src, int c)
Sets the values of the matrix in specified column.
Definition: floatmatrix.C:648
Class representing the a dynamic Input Record.
virtual FloatArray * giveCoordinates()
Definition: node.h:114
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
#define _IFT_PrescribedGenStrainShell7_centercoords
static void giveDualBase(FloatMatrix &base1, FloatMatrix &base2)
Definition: shell7base.C:283
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
DofManager * giveDofManager() const
Definition: dof.h:123
int giveSize() const
Definition: intarray.h:203
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.
Class implementing node in finite element mesh.
Definition: node.h:87
FloatArray initialGenEps
Initial generalized strain eps = [Xbar_{,alpha}, M_{,alpha}, M ] = [ (3 3) (3 3) 3] - 15 components...
#define IR_GIVE_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:69
Abstract class Dof represents Degree Of Freedom in finite element mesh.
Definition: dof.h:93
#define _IFT_PrescribedGenStrainShell7_initialgeneralizedstrain
virtual double evaluateAtTime(double t)
Returns the value of the function at given time.
Definition: function.C:76
virtual void giveInputRecord(DynamicInputRecord &input)
Setups the input record string of receiver.
virtual double give(Dof *dof, ValueModeType mode, double time)
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