1 /*! 2 * \file mtest/src/AsterSmallStrainBehaviour.cxx 3 * \brief 4 * \author Thomas Helfer 5 * \brief 07 avril 2013 6 * \copyright Copyright (C) 2006-2018 CEA/DEN, EDF R&D. All rights 7 * reserved. 8 * This project is publicly released under either the GNU GPL Licence 9 * or the CECILL-A licence. A copy of thoses licences are delivered 10 * with the sources of TFEL. CEA or EDF may also distribute this 11 * project under specific licensing conditions. 12 */ 13 14 #include<cmath> 15 #include<ostream> 16 #include<algorithm> 17 18 #include"TFEL/Raise.hxx" 19 #include"TFEL/Math/tmatrix.hxx" 20 #include"TFEL/Math/stensor.hxx" 21 #include"TFEL/Math/st2tost2.hxx" 22 #include"MFront/MFrontLogStream.hxx" 23 #include"MFront/Aster/Aster.hxx" 24 #include"MFront/Aster/AsterComputeStiffnessTensor.hxx" 25 26 #include"MTest/CurrentState.hxx" 27 #include"MTest/BehaviourWorkSpace.hxx" 28 #include"MTest/UmatNormaliseTangentOperator.hxx" 29 #include"MTest/AsterSmallStrainBehaviour.hxx" 30 31 namespace mtest 32 { 33 AsterSmallStrainBehaviour(const Hypothesis h,const std::string & l,const std::string & b)34 AsterSmallStrainBehaviour::AsterSmallStrainBehaviour(const Hypothesis h, 35 const std::string& l, 36 const std::string& b) 37 : AsterStandardBehaviour(h,l,b) 38 {} // end of AsterSmallStrainBehaviour::AsterSmallStrainBehaviour 39 40 void getGradientsDefaultInitialValues(tfel::math::vector<real> & v) const41 AsterSmallStrainBehaviour::getGradientsDefaultInitialValues(tfel::math::vector<real>& v) const 42 { 43 std::fill(v.begin(),v.end(),real(0)); 44 } 45 46 std::pair<bool,real> call_behaviour(tfel::math::matrix<real> & Kt,CurrentState & s,BehaviourWorkSpace & wk,const real dt,const StiffnessMatrixType ktype,const bool b) const47 AsterSmallStrainBehaviour::call_behaviour(tfel::math::matrix<real>& Kt, 48 CurrentState& s, 49 BehaviourWorkSpace& wk, 50 const real dt, 51 const StiffnessMatrixType ktype, 52 const bool b) const 53 { 54 using namespace std; 55 using namespace tfel::math; 56 using namespace aster; 57 using tfel::math::vector; 58 constexpr const auto sqrt2 = Cste<real>::sqrt2; 59 unsigned short dimension; 60 AsterInt ntens; 61 AsterInt nprops = s.mprops1.size() == 0 ? 1 : static_cast<AsterInt>(s.mprops1.size()); 62 AsterInt nstatv; 63 AsterInt nummod; 64 const auto h = this->getHypothesis(); 65 if (h==ModellingHypothesis::AXISYMMETRICAL){ 66 ntens = 4; 67 dimension = 2u; 68 nummod = 4u; 69 } else if (h==ModellingHypothesis::PLANESTRESS){ 70 ntens = 4; 71 dimension = 2u; 72 nummod = 5u; 73 } else if (h==ModellingHypothesis::PLANESTRAIN){ 74 ntens = 4; 75 dimension = 2u; 76 nummod = 6u; 77 } else if (h==ModellingHypothesis::TRIDIMENSIONAL){ 78 ntens = 6; 79 dimension = 3u; 80 nummod = 3u; 81 } else { 82 tfel::raise("AsterSmallStrainBehaviour::call_behaviour: " 83 "unsupported hypothesis"); 84 } 85 fill(wk.D.begin(),wk.D.end(),0.); 86 // choosing the type of stiffness matrix 87 StandardBehaviourBase::initializeTangentOperator(wk.D,ktype,b); 88 // using a local copy of material properties to handle the 89 // case where s.mprops1 is empty 90 copy(s.mprops1.begin(),s.mprops1.end(),wk.mps.begin()); 91 if(s.mprops1.empty()){ 92 wk.mps[0] = real(0); 93 } 94 // using a local copy of internal state variables to handle the 95 // case where iv0 is empty 96 copy(s.iv0.begin(),s.iv0.end(),wk.ivs.begin()); 97 if(s.iv0.empty()){ 98 wk.ivs[0] = real(0); 99 } 100 nstatv = static_cast<AsterInt>(wk.ivs.size()); 101 // rotation matrix 102 tmatrix<3u,3u,real> drot = transpose(s.r); 103 stensor<3u,real> ue0(real(0)); 104 stensor<3u,real> ude(real(0)); 105 copy(s.e0.begin(),s.e0.end(),ue0.begin()); 106 for(decltype(s.e1.size()) i=0;i!=s.e1.size();++i){ 107 ude(i) = s.e1(i)-s.e0(i); 108 } 109 copy(s.s0.begin(),s.s0.end(),s.s1.begin()); 110 // thermal strain 111 for(AsterInt i=0;i!=static_cast<unsigned short>(ntens);++i){ 112 ue0(i) -= s.e_th0(i); 113 ude(i) -= s.e_th1(i)-s.e_th0(i); 114 } 115 // aster conventions 116 for(AsterInt i=3;i!=static_cast<unsigned short>(ntens);++i){ 117 s.s1(i) /= sqrt2; 118 ue0(i) *= sqrt2; 119 ude(i) *= sqrt2; 120 } 121 AsterReal ndt(1.); 122 (this->fct)(&(s.s1(0)),&(wk.ivs(0)),&(wk.D(0,0)), 123 &ue0(0),&ude(0),&dt, 124 &(s.esv0(0)),&(s.desv(0)), 125 &(s.esv0(0))+1,&(s.desv(0))+1, 126 &ntens,&nstatv,&(wk.mps(0)), 127 &nprops,&drot(0,0),&ndt,&nummod); 128 if(ndt<1){ 129 if(mfront::getVerboseMode()>=mfront::VERBOSE_LEVEL1){ 130 if(this->emsg!=nullptr){ 131 mfront::getLogStream() << this->emsg() << std::endl; 132 } 133 } 134 return {false,ndt}; 135 } 136 if(ktype!=StiffnessMatrixType::NOSTIFFNESS){ 137 UmatNormaliseTangentOperator::exe(&Kt(0,0),wk.D,dimension); 138 } 139 if(b){ 140 if(!s.iv0.empty()){ 141 copy_n(wk.ivs.begin(), s.iv1.size(),s.iv1.begin()); 142 } 143 // turning things in standard conventions 144 for(AsterInt i=3;i!=static_cast<unsigned short>(ntens);++i){ 145 s.s1(i) *= sqrt2; 146 } 147 } 148 return {true,ndt}; 149 } 150 151 AsterSmallStrainBehaviour::~AsterSmallStrainBehaviour() = default; 152 153 } // end of namespace mtest 154