MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
constltp_impl.cc
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/base/constltp_impl.cc,v 1.55 2017/01/12 14:46:09 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@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 /* Legami costitutivi */
33 
34 #include "mbconfig.h" /* This goes first in every *.c,*.cc file */
35 
36 #include "myassert.h"
37 #include "mynewmem.h"
38 #include "drive_.h"
39 #include "constltp_impl.h"
40 
41 #include "symcltp.h"
42 #ifdef USE_GINAC
43 #include "ginaccltp.h"
44 #endif // USE_GINAC
45 #include "shockabsorber.h"
46 #include "constltp_ann.h"
47 
48 /* constitutive laws sponsored by Hutchinson CdR */
49 #include "constltp_nlp.h"
50 #include "constltp_nlsf.h"
51 
52 /* used by invariant constitutive law... */
53 #include "vehj.h"
54 
55 /* ... */
56 #include "tdclw.h"
57 #include "constltp_axw.h"
58 
59 /* constitutive law containers */
60 typedef std::map<std::string, ConstitutiveLawRead<doublereal, doublereal> *, ltstrcase> CL1DFuncMapType;
61 typedef std::map<std::string, ConstitutiveLawRead<Vec3, Mat3x3> *, ltstrcase> CL3DFuncMapType;
62 typedef std::map<std::string, ConstitutiveLawRead<Vec6, Mat6x6> *, ltstrcase> CL6DFuncMapType;
63 
67 
68 /* constitutive law parsing checkers */
70  bool IsWord(const std::string& s) const {
71  return CL1DFuncMap.find(std::string(s)) != CL1DFuncMap.end();
72  };
73 };
74 
76  bool IsWord(const std::string& s) const {
77  return CL3DFuncMap.find(std::string(s)) != CL3DFuncMap.end();
78  };
79 };
80 
82  bool IsWord(const std::string& s) const {
83  return CL6DFuncMap.find(std::string(s)) != CL6DFuncMap.end();
84  };
85 };
86 
90 
91 /* constitutive law registration functions: call to register one */
92 bool
94 {
95  pedantic_cout("registering constitutive law 1D \"" << name << "\""
96  << std::endl );
97  return CL1DFuncMap.insert(CL1DFuncMapType::value_type(name, rf)).second;
98 }
99 
100 bool
102 {
103  pedantic_cout("registering constitutive law 3D \"" << name << "\""
104  << std::endl );
105  return CL3DFuncMap.insert(CL3DFuncMapType::value_type(name, rf)).second;
106 }
107 
108 bool
110 {
111  pedantic_cout("registering constitutive law 6D \"" << name << "\""
112  << std::endl );
113  return CL6DFuncMap.insert(CL6DFuncMapType::value_type(name, rf)).second;
114 }
115 
116 /* function that reads a constitutive law */
119 {
120  const char *s, *sOrig = HP.IsWord(CL1DWordSet);
121  if (sOrig == 0) {
122  /* default to linear elastic? */
123  s = "linear" "elastic";
124  sOrig = "";
125 
126  } else {
127  s = sOrig;
128  }
129 
130  CL1DFuncMapType::iterator func = CL1DFuncMap.find(std::string(s));
131  if (func == CL1DFuncMap.end()) {
132  silent_cerr("unknown constitutive law 1D type "
133  "\"" << sOrig << "\" "
134  "at line " << HP.GetLineData() << std::endl);
136  }
137 
138  return func->second->Read(pDM, HP, CLType);
139 }
140 
143 {
144  const char *s, *sOrig = HP.IsWord(CL3DWordSet);
145  if (sOrig == 0) {
146 #if 0
147  s = "linear" "elastic";
148  sOrig = "";
149 #else
150  silent_cerr("unknown constitutive law 3D type "
151  "at line " << HP.GetLineData() << std::endl);
153 #endif
154 
155  } else {
156  s = sOrig;
157  }
158 
159  CL3DFuncMapType::iterator func = CL3DFuncMap.find(std::string(s));
160  if (func == CL3DFuncMap.end()) {
161  silent_cerr("unknown constitutive law 3D type "
162  "\"" << sOrig << "\" "
163  "at line " << HP.GetLineData() << std::endl);
165  }
166 
167  return func->second->Read(pDM, HP, CLType);
168 }
169 
172 {
173  const char *s, *sOrig = HP.IsWord(CL6DWordSet);
174  if (sOrig == 0) {
175 #if 0
176  s = "linear" "elastic";
177  sOrig = "";
178 #else
179  silent_cerr("unknown constitutive law 6D type "
180  "at line " << HP.GetLineData() << std::endl);
182 #endif
183 
184  } else {
185  s = sOrig;
186  }
187 
188  CL6DFuncMapType::iterator func = CL6DFuncMap.find(std::string(s));
189  if (func == CL6DFuncMap.end()) {
190  silent_cerr("unknown constitutive law 6D type "
191  "\"" << sOrig << "\" "
192  "at line " << HP.GetLineData() << std::endl);
194  }
195 
196  return func->second->Read(pDM, HP, CLType);
197 }
198 
199 /* specific functional object(s) */
200 struct CLArray1DR : public ConstitutiveLawRead<doublereal, doublereal> {
202  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
204 
205  unsigned n = HP.GetInt();
206  if (n <= 0) {
207  silent_cerr("constitutive law array 1D: invalid constitutive law number " << n
208  << " at line " << HP.GetLineData() << std::endl);
210  }
211 
212  if (n == 1) {
213  return ReadCL1D(pDM, HP, CLType);
214  }
215 
216  std::vector<ConstitutiveLaw<doublereal, doublereal> *> clv(n);
217  for (unsigned i = 0; i < n; i++) {
218  clv[i] = ReadCL1D(pDM, HP, CLType);
219  }
220 
222  SAFENEWWITHCONSTRUCTOR(pCL, L, L(clv));
223 
224  CLType = pCL->GetConstLawType();
225 
226  return pCL;
227  };
228 };
229 
230 struct CLArray3DR : public ConstitutiveLawRead<Vec3, Mat3x3> {
232  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
234 
235  unsigned n = HP.GetInt();
236  if (n <= 0) {
237  silent_cerr("constitutive law array 1D: invalid constitutive law number " << n
238  << " at line " << HP.GetLineData() << std::endl);
240  }
241 
242  if (n == 1) {
243  return ReadCL3D(pDM, HP, CLType);
244  }
245 
246  std::vector<ConstitutiveLaw<Vec3, Mat3x3> *> clv(n);
247  for (unsigned i = 0; i < n; i++) {
248  clv[i] = ReadCL3D(pDM, HP, CLType);
249  }
250 
252  SAFENEWWITHCONSTRUCTOR(pCL, L, L(clv));
253 
254  CLType = pCL->GetConstLawType();
255 
256  return pCL;
257  };
258 };
259 
260 struct CLArray6DR : public ConstitutiveLawRead<Vec6, Mat6x6> {
262  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
264 
265  unsigned n = HP.GetInt();
266  if (n <= 0) {
267  silent_cerr("constitutive law array 1D: invalid constitutive law number " << n
268  << " at line " << HP.GetLineData() << std::endl);
270  }
271 
272  if (n == 1) {
273  return ReadCL6D(pDM, HP, CLType);
274  }
275 
276  std::vector<ConstitutiveLaw<Vec6, Mat6x6> *> clv(n);
277  for (unsigned i = 0; i < n; i++) {
278  clv[i] = ReadCL6D(pDM, HP, CLType);
279  }
280 
282  SAFENEWWITHCONSTRUCTOR(pCL, L, L(clv));
283 
284  CLType = pCL->GetConstLawType();
285 
286  return pCL;
287  };
288 };
289 
290 template <class T, class Tder>
291 struct LinearElasticCLR : public ConstitutiveLawRead<T, Tder> {
292  virtual ConstitutiveLaw<T, Tder> *
293  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
294  ConstitutiveLaw<T, Tder>* pCL = 0;
295 
296  CLType = ConstLawType::ELASTIC;
297 
298  doublereal dS = HP.GetReal();
299  DEBUGCOUT("Linear Elastic Isotropic Constitutive Law, stiffness = "
300  << dS << std::endl);
301 
302  if (dS <= 0.) {
303  silent_cerr("warning, null or negative stiffness at line "
304  << HP.GetLineData() << std::endl);
305  }
306 
307  /* Prestress and prestrain */
308  T PreStress(mb_zero<T>());
309  GetPreStress(HP, PreStress);
310  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
311 
313  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dS));
314 
315  return pCL;
316  };
317 };
318 
319 template <class T, class Tder>
321  virtual ConstitutiveLaw<T, Tder> *
322  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
323  ConstitutiveLaw<T, Tder>* pCL = 0;
324 
325  CLType = ConstLawType::ELASTIC;
326 
327  DEBUGCOUT("Linear Elastic Generic Constitutive Law" << std::endl);
328  Tder S(mb_zero<Tder>());
329  S = HP.Get(S);
330 
331  /* Prestress and prestrain */
332  T PreStress(mb_zero<T>());
333  GetPreStress(HP, PreStress);
334  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
335 
337  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S));
338 
339  return pCL;
340  };
341 };
342 
343 template <class T, class Tder>
345  virtual ConstitutiveLaw<T, Tder> *
346  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
347  ConstitutiveLaw<T, Tder>* pCL = 0;
348 
349  CLType = ConstLawType::ELASTIC;
350 
351  DEBUGCOUT("Linear Elastic Generic Constitutive Law with Axial-Torsion Coupling" << std::endl);
352  Tder S(mb_zero<Tder>());
353  S = HP.Get(S);
354 
355  /* coefficiente di accoppiamento */
356  doublereal dCoupl = HP.GetReal();
357  DEBUGCOUT("coupling coefficient: " << dCoupl << std::endl);
358 
359  /* Prestress and prestrain */
360  T PreStress(mb_zero<T>());
361  GetPreStress(HP, PreStress);
362  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
363 
365  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S, dCoupl));
366 
367  return pCL;
368  };
369 };
370 
371 template <class T, class Tder>
373  virtual ConstitutiveLaw<T, Tder> *
374  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
375  ConstitutiveLaw<T, Tder>* pCL = 0;
376 
378 
379  Tder S(mb_zero<Tder>());
380  S = HP.Get(S);
381 
382  Tder SP(mb_zero<Tder>());
383  if (HP.IsKeyWord("proportional")) {
384  doublereal k = HP.GetReal();
385  SP = S*k;
386  } else {
387  SP = HP.Get(SP);
388  }
389 
390  /* coefficiente di accoppiamento */
391  doublereal dCoupl = HP.GetReal();
392  DEBUGCOUT("coupling coefficient: " << dCoupl << std::endl);
393 
394  /* Prestress and prestrain */
395  T PreStress(mb_zero<T>());
396  GetPreStress(HP, PreStress);
397  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
398 
400  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S, SP, dCoupl));
401 
402  return pCL;
403  };
404 };
405 
406 struct InverseSquareElasticCLR : public ConstitutiveLawRead<doublereal, doublereal> {
408  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
410 
411  CLType = ConstLawType::ELASTIC;
412 
413  DEBUGCOUT("Inverse Square Elastic Constitutive Law" << std::endl);
414  doublereal dA = HP.GetReal();
415  if (dA <= 0.) {
416  silent_cerr("warning, null or negative stiffness at line "
417  << HP.GetLineData() << std::endl);
418  }
419 
420  doublereal dL0 = HP.GetReal();
421  if (dL0 <= 0.) {
422  silent_cerr("warning, null or negative reference length at line "
423  << HP.GetLineData() << std::endl);
424  }
425 
426  /* Prestress and prestrain */
427  doublereal PreStress(mb_zero<doublereal>());
428  GetPreStress(HP, PreStress);
429  TplDriveCaller<doublereal>* pTplDC = GetPreStrain<doublereal>(pDM, HP);
430 
431 
434  InverseSquareConstitutiveLaw(pTplDC, PreStress, dA, dL0));
435 
436  return pCL;
437  };
438 };
439 
440 template <class T, class Tder>
441 struct LogElasticCLR : public ConstitutiveLawRead<T, Tder> {
442  virtual ConstitutiveLaw<T, Tder> *
443  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
444  ConstitutiveLaw<T, Tder>* pCL = 0;
445 
446  CLType = ConstLawType::ELASTIC;
447 
448  DEBUGCOUT("Logarithmic Elastic Constitutive Law" << std::endl);
449  doublereal dS = HP.GetReal();
450  if (dS <= 0.) {
451  silent_cerr("warning, null or negative stiffness at line "
452  << HP.GetLineData() << std::endl);
453  }
454 
455  /* Prestress and prestrain */
456  T PreStress(mb_zero<T>());
457  GetPreStress(HP, PreStress);
458  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
459 
460  typedef LogConstitutiveLaw<T, Tder> L;
461  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dS));
462 
463  return pCL;
464  };
465 };
466 
467 template <class T, class Tder>
468 struct DoubleLinearElasticCLR : public ConstitutiveLawRead<T, Tder> {
469  virtual ConstitutiveLaw<T, Tder> *
470  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
471  ConstitutiveLaw<T, Tder>* pCL = 0;
472 
473  CLType = ConstLawType::ELASTIC;
474 
475  doublereal dS = HP.GetReal();
476  DEBUGCOUT("stiffness = " << dS << std::endl);
477 
478  if (dS <= 0.) {
479  silent_cerr("warning, null or negative stiffness at line "
480  << HP.GetLineData() << std::endl);
481  }
482 
483  doublereal dUpp = HP.GetReal();
484  if (dUpp <= 0.) {
485  silent_cerr("warning, null or negative upper limit strain at line "
486  << HP.GetLineData() << std::endl);
487  }
488 
489  doublereal dLow = HP.GetReal();
490  if (dLow >= 0.) {
491  silent_cerr("warning, null or positive lower limit strain at line "
492  << HP.GetLineData() << std::endl);
493  }
494 
495  doublereal dSecondS = HP.GetReal();
496  if (dSecondS <= 0.) {
497  silent_cerr("warning, null or negative second stiffness at line "
498  << HP.GetLineData() << std::endl);
499  }
500 
501  /* Prestress and prestrain */
502  T PreStress(mb_zero<T>());
503  GetPreStress(HP, PreStress);
504  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
505 
508  L,
509  L(pTplDC, PreStress, dS, dUpp, dLow, dSecondS));
510 
511  return pCL;
512  };
513 };
514 
515 template <class T, class Tder>
516 struct IsotropicHardeningCLR : public ConstitutiveLawRead<T, Tder> {
517  virtual ConstitutiveLaw<T, Tder> *
518  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
519  ConstitutiveLaw<T, Tder>* pCL = 0;
520 
521  CLType = ConstLawType::ELASTIC;
522 
523  doublereal dS = HP.GetReal();
524  DEBUGCOUT("Stiffness = " << dS << std::endl);
525 
526  if (dS <= 0.) {
527  silent_cerr("warning, null or negative stiffness at line "
528  << HP.GetLineData() << std::endl);
529  }
530 
531  doublereal dE = HP.GetReal();
532  DEBUGCOUT("Reference strain = " << dE << std::endl);
533 
534  if (dE <= 0.) {
535  silent_cerr("error, null or negative reference strain at line "
536  << HP.GetLineData() << std::endl);
538  }
539 
540  doublereal dS0 = 0.;
541  if (HP.IsKeyWord("linear" "stiffness")) {
542  dS0 = HP.GetReal();
543  }
544 
545  /* Prestress and prestrain */
546  T PreStress(mb_zero<T>());
547  GetPreStress(HP, PreStress);
548  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
549 
551  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dS, dS0, dE));
552 
553  return pCL;
554  };
555 };
556 
557 template <class T, class Tder>
558 struct ContactElasticCLR : public ConstitutiveLawRead<T, Tder> {
559  virtual ConstitutiveLaw<T, Tder> *
560  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
561  ConstitutiveLaw<T, Tder>* pCL = 0;
562 
563  CLType = ConstLawType::ELASTIC;
564 
565  doublereal dK = HP.GetReal();
566  DEBUGCOUT("Stiffness = " << dK << std::endl);
567 
568  if (dK <= 0.) {
569  silent_cerr("warning, null or negative stiffness at line "
570  << HP.GetLineData() << std::endl);
571  }
572 
573  doublereal dGamma = HP.GetReal();
574  DEBUGCOUT("Exponent = " << dGamma << std::endl);
575 
576  if (dGamma < 1.) {
577  silent_cerr("error, exponent < 1. at line "
578  << HP.GetLineData() << std::endl);
580  }
581 
582  /* Prestress and prestrain */
583  T PreStress(mb_zero<T>());
584  GetPreStress(HP, PreStress);
585  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
586 
588  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dK, dGamma));
589 
590  return pCL;
591  };
592 };
593 
594 template <class T, class Tder>
595 struct SymbolicCLR : public ConstitutiveLawRead<T, Tder> {
596  virtual ConstitutiveLaw<T, Tder> *
597  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
598  ConstitutiveLaw<T, Tder>* pCL = 0;
599 
600  unsigned dim;
601  if (typeid(T) == typeid(doublereal)) {
602  dim = 1;
603 
604  } else if (typeid(T) == typeid(Vec3)) {
605  dim = 3;
606 
607  } else if (typeid(T) == typeid(Vec6)) {
608  dim = 6;
609 
610  } else {
611  silent_cerr("Invalid dimensionality "
612  "for symbolic constitutive law "
613  "at line " << HP.GetLineData()
614  << std::endl);
616  }
617 
618  std::vector<std::string> epsilon;
619  if (CLType & ConstLawType::ELASTIC) {
620  if (!HP.IsKeyWord("epsilon")) {
621  silent_cerr("keyword \"epsilon\" expected at line " << HP.GetLineData() << std::endl);
623  }
624 
625  epsilon.resize(dim);
626 
627  for (unsigned row = 0; row < dim; row++) {
628  const char *tmp = HP.GetStringWithDelims();
629 
630  if (tmp == 0) {
631  silent_cerr("unable to get \"epsilon\" "
632  "symbol #" << row << " "
633  "at line " << HP.GetLineData() << std::endl);
635  }
636  epsilon[row] = tmp;
637  }
638  }
639 
640  std::vector<std::string> epsilonPrime;
641  if (CLType & ConstLawType::VISCOUS) {
642  if (!HP.IsKeyWord("epsilon" "prime")) {
643  silent_cerr("keyword \"epsilon prime\" expected at line " << HP.GetLineData() << std::endl);
645  }
646 
647  epsilonPrime.resize(dim);
648 
649  for (unsigned row = 0; row < dim; row++) {
650  const char *tmp = HP.GetStringWithDelims();
651 
652  if (tmp == 0) {
653  silent_cerr("unable to get \"epsilonPrime\" "
654  "symbol #" << row << " "
655  "at line " << HP.GetLineData() << std::endl);
657  }
658  epsilonPrime[row] = tmp;
659  }
660  }
661 
662  if (!HP.IsKeyWord("expression")) {
663  silent_cerr("keyword \"expression\" expected at line " << HP.GetLineData() << std::endl);
665  }
666 
667  std::vector<std::string> expression(dim);
668  for (unsigned row = 0; row < dim; row++) {
669  const char *tmp = HP.GetStringWithDelims();
670  if (tmp == 0) {
671  silent_cerr("unable to get \"expression\" "
672  "#" << row << " "
673  "at line " << HP.GetLineData()
674  << std::endl);
676  }
677  expression[row] = tmp;
678  }
679 
680  /* Prestress and prestrain */
681  T PreStress(mb_zero<T>());
682  GetPreStress(HP, PreStress);
683 #ifdef USE_GINAC
684  TplDriveCaller<T>* pTplDC =
685 #endif /* ! USE_GINAC */
686  GetPreStrain<T>(pDM, HP);
687 
688  switch (CLType) {
689  case ConstLawType::ELASTIC: {
690 #ifdef USE_GINAC
692  SAFENEWWITHCONSTRUCTOR(pCL, L,
693  L(pTplDC, PreStress,
694  epsilon,
695  expression));
696 #else /* ! USE_GINAC */
697  silent_cerr("symbolic constitutive law not supported "
698  "at line " << HP.GetLineData() << std::endl);
700 #endif /* ! USE_GINAC */
701  break;
702  }
703 
704  case ConstLawType::VISCOUS: {
705 #ifdef USE_GINAC
707  SAFENEWWITHCONSTRUCTOR(pCL, L,
708  L(PreStress, epsilonPrime, expression));
709 #else /* ! USE_GINAC */
710  silent_cerr("symbolic constitutive law not supported "
711  "at line " << HP.GetLineData() << std::endl);
713 #endif /* ! USE_GINAC */
714  break;
715  }
716 
718 #ifdef USE_GINAC
720  SAFENEWWITHCONSTRUCTOR(pCL, L,
721  L(pTplDC, PreStress,
722  epsilon, epsilonPrime,
723  expression));
724 #else /* ! USE_GINAC */
725  silent_cerr("symbolic constitutive law not supported "
726  "at line " << HP.GetLineData() << std::endl);
728 #endif /* ! USE_GINAC */
729  break;
730  }
731 
732  default:
733  ASSERT(0);
735  }
736 
737  return pCL;
738  };
739 };
740 
741 template <class T, class Tder>
742 struct SymbolicElasticCLR : public SymbolicCLR<T, Tder> {
743  virtual ConstitutiveLaw<T, Tder> *
744  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
745  CLType = ConstLawType::ELASTIC;
746  return SymbolicCLR<T, Tder>::Read(pDM, HP, CLType);
747  };
748 };
749 
750 template <class T, class Tder>
751 struct SymbolicViscousCLR : public SymbolicCLR<T, Tder> {
752  virtual ConstitutiveLaw<T, Tder> *
753  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
754  CLType = ConstLawType::VISCOUS;
755  return SymbolicCLR<T, Tder>::Read(pDM, HP, CLType);
756  };
757 };
758 
759 template <class T, class Tder>
760 struct SymbolicViscoElasticCLR : public SymbolicCLR<T, Tder> {
761  virtual ConstitutiveLaw<T, Tder> *
762  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
764  return SymbolicCLR<T, Tder>::Read(pDM, HP, CLType);
765  };
766 };
767 
768 template <class T, class Tder>
769 struct LinearViscousCLR : public ConstitutiveLawRead<T, Tder> {
770  virtual ConstitutiveLaw<T, Tder> *
771  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
772  ConstitutiveLaw<T, Tder>* pCL = 0;
773 
774  CLType = ConstLawType::VISCOUS;
775 
776  doublereal dSP = HP.GetReal();
777  DEBUGCOUT("stiffness prime = " << dSP << std::endl);
778 
779  if (dSP <= 0.) {
780  silent_cerr("warning, null or negative stiffness prime at line "
781  << HP.GetLineData() << std::endl);
782  }
783 
784  /* Prestress (no prestrain) */
785  T PreStress(mb_zero<T>());
786  GetPreStress(HP, PreStress);
787 
789  SAFENEWWITHCONSTRUCTOR(pCL, L, L(PreStress, dSP));
790 
791  return pCL;
792  };
793 };
794 
795 template <class T, class Tder>
797  virtual ConstitutiveLaw<T, Tder> *
798  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
799  ConstitutiveLaw<T, Tder>* pCL = 0;
800 
801  CLType = ConstLawType::VISCOUS;
802 
803  Tder SP(mb_zero<Tder>());
804  SP = HP.Get(SP);
805 
806  /* Prestress (no prestrain) */
807  T PreStress(mb_zero<T>());
808  GetPreStress(HP, PreStress);
809 
811  SAFENEWWITHCONSTRUCTOR(pCL, L, L(PreStress, SP));
812 
813  return pCL;
814  };
815 };
816 
817 template <class T, class Tder>
818 struct LinearViscoElasticCLR : public ConstitutiveLawRead<T, Tder> {
819  virtual ConstitutiveLaw<T, Tder> *
820  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
821  ConstitutiveLaw<T, Tder>* pCL = 0;
822 
824 
825  doublereal dS = HP.GetReal();
826  DEBUGCOUT("Stiffness = " << dS << std::endl);
827 
828  if (dS <= 0.) {
829  silent_cerr("warning, null or negative stiffness at line "
830  << HP.GetLineData() << std::endl);
831  }
832 
833  doublereal dSP = 0.;
834  if (HP.IsKeyWord("proportional")) {
835  doublereal k = HP.GetReal();
836  dSP = k*dS;
837  } else {
838  dSP = HP.GetReal();
839  }
840  DEBUGCOUT("stiffness prime = " << dSP << std::endl);
841 
842  if (dSP <= 0.) {
843  silent_cerr("warning, null or negative stiffness prime at line "
844  << HP.GetLineData() << std::endl);
845  }
846 
847  /* Prestress and prestrain */
848  T PreStress(mb_zero<T>());
849  GetPreStress(HP, PreStress);
850  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
851 
852 #if 0 // TODO: implement a "null" constitutive law that does nothing
853  if (dS == 0. && dSP == 0.) {
854 
855  } else
856 #endif
857  if (dS == 0.) {
858  silent_cerr("warning, null stiffness, "
859  "using linear viscous constitutive law instead"
860  << std::endl);
861 
863  SAFENEWWITHCONSTRUCTOR(pCL, L, L(PreStress, dSP));
864  if (pTplDC) {
865  delete pTplDC;
866  }
867 
868  } else if (dSP == 0.) {
869  silent_cerr("warning, null stiffness prime, "
870  "using linear elastic constitutive law instead"
871  << std::endl);
872 
874  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dS));
875 
876  } else {
878  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, dS, dSP));
879  }
880 
881  return pCL;
882  };
883 };
884 
885 template <class T, class Tder>
887  virtual ConstitutiveLaw<T, Tder> *
888  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
889  ConstitutiveLaw<T, Tder>* pCL = 0;
890 
892 
893  Tder S(mb_zero<Tder>());
894  S = HP.Get(S);
895 
896  Tder SP(mb_zero<Tder>());
897  if (HP.IsKeyWord("proportional")) {
898  doublereal k = HP.GetReal();
899  SP = S*k;
900 
901  } else {
902  SP = HP.Get(SP);
903  }
904 
905  /* Prestress and prestrain */
906  T PreStress(mb_zero<T>());
907  GetPreStress(HP, PreStress);
908  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
909 
910 #if 0 // TODO: implement a "null" constitutive law that does nothing
911  if (IsNull(S) && IsNull(SP)) {
912 
913  } else
914 #endif
915  if (IsNull(S)) {
916  silent_cerr("warning, null stiffness, "
917  "using linear viscous generic constitutive law instead"
918  << std::endl);
919 
921  SAFENEWWITHCONSTRUCTOR(pCL, L, L(PreStress, SP));
922  if (pTplDC) {
923  delete pTplDC;
924  }
925 
926  } else if (IsNull(SP)) {
927  silent_cerr("warning, null stiffness prime, "
928  "using linear elastic generic constitutive law instead"
929  << std::endl);
930 
932  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S));
933 
934  } else {
936  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S, SP));
937  }
938 
939  return pCL;
940  };
941 };
942 
943 template <class T, class Tder>
945  virtual ConstitutiveLaw<T, Tder> *
946  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
947  ConstitutiveLaw<T, Tder>* pCL = 0;
948 
950 
951  Tder S(mb_zero<Tder>());
952  S = HP.Get(S);
953 
954  DriveCaller *pdc = HP.GetDriveCaller();
955 
956  Tder SP(mb_zero<Tder>());
957  if (HP.IsKeyWord("proportional")) {
958  doublereal k = HP.GetReal();
959  SP = S*k;
960 
961  } else {
962  SP = HP.Get(SP);
963  }
964 
965  DriveCaller *pdcp = HP.GetDriveCaller();
966 
967  /* Prestress and prestrain */
968  T PreStress(mb_zero<T>());
969  GetPreStress(HP, PreStress);
970  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
971 
972 #if 0 // TODO: implement a "null" constitutive law that does nothing
973  if (IsNull(S) && IsNull(SP)) {
974 
975  } else if (IsNull(S)) {
976  silent_cerr("warning, null stiffness, "
977  "using linear viscous generic constitutive law instead"
978  << std::endl);
979 
980  SAFEDELETE(pdc);
981 
982  typedef LTVViscousGenericConstitutiveLaw<T, Tder> L;
983  SAFENEWWITHCONSTRUCTOR(pCL, L, L(PreStress, SP, pdcp));
984  if (pTplDC) {
985  delete pTplDC;
986  }
987 
988  } else
989  if (IsNull(SP)) {
990  silent_cerr("warning, null stiffness prime, "
991  "using linear elastic generic constitutive law instead"
992  << std::endl);
993 
994  SAFEDELETE(pdcp);
995 
996  typedef LTVElasticGenericConstitutiveLaw<T, Tder> L;
997  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S, pdc));
998 
999  } else
1000 #endif
1001  {
1003  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S, pdc, SP, pdcp));
1004  }
1005 
1006  return pCL;
1007  };
1008 };
1009 
1010 template <class T, class Tder>
1012  virtual ConstitutiveLaw<T, Tder> *
1013  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1014  ConstitutiveLaw<T, Tder>* pCL = 0;
1015 
1016  CLType = ConstLawType::ELASTIC;
1017 
1018  T S1(mb_zero<T>());
1019  S1 = HP.Get(S1);
1020 
1021  T S2(mb_zero<T>());
1022  S2 = HP.Get(S2);
1023 
1024  T S3(mb_zero<T>());
1025  S3 = HP.Get(S3);
1026 
1027  /* Prestress and prestrain */
1028  T PreStress(mb_zero<T>());
1029  GetPreStress(HP, PreStress);
1030  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1031 
1033  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S1, S2, S3));
1034 
1035  return pCL;
1036  };
1037 };
1038 
1039 template <class T, class Tder>
1041  virtual ConstitutiveLaw<T, Tder> *
1042  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1043  ConstitutiveLaw<T, Tder>* pCL = 0;
1044 
1045  CLType = ConstLawType::VISCOELASTIC;
1046 
1047  T S1(mb_zero<T>());
1048  S1 = HP.Get(S1);
1049 
1050  T S2(mb_zero<T>());
1051  S2 = HP.Get(S2);
1052 
1053  T S3(mb_zero<T>());
1054  S3 = HP.Get(S3);
1055 
1056  Tder SP(mb_zero<Tder>());
1057  SP = HP.Get(SP);
1058 
1059  /* Prestress and prestrain */
1060  T PreStress(mb_zero<T>());
1061  GetPreStress(HP, PreStress);
1062  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1063 
1065  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pTplDC, PreStress, S1, S2, S3, SP));
1066 
1067  return pCL;
1068  };
1069 };
1070 
1071 template <class T, class Tder>
1073  virtual ConstitutiveLaw<T, Tder> *
1074  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1075  ConstitutiveLaw<T, Tder>* pCL = 0;
1076 
1077  CLType = ConstLawType::VISCOELASTIC;
1078 
1079  doublereal dS = HP.GetReal();
1080  DEBUGCOUT("stiffness = " << dS << std::endl);
1081 
1082  if (dS <= 0.) {
1083  silent_cerr("warning, null or negative stiffness at line "
1084  << HP.GetLineData() << std::endl);
1085  }
1086 
1087  doublereal dUpp = HP.GetReal();
1088  if (dUpp <= 0.) {
1089  silent_cerr("warning, null or negative upper limit strain at line "
1090  << HP.GetLineData() << std::endl);
1091  }
1092 
1093  doublereal dLow = HP.GetReal();
1094  if (dLow >= 0.) {
1095  silent_cerr("warning, null or positive lower limit strain at line "
1096  << HP.GetLineData() << std::endl);
1097  }
1098 
1099  doublereal dSecondS = HP.GetReal();
1100  if (dSecondS <= 0.) {
1101  silent_cerr("warning, null or negative second stiffness at line "
1102  << HP.GetLineData() << std::endl);
1103  }
1104 
1105  doublereal dSP = HP.GetReal();
1106  DEBUGCOUT("stiffness prime = " << dSP << std::endl);
1107 
1108  if (dSP <= 0.) {
1109  silent_cerr("warning, null or negative stiffness prime at line "
1110  << HP.GetLineData() << std::endl);
1111  }
1112 
1113  doublereal dSecondSP = dSP;
1114  if (HP.IsKeyWord("second" "damping")) {
1115  dSecondSP = HP.GetReal();
1116  DEBUGCOUT("second stiffness prime = " << dSecondSP << std::endl);
1117 
1118  if (dSecondSP <= 0.) {
1119  silent_cerr("warning, null or negative second stiffness prime at line "
1120  << HP.GetLineData() << std::endl);
1121  }
1122  }
1123 
1124  /* Prestress and prestrain */
1125  T PreStress(mb_zero<T>());
1126  GetPreStress(HP, PreStress);
1127  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1128 
1131  L,
1132  L(pTplDC, PreStress,
1133  dS, dUpp, dLow, dSecondS, dSP, dSecondSP));
1134 
1135  return pCL;
1136  };
1137 };
1138 
1139 template <class T, class Tder>
1141  virtual ConstitutiveLaw<T, Tder> *
1142  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1143  ConstitutiveLaw<T, Tder>* pCL = 0;
1144 
1145  CLType = ConstLawType::VISCOELASTIC;
1146 
1147  doublereal dS = HP.GetReal();
1148  DEBUGCOUT("Visco-Elastic Turbulent Rod Joint, stiffness = "
1149  << dS << std::endl);
1150 
1151  if (dS <= 0.) {
1152  silent_cerr("warning, null or negative stiffness at line "
1153  << HP.GetLineData() << std::endl);
1154  }
1155 
1156  doublereal dParabStiff = HP.GetReal();
1157  DEBUGCOUT("stiffness prime = " << dParabStiff << std::endl);
1158 
1159  if (dParabStiff <= 0.) {
1160  silent_cerr("warning, null or negative derivative stiffness at line "
1161  << HP.GetLineData() << std::endl);
1162  }
1163 
1164  doublereal dTreshold = 0.;
1165  if (HP.IsArg()) {
1166  dTreshold = HP.GetReal(dTreshold);
1167 
1168  /*
1169  * Il legame costitutivo ha la forma seguente:
1170  * F = Kp*e + Kd*(de/dt)
1171  * con Kp costante e Kd dato dalla seguente legge:
1172  * Kd = cost2 per fabs(de/dt) < Treshold
1173  * Kd = 2*cost1*fabs(de/dt) per fabs(de/dt) > Treshold
1174  * se non viene inserito il valore di treshold, lo si
1175  * assume = 0. e quindi il legame e' sempre del secondo
1176  * tipo. Altrimenti, se non viene inserita la seconda
1177  * costante cost2, si assume che vi sia raccordo tra
1178  * i due tipi di legge, ovvero cost2 = cost1*Treshold
1179  * altrimenti e' possibile avere un comportamento,
1180  * che in prima approssimazione e' valido
1181  * per numerosi fluidi, in cui vi e' un salto tra i due
1182  * tipi di legge costitutiva. */
1183  }
1184 
1185  doublereal dSP = dTreshold*dParabStiff;
1186  if (HP.IsArg()) {
1187  dSP = HP.GetReal(dSP);
1188  }
1189 
1190  /* Prestress and prestrain */
1191  T PreStress(mb_zero<T>());
1192  GetPreStress(HP, PreStress);
1193  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1194 
1197  L,
1198  L(pTplDC, PreStress,
1199  dS, dSP, dTreshold, dParabStiff));
1200 
1201  return pCL;
1202  };
1203 };
1204 
1205 template <class T, class Tder>
1206 struct LinearBiStopCLR : public ConstitutiveLawRead<T, Tder> {
1207  virtual ConstitutiveLaw<T, Tder> *
1208  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1209  ConstitutiveLaw<T, Tder>* pCL = 0;
1210 
1211  typedef BiStopCLWrapper<T, Tder> L;
1214 
1215  DEBUGCOUT("Linear Viscoelastic Bi Stop Constitutive Law" << std::endl);
1216  doublereal dS = HP.GetReal();
1217  if (dS <= 0.) {
1218  silent_cerr("warning, null or negative stiffness at line "
1219  << HP.GetLineData() << std::endl);
1220  }
1221 
1222  doublereal dSp = 0.;
1223  if (CLType == ConstLawType::VISCOELASTIC) {
1224  dSp = HP.GetReal();
1225  if (dSp <= 0.) {
1226  silent_cerr("warning, null or negative stiffness prime at line "
1227  << HP.GetLineData() << std::endl);
1228  }
1229  }
1230 
1231  bool s(false);
1232  if (HP.IsKeyWord("initial" "status")) {
1233  if (HP.IsKeyWord("active")) {
1234  s = true;
1235 
1236  } else if (HP.IsKeyWord("inactive")) {
1237  s = false;
1238 
1239  } else {
1240  s = HP.GetBool();
1241  }
1242  }
1243 
1244  const DriveCaller *pA = HP.GetDriveCaller();
1245  const DriveCaller *pD = HP.GetDriveCaller();
1246 
1247  /* Prestress and prestrain */
1248  T PreStress(mb_zero<T>());
1249  GetPreStress(HP, PreStress);
1250  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1251 
1252  ConstitutiveLaw<T, Tder> *pWrappedCL = 0;
1253  if (CLType == ConstLawType::VISCOELASTIC && dSp != 0.) {
1254  SAFENEWWITHCONSTRUCTOR(pWrappedCL, LVECL, LVECL(pTplDC, PreStress, dS, dSp));
1255 
1256  } else {
1257  SAFENEWWITHCONSTRUCTOR(pWrappedCL, LECL, LECL(pTplDC, PreStress, dS));
1258  }
1259 
1260  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pWrappedCL, s, pA, pD));
1261 
1262  return pCL;
1263  };
1264 };
1265 
1266 template <class T, class Tder>
1267 struct LinearElasticBiStopCLR : public LinearBiStopCLR<T, Tder> {
1268  virtual ConstitutiveLaw<T, Tder> *
1269  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1270  CLType = ConstLawType::ELASTIC;
1271  return LinearBiStopCLR<T, Tder>::Read(pDM, HP, CLType);
1272  };
1273 };
1274 
1275 template <class T, class Tder>
1277  virtual ConstitutiveLaw<T, Tder> *
1278  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1279  CLType = ConstLawType::VISCOELASTIC;
1280  return LinearBiStopCLR<T, Tder>::Read(pDM, HP, CLType);
1281  };
1282 };
1283 
1284 static void
1285 ReadBiStopBase(MBDynParser& HP, bool& bStatus, const DriveCaller *& pA, const DriveCaller *& pD)
1286 {
1287  if (HP.IsKeyWord("initial" "status")) {
1288  if (HP.IsKeyWord("active")) {
1289  bStatus = true;
1290 
1291  } else if (HP.IsKeyWord("inactive")) {
1292  bStatus = false;
1293 
1294  } else {
1295  bStatus = HP.GetBool();
1296  }
1297  }
1298 
1299  pA = HP.GetDriveCaller();
1300  pD = HP.GetDriveCaller();
1301 }
1302 
1303 struct BiStopCLW1DR : public ConstitutiveLawRead<doublereal, doublereal> {
1305  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1307 
1309 
1310  const DriveCaller *pA = 0;
1311  const DriveCaller *pD = 0;
1312  bool bStatus(false);
1313  ReadBiStopBase(HP, bStatus, pA, pD);
1314 
1315  ConstitutiveLaw<doublereal, doublereal> *pWrappedCL = ReadCL1D(pDM, HP, CLType);
1316  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pWrappedCL, bStatus, pA, pD));
1317 
1318  return pCL;
1319  };
1320 };
1321 
1322 struct BiStopCLW3DR : public ConstitutiveLawRead<Vec3, Mat3x3> {
1324  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1326 
1328 
1329  const DriveCaller *pA = 0;
1330  const DriveCaller *pD = 0;
1331  bool bStatus(false);
1332  ReadBiStopBase(HP, bStatus, pA, pD);
1333 
1334  ConstitutiveLaw<Vec3, Mat3x3> *pWrappedCL = ReadCL3D(pDM, HP, CLType);
1335  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pWrappedCL, bStatus, pA, pD));
1336 
1337  return pCL;
1338  };
1339 };
1340 
1341 struct BiStopCLW6DR : public ConstitutiveLawRead<Vec6, Mat6x6> {
1343  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1345 
1347 
1348  const DriveCaller *pA = 0;
1349  const DriveCaller *pD = 0;
1350  bool bStatus(false);
1351  ReadBiStopBase(HP, bStatus, pA, pD);
1352 
1353  ConstitutiveLaw<Vec6, Mat6x6> *pWrappedCL = ReadCL6D(pDM, HP, CLType);
1354  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pWrappedCL, bStatus, pA, pD));
1355 
1356  return pCL;
1357  };
1358 };
1359 
1360 /*
1361  * Shock absorber per Stefy:
1362  *
1363  * ``Riprogettazione dell'ammortizzatore del carrello anteriore
1364  * di un velivolo di aviazione generale'',
1365  * S. Carlucci e S. Gualdi,
1366  * A.A. 1997-98
1367  */
1368 template <class T, class Tder>
1369 struct ShockAbsorberCLR : public ConstitutiveLawRead<T, Tder> {
1370  virtual ConstitutiveLaw<T, Tder> *
1371  Read(const DataManager* pDM, MBDynParser& HP, ConstLawType::Type& CLType) {
1372  ConstitutiveLaw<T, Tder>* pCL = 0;
1373 
1374  CLType = ConstLawType::VISCOELASTIC;
1375 
1376  TplDriveCaller<T>* pTplDC = GetPreStrain<T>(pDM, HP);
1377 
1379  SAFENEWWITHCONSTRUCTOR(pCL, L, L(pDM, pTplDC, HP));
1380 
1381  return pCL;
1382  };
1383 };
1384 
1385 static unsigned done = 0;
1386 
1387 /* initialization function */
1388 void
1389 InitCL(void)
1390 {
1391  if (::done++ > 0) {
1392  return;
1393  }
1394 
1395  /* constitutive law array */
1396  SetCL1D("array", new CLArray1DR);
1397  SetCL3D("array", new CLArray3DR);
1398  SetCL6D("array", new CLArray6DR);
1399 
1400  /* linear elastic */
1401  SetCL1D("linear" "elastic", new LinearElasticCLR<doublereal, doublereal>);
1402  SetCL3D("linear" "elastic", new LinearElasticCLR<Vec3, Mat3x3>);
1403  SetCL6D("linear" "elastic", new LinearElasticCLR<Vec6, Mat6x6>);
1404 
1405  /* linear elastic isotropic */
1406  SetCL1D("linear" "elastic" "isotropic", new LinearElasticCLR<doublereal, doublereal>);
1407  SetCL3D("linear" "elastic" "isotropic", new LinearElasticCLR<Vec3, Mat3x3>);
1408  SetCL6D("linear" "elastic" "isotropic", new LinearElasticCLR<Vec6, Mat6x6>);
1409 
1410  /* linear elastic generic */
1411  SetCL1D("linear" "elastic" "generic", new LinearElasticGenericCLR<doublereal, doublereal>);
1412  SetCL3D("linear" "elastic" "generic", new LinearElasticGenericCLR<Vec3, Mat3x3>);
1413  SetCL6D("linear" "elastic" "generic", new LinearElasticGenericCLR<Vec6, Mat6x6>);
1414 
1415  /* linear (visco)elastic generic axial torsion coupling*/
1416  SetCL6D("linear" "elastic" "generic" "axial" "torsion" "coupling",
1418  SetCL6D("linear" "viscoelastic" "generic" "axial" "torsion" "coupling",
1420 
1421  /* inverse square elastic */
1422  SetCL1D("inverse" "square" "elastic", new InverseSquareElasticCLR);
1423 
1424  /* log elastic */
1425  SetCL1D("log" "elastic", new LogElasticCLR<doublereal, doublereal>);
1426 
1427  /* double linear elastic */
1428  SetCL1D("double" "linear" "elastic", new DoubleLinearElasticCLR<doublereal, doublereal>);
1429  SetCL3D("double" "linear" "elastic", new DoubleLinearElasticCLR<Vec3, Mat3x3>);
1430 
1431  /* isotropic hardening elastic */
1432  SetCL1D("isotropic" "hardening" "elastic", new IsotropicHardeningCLR<doublereal, doublereal>);
1433  SetCL3D("isotropic" "hardening" "elastic", new IsotropicHardeningCLR<Vec3, Mat3x3>);
1434  SetCL6D("isotropic" "hardening" "elastic", new IsotropicHardeningCLR<Vec6, Mat6x6>);
1435 
1436  /* contact elastic */
1437  SetCL1D("contact" "elastic", new ContactElasticCLR<doublereal, doublereal>);
1438  SetCL3D("contact" "elastic", new ContactElasticCLR<Vec3, Mat3x3>);
1439 
1440  /* symbolic (elastic, viscous, viscoelastic) */
1442  SetCL3D("symbolic", new SymbolicCLR<Vec3, Mat3x3>);
1443  SetCL6D("symbolic", new SymbolicCLR<Vec6, Mat6x6>);
1444 
1445  SetCL1D("symbolic" "elastic", new SymbolicElasticCLR<doublereal, doublereal>);
1446  SetCL3D("symbolic" "elastic", new SymbolicElasticCLR<Vec3, Mat3x3>);
1447  SetCL6D("symbolic" "elastic", new SymbolicElasticCLR<Vec6, Mat6x6>);
1448 
1449  SetCL1D("symbolic" "viscous", new SymbolicViscousCLR<doublereal, doublereal>);
1450  SetCL3D("symbolic" "viscous", new SymbolicViscousCLR<Vec3, Mat3x3>);
1451  SetCL6D("symbolic" "viscous", new SymbolicViscousCLR<Vec6, Mat6x6>);
1452 
1453  SetCL1D("symbolic" "viscoelastic", new SymbolicViscoElasticCLR<doublereal, doublereal>);
1454  SetCL3D("symbolic" "viscoelastic", new SymbolicViscoElasticCLR<Vec3, Mat3x3>);
1455  SetCL6D("symbolic" "viscoelastic", new SymbolicViscoElasticCLR<Vec6, Mat6x6>);
1456 
1457  /* linear viscous */
1458  SetCL1D("linear" "viscous", new LinearViscousCLR<doublereal, doublereal>);
1459  SetCL3D("linear" "viscous", new LinearViscousCLR<Vec3, Mat3x3>);
1460  SetCL6D("linear" "viscous", new LinearViscousCLR<Vec6, Mat6x6>);
1461 
1462  /* linear viscous isotropic */
1463  SetCL1D("linear" "viscous" "isotropic", new LinearViscousCLR<doublereal, doublereal>);
1464  SetCL3D("linear" "viscous" "isotropic", new LinearViscousCLR<Vec3, Mat3x3>);
1465  SetCL6D("linear" "viscous" "isotropic", new LinearViscousCLR<Vec6, Mat6x6>);
1466 
1467  /* linear viscous generic */
1468  SetCL1D("linear" "viscous" "generic", new LinearViscousGenericCLR<doublereal, doublereal>);
1469  SetCL3D("linear" "viscous" "generic", new LinearViscousGenericCLR<Vec3, Mat3x3>);
1470  SetCL6D("linear" "viscous" "generic", new LinearViscousGenericCLR<Vec6, Mat6x6>);
1471 
1472  /* linear viscoelastic */
1473  SetCL1D("linear" "viscoelastic", new LinearViscoElasticCLR<doublereal, doublereal>);
1474  SetCL3D("linear" "viscoelastic", new LinearViscoElasticCLR<Vec3, Mat3x3>);
1475  SetCL6D("linear" "viscoelastic", new LinearViscoElasticCLR<Vec6, Mat6x6>);
1476 
1477  /* linear viscoelastic isotropic */
1478  SetCL1D("linear" "viscoelastic" "isotropic", new LinearViscoElasticCLR<doublereal, doublereal>);
1479  SetCL3D("linear" "viscoelastic" "isotropic", new LinearViscoElasticCLR<Vec3, Mat3x3>);
1480  SetCL6D("linear" "viscoelastic" "isotropic", new LinearViscoElasticCLR<Vec6, Mat6x6>);
1481 
1482  /* linear viscoelastic generic */
1483  SetCL1D("linear" "viscoelastic" "generic", new LinearViscoElasticGenericCLR<doublereal, doublereal>);
1484  SetCL3D("linear" "viscoelastic" "generic", new LinearViscoElasticGenericCLR<Vec3, Mat3x3>);
1485  SetCL6D("linear" "viscoelastic" "generic", new LinearViscoElasticGenericCLR<Vec6, Mat6x6>);
1486 
1487  /* linear time variant viscoelastic generic */
1488  SetCL1D("linear" "time" "variant" "viscoelastic" "generic", new LTVViscoElasticGenericCLR<doublereal, doublereal>);
1489  SetCL3D("linear" "time" "variant" "viscoelastic" "generic", new LTVViscoElasticGenericCLR<Vec3, Mat3x3>);
1490  SetCL6D("linear" "time" "variant" "viscoelastic" "generic", new LTVViscoElasticGenericCLR<Vec6, Mat6x6>);
1491 
1492  /* cubic elastic generic */
1493  SetCL1D("cubic" "elastic" "generic", new CubicElasticGenericCLR<doublereal, doublereal>);
1494  SetCL3D("cubic" "elastic" "generic", new CubicElasticGenericCLR<Vec3, Mat3x3>);
1495 
1496  /* cubic viscoelastic generic */
1497  SetCL1D("cubic" "viscoelastic" "generic", new CubicViscoElasticGenericCLR<doublereal, doublereal>);
1498  SetCL3D("cubic" "viscoelastic" "generic", new CubicViscoElasticGenericCLR<Vec3, Mat3x3>);
1499 
1500  /* double linear viscoelastic */
1501  SetCL1D("double" "linear" "viscoelastic", new DoubleLinearViscoElasticCLR<doublereal, doublereal>);
1502  SetCL3D("double" "linear" "viscoelastic", new DoubleLinearViscoElasticCLR<Vec3, Mat3x3>);
1503 
1504  /* turbulent viscoelastic */
1505  SetCL1D("turbulent" "viscoelastic", new TurbulentViscoElasticCLR<doublereal, doublereal>);
1506 
1507  /* linear elastic bistop */
1508  SetCL1D("linear" "elastic" "bistop", new LinearElasticBiStopCLR<doublereal, doublereal>);
1509  SetCL3D("linear" "elastic" "bistop", new LinearElasticBiStopCLR<Vec3, Mat3x3>);
1510  SetCL6D("linear" "elastic" "bistop", new LinearElasticBiStopCLR<Vec6, Mat6x6>);
1511 
1512  /* linear viscoelastic bistop */
1513  SetCL1D("linear" "viscoelastic" "bistop", new LinearViscoElasticBiStopCLR<doublereal, doublereal>);
1514  SetCL3D("linear" "viscoelastic" "bistop", new LinearViscoElasticBiStopCLR<Vec3, Mat3x3>);
1515  SetCL6D("linear" "viscoelastic" "bistop", new LinearViscoElasticBiStopCLR<Vec6, Mat6x6>);
1516 
1517  /* bistop wrapper */
1518  SetCL1D("bistop", new BiStopCLW1DR);
1519  SetCL3D("bistop", new BiStopCLW3DR);
1520  SetCL6D("bistop", new BiStopCLW6DR);
1521 
1522  /* shock absorber */
1523  SetCL1D("shock" "absorber", new ShockAbsorberCLR<doublereal, doublereal>);
1524 
1525  /* Artificial Neural Network */
1526  SetCL1D("ann" "elastic", new AnnElasticCLR<doublereal, doublereal>);
1527  SetCL1D("ann" "viscoelastic", new AnnViscoElasticCLR<doublereal, doublereal>);
1528 
1529  /* constitutive laws sponsored by Hutchinson CdR */
1530  NLP_init();
1531  NLSF_init();
1532 
1533  /* invariant constitutive law */
1534  SetCL3D("invariant" "angular", new InvAngularCLR);
1535  SetCL3D("axial" "wrapper", new AxialCLR);
1536 
1537  /* ... */
1538  TDCLW_init();
1539 
1540  /* NOTE: add here initialization of new built-in constitutive laws;
1541  * alternative ways to register new custom constitutive laws are:
1542  * - call SetCL*D() from anywhere in the code
1543  * - write a module that calls SetCL*D() from inside a function
1544  * called module_init(), and load it using "module load".
1545  */
1546 }
1547 
1548 void
1550 {
1551  if (::done == 0) {
1552  silent_cerr("DestroyCL() called once too many" << std::endl);
1554  }
1555 
1556  if (--::done > 0) {
1557  return;
1558  }
1559 
1560  /* free stuff */
1561  for (CL1DFuncMapType::iterator i = CL1DFuncMap.begin(); i != CL1DFuncMap.end(); ++i) {
1562  delete i->second;
1563  }
1564  CL1DFuncMap.clear();
1565 
1566  for (CL3DFuncMapType::iterator i = CL3DFuncMap.begin(); i != CL3DFuncMap.end(); ++i) {
1567  delete i->second;
1568  }
1569  CL3DFuncMap.clear();
1570 
1571  for (CL6DFuncMapType::iterator i = CL6DFuncMap.begin(); i != CL6DFuncMap.end(); ++i) {
1572  delete i->second;
1573  }
1574  CL6DFuncMap.clear();
1575 }
static CL1DWordSetType CL1DWordSet
virtual ConstitutiveLaw< Vec3, Mat3x3 > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
static CL3DFuncMapType CL3DFuncMap
Definition: matvec3.h:98
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual integer GetInt(integer iDefval=0)
Definition: parser.cc:1050
bool SetCL3D(const char *name, ConstitutiveLawRead< Vec3, Mat3x3 > *rf)
virtual ConstitutiveLaw< Vec6, Mat6x6 > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual const char * IsWord(const HighParser::WordSet &ws)
Definition: parser.cc:977
bool IsNull(const doublereal &d)
Definition: matvec3.cc:1124
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
ConstitutiveLaw< Vec6, Mat6x6 > * ReadCL6D(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
void TDCLW_init(void)
Definition: tdclw.cc:236
bool IsWord(const std::string &s) const
virtual bool GetBool(bool bDefval=false)
Definition: parser.cc:1010
static CL6DFuncMapType CL6DFuncMap
virtual ConstitutiveLaw< doublereal, doublereal > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
bool SetCL6D(const char *name, ConstitutiveLawRead< Vec6, Mat6x6 > *rf)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
ConstitutiveLaw< doublereal, doublereal > * ReadCL1D(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< doublereal, doublereal > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
bool SetCL1D(const char *name, ConstitutiveLawRead< doublereal, doublereal > *rf)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
static CL6DWordSetType CL6DWordSet
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
static CL3DWordSetType CL3DWordSet
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
void func(const T &u, const T &v, const T &w, doublereal e, T &f)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< Vec6, Mat6x6 > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
static unsigned done
int NLSF_init(void)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
Definition: matvec6.h:37
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstLawType::Type GetConstLawType(void) const =0
virtual bool IsKeyWord(const char *sKeyWord)
Definition: parser.cc:910
void InitCL(void)
#define DEBUGCOUT(msg)
Definition: myassert.h:232
virtual const char * GetStringWithDelims(enum Delims Del=DEFAULTDELIM, bool escape=true)
Definition: parser.cc:1228
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
bool IsWord(const std::string &s) const
std::map< std::string, ConstitutiveLawRead< doublereal, doublereal > *, ltstrcase > CL1DFuncMapType
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
int NLP_init(void)
#define ASSERT(expression)
Definition: colamd.c:977
#define SAFENEWWITHCONSTRUCTOR(pnt, item, constructor)
Definition: mynewmem.h:698
void GetPreStress(MBDynParser &HP, T &PreStress)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
std::map< std::string, ConstitutiveLawRead< Vec3, Mat3x3 > *, ltstrcase > CL3DFuncMapType
static CL1DFuncMapType CL1DFuncMap
const doublereal dS
Definition: beamslider.cc:71
virtual bool IsArg(void)
Definition: parser.cc:807
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
DriveCaller * GetDriveCaller(bool bDeferred=false)
Definition: mbpar.cc:2033
ConstitutiveLaw< Vec3, Mat3x3 > * ReadCL3D(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual doublereal Get(const doublereal &d)
Definition: mbpar.cc:2213
double doublereal
Definition: colamd.c:52
const doublereal & mb_zero< doublereal >(void)
Definition: matvec3.h:1554
virtual ConstitutiveLaw< doublereal, doublereal > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
virtual HighParser::ErrOut GetLineData(void) const
Definition: parsinc.cc:697
void DestroyCL(void)
bool IsWord(const std::string &s) const
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
static void ReadBiStopBase(MBDynParser &HP, bool &bStatus, const DriveCaller *&pA, const DriveCaller *&pD)
virtual ConstitutiveLaw< Vec3, Mat3x3 > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
#define SAFEDELETE(pnt)
Definition: mynewmem.h:710
virtual ConstitutiveLaw< T, Tder > * Read(const DataManager *pDM, MBDynParser &HP, ConstLawType::Type &CLType)
std::map< std::string, ConstitutiveLawRead< Vec6, Mat6x6 > *, ltstrcase > CL6DFuncMapType
virtual doublereal GetReal(const doublereal &dDefval=0.0)
Definition: parser.cc:1056