MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
shelleas.h
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/struct/shelleas.h,v 1.13 2017/01/12 14:46:44 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 2010-2017
7  *
8  * Marco Morandini <morandini@aero.polimi.it>
9  * Riccardo Vescovini <vescovini@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 /*
33  * Inspired by
34  * Wojciech Witkowski
35  * "4-Node combined shell element with semi-EAS-ANS strain interpolations
36  * in 6-parameter shell theories with drilling degrees of freedom"
37  * Comput Mech (2009) 43:307­319 DOI 10.1007/s00466-008-0307-x
38  */
39 
40 #ifndef SHELLEAS_H
41 #define SHELLEAS_H
42 
43 #include "myassert.h"
44 #include "except.h"
45 
46 #include "strnode.h"
47 #include "elem.h"
48 
49 #include "constltp.h"
50 
51 /* da spostare in .cc */
52 #include "Rot.hh"
53 #include "joint.h"
54 
55 #include "shell.h"
56 
57 // Shell4EAS - begin
58 
59 class Shell4EAS
60 : virtual public Elem,
61 public Shell
62 {
63 #if 0
64 protected:
65  static const unsigned int iNumPrivData =
66  +3 // 0 ( 1 -> 3) - strain
67  +3 // 3 ( 4 -> 6) - curvature
68  +3 // 6 ( 7 -> 9) - force
69  +3 // 9 (10 -> 12) - moment
70  +3 // 12 (13 -> 15) - position
71  +3 // 15 (16 -> 18) - orientation vector
72  +3 // 18 (19 -> 21) - angular velocity
73  +3 // 21 (22 -> 24) - strain rate
74  +3 // 24 (25 -> 27) - curvature rate
75  ;
76 
77  static unsigned int iGetPrivDataIdx_int(const char *s,
78  ConstLawType::Type type);
79 #endif
80 private:
81 public:
82  // numbered according to
83  //
84  // ^
85  // 4 o-----+-----o 3
86  // | 1_2 | 2_2 |
87  // --+-----+-----+->
88  // | 1_1 | 2_1 |
89  // 1 o-----+-----o 2
90  //
92  IP_1_1 = 0,
93  IP_1_2 = 1,
94  IP_1_3 = 2,
95  IP_2_1 = 3,
96  IP_2_2 = 4,
97  IP_2_3 = 5,
98  IP_3_1 = 6,
99  IP_3_2 = 7,
100  IP_3_3 = 8,
101 
102  NUMIP = 4
103  };
104 
105  static doublereal xi_i[NUMIP][2];
107 
108  // numbered according to the side they are defined on
110  SSEP_1 = 0,
111  SSEP_2 = 1,
112  SSEP_3 = 2,
113  SSEP_4 = 3,
114  };
115 
116  enum NodeName {
117  NODE1 = 0,
118  NODE2 = 1,
119  NODE3 = 2,
120  NODE4 = 3,
121 
123  };
124 
126 
127  static doublereal xi_0[2];
128 
130  STRAIN = 0,
131  CURVAT = 1,
132 
134  };
135 
136 protected:
137  // Pointers to nodes
139 
140 #if 0
141  // Node offsets - TODO: offsets
142  const Vec3 f[NUMNODES];
143  Vec3 fRef[NUMNODES];
144 #endif
145 
146  // nodal positions (0: initial; otherwise current)
149  // current nodal orientation
152  // Euler vector of Ra
154 
155  // Average orientation matrix
157  // Average orientation matrix
158  // .. in reference configuration
160  // .. in current configuration
162 
163  // Orientation matrix
164  // .. in reference configuration
167  // .. in current configuration
170 
174 
175 // // Orientation matrix of the shear strain evaluation points
176 // Mat3x3 R[NUMSSEP];
177 // Mat3x3 RRef[NUMSSEP];
178 // Mat3x3 RPrev[NUMSSEP];
179 
180  // rotation tensors
182 
183  // Orientation tensor derivative axial vector
186 // Mat3x3 T_1_i[NUMIP];
187 // Mat3x3 T_2_i[NUMIP];
188 
189  // linear deformation vectors
190  // .. in reference configuration
193  // .. in current configuration
196 
197  // angular deformation vectors
198  // .. in reference configuration
201  // .. in current configuration
204 
205 
206 #if 0
207  // Angular velocity of the sections - TODO: viscoelastic
208  Vec3 Omega[NUMSEZ];
209  Vec3 OmegaRef[NUMSEZ];
210 #endif
211 
212  // Temporary data - TODO
213 #if 0
214  Vec6 Az[NUMSEZ];
215  Vec6 AzRef[NUMSEZ];
216  Vec6 AzLoc[NUMSEZ];
217  Vec6 DefLoc[NUMSEZ];
218  Vec6 DefLocRef[NUMSEZ];
219  Vec6 DefLocPrev[NUMSEZ];
220 #endif
221 
222 protected:
227 
229 
231 
234 
238 
239 #ifdef USE_CL_IN_SHELL
240  // Constitutive law handlers at integration points
242 #else // ! USE_CL_IN_SHELL
246 #endif // ! USE_CL_IN_SHELL
247 
248  // Reference constitutive law tangent matrices
250 
251  //stress
253 
254  // Is first residual
255  bool bFirstRes;
256 
257  // for derived elements that add external contributions
258  // to internal forces
259  virtual void
260  AddInternalForces(Vec6& /* AzLoc */ , unsigned int /* iSez */ ) {
261  NO_OP;
262  };
263 
264 private:
267  void InterpolateOrientation();
268 // void ComputeIPSEPRotations();
269  void ComputeIPCurvature();
270 
271 public:
272  Shell4EAS(unsigned int uL,
273  const DofOwner* pDO,
274  const StructNode* pN1, const StructNode* pN2,
275  const StructNode* pN3, const StructNode* pN4,
276 #if 0 // TODO: offset
277  const Vec3& f1, const Vec3& f2,
278  const Vec3& f3, const Vec3& f4,
279 #endif
280  const Mat3x3& R1, const Mat3x3& R2,
281  const Mat3x3& R3, const Mat3x3& R4,
282 #ifdef USE_CL_IN_SHELL
283  const ConstitutiveLaw<vh, fmh>** pDTmp,
284 #else // ! USE_CL_IN_SHELL
285  const fmh& pDTmp,
286  const vh& PreStress,
287 #endif // ! USE_CL_IN_SHELL
288  flag fOut);
289 
290  virtual ~Shell4EAS(void);
291 
292  // Shell type
293  virtual Shell::Type GetShellType(void) const {
294  return Shell::ELASTIC;
295  };
296 
297  virtual unsigned int iGetNumDof(void) const {
298  //return 14;
299  return 13;
300  };
301 
302  virtual DofOrder::Order GetDofType(unsigned int i) const {
303  ASSERT(i >= 0 && i < iGetNumDof());
304  return DofOrder::ALGEBRAIC;
305  };
306 
307  virtual DofOrder::Order GetEqType(unsigned int i) const {
308  ASSERT(i >= 0 && i < iGetNumDof());
309  return DofOrder::ALGEBRAIC;
310  };
311 
312  // Contribution to restart file
313  virtual std::ostream& Restart(std::ostream& out) const;
314 
315 #if 0
316  virtual void
317  AfterConvergence(const VectorHandler& X, const VectorHandler& XP);
318 
319  // Inverse dynamics
320  virtual void
321  AfterConvergence(const VectorHandler& X, const VectorHandler& XP,
322  const VectorHandler& XPP);
323 #endif
324 
325  // Workspace size
326  virtual void
327  WorkSpaceDim(integer* piNumRows, integer* piNumCols) const {
328  *piNumRows = 6*4 + iGetNumDof();
329  *piNumCols = 6*4 + iGetNumDof();
330  };
331 
332  // Initial settings
333  void
334  SetValue(DataManager *pDM,
335  VectorHandler& /* X */ , VectorHandler& /* XP */ ,
336  SimulationEntity::Hints *ph = 0);
337 
338 #if 0
339  // Prepares reference parameters after prediction
340  virtual void
341  AfterPredict(VectorHandler& /* X */ , VectorHandler& /* XP */ );
342 #endif
343 
344  // Residual assembly
345  virtual SubVectorHandler& AssRes(SubVectorHandler& WorkVec,
346  doublereal dCoef,
347  const VectorHandler& XCurr,
348  const VectorHandler& XPrimeCurr);
349 
350 #if 0
351  // Inverse dynamics
352  virtual SubVectorHandler&
353  AssRes(SubVectorHandler& WorkVec,
354  const VectorHandler& XCurr ,
355  const VectorHandler& XPrimeCurr ,
356  const VectorHandler& XPrimePrimeCurr ,
358 #endif
359 
360  // Jacobian matrix assembly
361  virtual VariableSubMatrixHandler&
363  doublereal dCoef,
364  const VectorHandler& XCurr,
365  const VectorHandler& XPrimeCurr);
366 
367 #if 0
368  // Matrix assembly for eigenvalues
369  void
371  VariableSubMatrixHandler& WorkMatB,
372  const VectorHandler& XCurr,
373  const VectorHandler& XPrimeCurr);
374 #endif
375 
376  virtual void Output(OutputHandler& OH) const;
377 
378  virtual unsigned int iGetInitialNumDof(void) const {
379  return 0;
380  };
381 
382  virtual void
383  InitialWorkSpaceDim(integer* piNumRows, integer* piNumCols) const {
384  *piNumRows = 6*4;
385  *piNumCols = 6*4;
386  };
387 
388  virtual void SetInitialValue(VectorHandler& /* X */ ) {
389  NO_OP;
390  };
391 
392  virtual VariableSubMatrixHandler&
394  const VectorHandler& XCurr);
395 
396  virtual SubVectorHandler&
397  InitialAssRes(SubVectorHandler& WorkVec, const VectorHandler& XCurr);
398 
399 #if 0
400  // Access to private data
401  virtual unsigned int iGetNumPrivData(void) const;
402  virtual unsigned int iGetPrivDataIdx(const char *s) const;
403  virtual doublereal dGetPrivData(unsigned int i) const;
404 #endif
405 
406  // Access to nodes
407  virtual const StructNode* pGetNode(unsigned int i) const;
408 
409  /******** PER IL SOLUTORE PARALLELO *********/
410  /* Fornisce il tipo e la label dei nodi che sono connessi all'elemento
411  * utile per l'assemblaggio della matrice di connessione fra i dofs */
412  virtual void
413  GetConnectedNodes(std::vector<const Node *>& connectedNodes) const {
414  connectedNodes.resize(NUMNODES);
415  for (int i = 0; i < NUMNODES; i++) {
416  connectedNodes[i] = pNode[i];
417  }
418  };
419  /**************************************************/
420 };
421 
422 // Shell4EAS - end
423 
424 #endif // SHELLEAS_H
425 
virtual DofOrder::Order GetDofType(unsigned int i) const
Definition: shelleas.h:302
Mat3x3 T_overline
Definition: shelleas.h:161
void ComputeInitialNodeAndIptOrientation()
Definition: shelleas.cc:123
Shell4EAS(unsigned int uL, const DofOwner *pDO, const StructNode *pN1, const StructNode *pN2, const StructNode *pN3, const StructNode *pN4, const Mat3x3 &R1, const Mat3x3 &R2, const Mat3x3 &R3, const Mat3x3 &R4, const fmh &pDTmp, const vh &PreStress, flag fOut)
Definition: shelleas.cc:215
vfmh P_i
Definition: shelleas.h:230
virtual std::ostream & Restart(std::ostream &out) const
Definition: shelleas.cc:1205
void SetValue(DataManager *pDM, VectorHandler &, VectorHandler &, SimulationEntity::Hints *ph=0)
Definition: shelleas.cc:1212
fmh S_alpha_beta_0
Definition: shelleas.h:223
vh epsilon_hat
Definition: shelleas.h:236
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: shelleas.cc:545
long int flag
Definition: mbdyn.h:43
Vec3 y_i_1[NUMIP]
Definition: shelleas.h:232
bool bFirstRes
Definition: shelleas.h:255
Definition: matvec3.h:98
Mat3x3 T_0_i[NUMIP]
Definition: shelleas.h:166
Vec3 k_tilde_1_i[NUMIP]
Definition: shelleas.h:202
static doublereal xi_i[NUMIP][2]
Definition: shelleas.h:105
Vec3 y_i_2[NUMIP]
Definition: shelleas.h:233
bool bPreStress
Definition: shelleas.h:244
vfmh L_alpha_beta_i
Definition: shelleas.h:226
virtual void Output(OutputHandler &OH) const
Definition: shelleas.cc:1251
Mat3x3 Phi_Delta_i[NUMIP][NUMNODES]
Definition: shelleas.h:171
std::vector< vh > vvh
Definition: shell.h:81
vvh stress_i
Definition: shelleas.h:252
virtual unsigned int iGetPrivDataIdx(const char *s) const
Definition: simentity.cc:142
static const unsigned int iNumPrivData
Definition: beam.cc:249
vvh FRef
Definition: shelleas.h:243
Mat3x3 Kappa_delta_i_1[NUMIP][NUMNODES]
Definition: shelleas.h:172
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: shelleas.h:383
Vec3 k_tilde_2_i[NUMIP]
Definition: shelleas.h:203
Mat3x3 T0_overline
Definition: shelleas.h:159
Vec3 phi_tilde_i[NUMIP]
Definition: shelleas.h:156
#define NO_OP
Definition: myassert.h:74
void ComputeIPCurvature()
Definition: shelleas.cc:191
std::vector< Hint * > Hints
Definition: simentity.h:89
static doublereal xi_0[2]
Definition: shelleas.h:127
Vec3 eps_tilde_1_0_i[NUMIP]
Definition: shelleas.h:191
static doublereal xi_n[NUMNODES][2]
Definition: shelleas.h:125
virtual ~Shell4EAS(void)
Definition: shelleas.cc:410
doublereal alpha_i[NUMIP]
Definition: shelleas.h:225
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const
Definition: shelleas.h:413
ShearStrainEvaluationPoint
Definition: shelleas.h:109
virtual DofOrder::Order GetEqType(unsigned int i) const
Definition: shelleas.h:307
virtual void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: shelleas.h:327
Mat3x3 T_i[NUMIP]
Definition: shelleas.h:169
virtual unsigned int iGetNumDof(void) const
Definition: shelleas.h:297
vfmh DRef
Definition: shelleas.h:249
Vec3 phi_tilde_n[NUMNODES]
Definition: shelleas.h:153
Definition: matvec6.h:37
Vec3 k_tilde_1_0_i[NUMIP]
Definition: shelleas.h:199
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: shelleas.cc:1228
virtual void AddInternalForces(Vec6 &, unsigned int)
Definition: shelleas.h:260
Mat3x3 Kappa_delta_i_2[NUMIP][NUMNODES]
Definition: shelleas.h:173
Vec3 xa[NUMNODES]
Definition: shelleas.h:148
IntegrationPoint
Definition: shelleas.h:91
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: shelleas.cc:422
static doublereal w_i[NUMIP]
Definition: shelleas.h:106
doublereal alpha_0
Definition: shelleas.h:224
Vec3 k_1_i[NUMIP]
Definition: shelleas.h:184
#define ASSERT(expression)
Definition: colamd.c:977
vh epsilon
Definition: shelleas.h:237
vh PreStress
Definition: shelleas.h:245
virtual doublereal dGetPrivData(unsigned int i) const
Definition: simentity.cc:149
virtual unsigned int iGetNumPrivData(void) const
Definition: simentity.cc:136
Mat3x3 iTa_i[NUMIP]
Definition: shelleas.h:151
Mat3x3 iTa[NUMNODES]
Definition: shelleas.h:150
Vec3 eps_tilde_2_0_i[NUMIP]
Definition: shelleas.h:192
vfmh B_overline_i
Definition: shelleas.h:228
Definition: elem.h:75
Vec3 xa_0[NUMNODES]
Definition: shelleas.h:147
Type
Definition: shell.h:71
virtual void AfterConvergence(const VectorHandler &X, const VectorHandler &XP)
Definition: simentity.cc:120
virtual void AssMats(VariableSubMatrixHandler &WorkMatA, VariableSubMatrixHandler &WorkMatB, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: elem.cc:55
Mat3x3 Q_i[NUMIP]
Definition: shelleas.h:181
FullMatrixHandler fmh
Definition: shell.h:82
virtual const StructNode * pGetNode(unsigned int i) const
Definition: shelleas.cc:1236
std::vector< fmh > vfmh
Definition: shell.h:83
void UpdateNodalAndAveragePosAndOrientation()
Definition: shelleas.cc:101
virtual void AfterPredict(VectorHandler &X, VectorHandler &XP)
Definition: simentity.cc:91
Mat3x3 T_0_0
Definition: shelleas.h:165
double doublereal
Definition: colamd.c:52
long int integer
Definition: colamd.c:51
virtual unsigned int iGetInitialNumDof(void) const
Definition: shelleas.h:378
Vec3 eps_tilde_1_i[NUMIP]
Definition: shelleas.h:194
const StructNode * pNode[NUMNODES]
Definition: shelleas.h:138
Mat3x3 T_0
Definition: shelleas.h:168
Vec3 k_2_i[NUMIP]
Definition: shelleas.h:185
virtual Shell::Type GetShellType(void) const
Definition: shelleas.h:293
Vec3 eps_tilde_2_i[NUMIP]
Definition: shelleas.h:195
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: shelleas.cc:1220
void InterpolateOrientation()
Definition: shelleas.cc:157
Definition: shell.h:63
Vec3 k_tilde_2_0_i[NUMIP]
Definition: shelleas.h:200
virtual void SetInitialValue(VectorHandler &)
Definition: shelleas.h:388
MyVectorHandler vh
Definition: shell.h:80