1 /*! 2 * \file mfront/src/Cazacu2004OrthotropicStressCriterion.cxx 3 * \brief 4 * \author Thomas Helfer 5 * \date 15/03/2018 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 "TFEL/Raise.hxx" 15 #include "MFront/BehaviourBrick/BrickUtilities.hxx" 16 #include "MFront/BehaviourBrick/StressPotential.hxx" 17 #include "MFront/BehaviourBrick/OptionDescription.hxx" 18 #include "MFront/BehaviourBrick/Cazacu2004OrthotropicStressCriterion.hxx" 19 20 namespace mfront { 21 22 namespace bbrick { 23 24 std::vector<OptionDescription> getOptions() const25 Cazacu2004OrthotropicStressCriterion::getOptions() const { 26 auto opts = StressCriterionBase::getOptions(); 27 opts.emplace_back("a", "coefficients of J20", 28 OptionDescription::ARRAYOFMATERIALPROPERTIES); 29 opts.emplace_back("b", "coefficients of J30", 30 OptionDescription::ARRAYOFMATERIALPROPERTIES); 31 opts.emplace_back("c", "", OptionDescription::MATERIALPROPERTY); 32 return opts; 33 } // end of Cazacu2004OrthotropicStressCriterion::getOptions() 34 35 std::vector<mfront::BehaviourSymmetryType> getSupportedBehaviourSymmetries() const36 Cazacu2004OrthotropicStressCriterion::getSupportedBehaviourSymmetries() 37 const { 38 return {mfront::ORTHOTROPIC}; 39 } // end of 40 // Cazacu2004OrthotropicStressCriterion::getSupportedBehaviourSymmetries() 41 initialize(BehaviourDescription & bd,AbstractBehaviourDSL & dsl,const std::string & id,const DataMap & d,const Role r)42 void Cazacu2004OrthotropicStressCriterion::initialize( 43 BehaviourDescription& bd, 44 AbstractBehaviourDSL& dsl, 45 const std::string& id, 46 const DataMap& d, 47 const Role r) { 48 using ConstantMaterialProperty = 49 mfront::BehaviourDescription::ConstantMaterialProperty; 50 tfel::raise_if(bd.getSymmetryType() != mfront::ORTHOTROPIC, 51 "Cazacu2004OrthotropicStressCriterion::initialize: " 52 "the behaviour must be orthotropic"); 53 StressCriterionBase::initialize(bd, dsl, id, d, r); 54 const auto an = StressCriterion::getVariableId("a", id, r); 55 const auto bn = StressCriterion::getVariableId("b", id, r); 56 const auto cn = StressCriterion::getVariableId("c", id, r); 57 tfel::raise_if(d.count("c") == 0, 58 "Cazacu2004OrthotropicStressCriterion::initialize: " 59 "material property 'c' is not defined"); 60 bd.appendToIncludes( 61 "#include\"TFEL/Material/Cazacu2004OrthotropicYieldCriterion.hxx\""); 62 // J2O coefficients 63 if (d.count("a") == 0) { 64 tfel::raise( 65 "Cazacu2004OrthotropicStressCriterion::initialize: entry '" + 66 std::string("a") + "' is not defined"); 67 } 68 this->a = getArrayOfBehaviourDescriptionMaterialProperties<6u>(dsl, "a", 69 d.at("a")); 70 if (this->a[0].is<ConstantMaterialProperty>()) { 71 std::vector<double> values(6u); 72 for (unsigned short i = 0; i != 6; ++i) { 73 tfel::raise_if(!this->a[i].is<ConstantMaterialProperty>(), 74 "Cazacu2004OrthotropicStressCriterion::initialize: " 75 "if one component of '" + 76 std::string("a") + 77 "' is a constant value, all components must be " 78 "constant values"); 79 values[i] = this->a[i].get<ConstantMaterialProperty>().value; 80 } 81 const auto ean = [&r] { 82 if (r == StressCriterion::STRESSANDFLOWCRITERION) { 83 return "J2OCoefficients"; 84 } else if (r == StressCriterion::STRESSCRITERION) { 85 return "StressCriterionJ2OCoefficients"; 86 } 87 return "FlowCriterionJ2OCoefficients"; 88 }(); 89 addParameter(bd, an, ean, 6, values); 90 } else { 91 addLocalVariable(bd, "real", an, 6); 92 } 93 // J3O coefficients 94 if (d.count("b") == 0) { 95 tfel::raise( 96 "Cazacu2004OrthotropicStressCriterion::initialize: entry '" + 97 std::string("b") + "' is not defined"); 98 } 99 this->b = getArrayOfBehaviourDescriptionMaterialProperties<11u>( 100 dsl, "b", d.at("b")); 101 if (this->b[0].is<ConstantMaterialProperty>()) { 102 std::vector<double> values(11u); 103 for (unsigned short i = 0; i != 11; ++i) { 104 tfel::raise_if(!this->b[i].is<ConstantMaterialProperty>(), 105 "Cazacu2004OrthotropicStressCriterion::initialize: " 106 "if one component of '" + 107 std::string("b") + 108 "' is a constant value, all components must be " 109 "constant values"); 110 values[i] = this->b[i].get<ConstantMaterialProperty>().value; 111 } 112 const auto ebn = [&r] { 113 if (r == StressCriterion::STRESSANDFLOWCRITERION) { 114 return "J3OCoefficients"; 115 } else if (r == StressCriterion::STRESSCRITERION) { 116 return "StressCriterionJ3OCoefficients"; 117 } 118 return "FlowCriterionJ3OCoefficients"; 119 }(); 120 addParameter(bd, bn, ebn, 11, values); 121 } else { 122 addLocalVariable(bd, "real", bn, 11); 123 } 124 this->cp = getBehaviourDescriptionMaterialProperty(dsl, "c", d.at("c")); 125 declareParameterOrLocalVariable(bd, this->cp, "real", cn); 126 } // end of Cazacu2004OrthotropicStressCriterion::initialize 127 endTreatment(BehaviourDescription & bd,const AbstractBehaviourDSL & dsl,const std::string & id,const Role r)128 void Cazacu2004OrthotropicStressCriterion::endTreatment( 129 BehaviourDescription& bd, 130 const AbstractBehaviourDSL& dsl, 131 const std::string& id, 132 const Role r) { 133 constexpr const auto uh = 134 tfel::material::ModellingHypothesis::UNDEFINEDHYPOTHESIS; 135 const auto an = StressCriterion::getVariableId("a", id, r); 136 const auto bn = StressCriterion::getVariableId("b", id, r); 137 const auto cn = StressCriterion::getVariableId("c", id, r); 138 if (bd.getOrthotropicAxesConvention() == 139 tfel::material::OrthotropicAxesConvention::PIPE) { 140 for (const auto mh : bd.getModellingHypotheses()) { 141 tfel::raise_if( 142 ((mh != ModellingHypothesis::TRIDIMENSIONAL) && 143 (mh != ModellingHypothesis::AXISYMMETRICAL) && 144 (mh != 145 ModellingHypothesis::AXISYMMETRICALGENERALISEDPLANESTRAIN) && 146 (mh != 147 ModellingHypothesis::AXISYMMETRICALGENERALISEDPLANESTRESS)), 148 "Cazacu2004OrthotropicStressCriterion::endTreatment: " 149 "with the Pipe orthotropic axis convention, " 150 "only the following modelling hypotheses are supported:\n" 151 "- Tridimensional\n" 152 "- Axisymmetrical\n" 153 "- AxisymmetricalGeneralisedPlaneStrain\n" 154 "- AxisymmetricalGeneralisedPlaneStress."); 155 } 156 } 157 if (bd.getOrthotropicAxesConvention() == 158 tfel::material::OrthotropicAxesConvention::PLATE) { 159 for (const auto mh : bd.getModellingHypotheses()) { 160 tfel::raise_if( 161 ((mh != ModellingHypothesis::TRIDIMENSIONAL) && 162 (mh != ModellingHypothesis::PLANESTRESS)), 163 "Cazacu2004OrthotropicStressCriterion::endTreatment: " 164 "with the Plate orthotropic axis convention, " 165 "only the following modelling hypotheses are supported:\n" 166 "- Tridimensional\n" 167 "- Plane stress"); 168 } 169 } else { 170 tfel::raise_if( 171 bd.getOrthotropicAxesConvention() != 172 tfel::material::OrthotropicAxesConvention::DEFAULT, 173 "Cazacu2004OrthotropicStressCriterion::endTreatment: " 174 "internal error, unsupported orthotropic axes convention"); 175 for (const auto mh : bd.getModellingHypotheses()) { 176 tfel::raise_if(mh != ModellingHypothesis::TRIDIMENSIONAL, 177 "Cazacu2004OrthotropicStressCriterion::endTreatment: " 178 "an orthotropic axes convention must be choosen when " 179 "defining a stress free expansion in behaviours " 180 "which shall be valid in other modelling hypothesis " 181 "than 'Tridimensional'.\n" 182 "Either restrict the validity of the behaviour to " 183 "'Tridimensional' (see @ModellingHypothesis) or " 184 "choose and orthotropic axes convention as on option " 185 "to the @OrthotropicBehaviour keyword"); 186 } 187 } 188 auto c = std::string{}; 189 c += generateMaterialPropertiesInitializationCode(dsl, bd, an, this->a); 190 c += generateMaterialPropertiesInitializationCode(dsl, bd, bn, this->b); 191 c += generateMaterialPropertyInitializationCode(dsl, bd, cn, this->cp); 192 if (!c.empty()) { 193 CodeBlock i; 194 i.code = c; 195 bd.setCode(uh, BehaviourData::BeforeInitializeLocalVariables, i, 196 BehaviourData::CREATEORAPPEND, BehaviourData::AT_BEGINNING); 197 } 198 } // end of Cazacu2004OrthotropicStressCriterion::endTreatment 199 computeElasticPrediction(const std::string & id,const BehaviourDescription &,const StressPotential &) const200 std::string Cazacu2004OrthotropicStressCriterion::computeElasticPrediction( 201 const std::string& id, 202 const BehaviourDescription&, 203 const StressPotential&) const { 204 const auto an = StressCriterion::getVariableId( 205 "a", id, StressCriterion::STRESSCRITERION); 206 const auto bn = StressCriterion::getVariableId( 207 "b", id, StressCriterion::STRESSCRITERION); 208 const auto cn = StressCriterion::getVariableId( 209 "c", id, StressCriterion::STRESSCRITERION); 210 return "const auto seqel" + id + 211 " = computeCazacu2004OrthotropicStressCriterion(sel" + id + 212 ",this->" + an + ",this->" + bn + ",this->" + cn + ");\n"; 213 } // end of Cazacu2004OrthotropicStressCriterion::computeElasticPrediction 214 computeCriterion(const std::string & id,const BehaviourDescription &,const StressPotential &) const215 std::string Cazacu2004OrthotropicStressCriterion::computeCriterion( 216 const std::string& id, 217 const BehaviourDescription&, 218 const StressPotential&) const { 219 const auto an = StressCriterion::getVariableId( 220 "a", id, StressCriterion::STRESSCRITERION); 221 const auto bn = StressCriterion::getVariableId( 222 "b", id, StressCriterion::STRESSCRITERION); 223 const auto cn = StressCriterion::getVariableId( 224 "c", id, StressCriterion::STRESSCRITERION); 225 return "const auto seq" + id + 226 " = computeCazacu2004OrthotropicStressCriterion(s" + id + 227 ",this->" + an + ",this->" + bn + ",this->" + cn + ");\n"; 228 } // end of Cazacu2004OrthotropicStressCriterion::computeCriterion 229 computeNormal(const std::string & id,const BehaviourDescription & bd,const StressPotential & sp,const Role r) const230 std::string Cazacu2004OrthotropicStressCriterion::computeNormal( 231 const std::string& id, 232 const BehaviourDescription& bd, 233 const StressPotential& sp, 234 const Role r) const { 235 const auto an = StressCriterion::getVariableId( 236 "a", id, StressCriterion::STRESSCRITERION); 237 const auto bn = StressCriterion::getVariableId( 238 "b", id, StressCriterion::STRESSCRITERION); 239 const auto cn = StressCriterion::getVariableId("c", id, r); 240 auto c = std::string{}; 241 if ((r == STRESSCRITERION) || (r == STRESSANDFLOWCRITERION)) { 242 #if __cplusplus >= 201703L 243 c += "const auto [seq" + id + ",dseq" + id + "_ds" + id + "] = "; 244 c += "computeCazacu2004OrthotropicStressCriterionNormal(s" + id + 245 ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 246 sp.getEquivalentStressLowerBound(bd) + ");\n"; 247 #else /* __cplusplus >= 201703L */ 248 c += "stress seq" + id + ";\n"; 249 c += "Stensor dseq" + id + "_ds" + id + ";\n"; 250 c += "std::tie(seq" + id + ",dseq" + id + "_ds" + id + ") = "; 251 c += "computeCazacu2004OrthotropicStressCriterionNormal(s" + id + 252 ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 253 sp.getEquivalentStressLowerBound(bd) + ");\n"; 254 #endif /* __cplusplus >= 201703L */ 255 } 256 if (r == STRESSANDFLOWCRITERION) { 257 c += "const auto& n" + id + " = dseq" + id + "_ds" + id + ";\n"; 258 } 259 if (r == FLOWCRITERION) { 260 #if __cplusplus >= 201703L 261 c += "const auto [seqf" + id + ", n" + id + "] = "; 262 c += "computeCazacu2004OrthotropicStressCriterionNormal(s" + id + 263 ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 264 sp.getEquivalentStressLowerBound(bd) + ");\n"; 265 #else /* __cplusplus >= 201703L */ 266 c += "stress seqf" + id + ";\n"; 267 c += "Stensor n" + id + ";\n"; 268 c += "std::tie(seqf" + id + ",n" + id + ") = "; 269 c += "computeCazacu2004OrthotropicStressCriterionNormal(s" + id + 270 ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 271 sp.getEquivalentStressLowerBound(bd) + ");\n"; 272 #endif /* __cplusplus >= 201703L */ 273 } 274 return c; 275 } // end of Cazacu2004OrthotropicStressCriterion::computeNormal 276 computeNormalDerivative(const std::string & id,const BehaviourDescription & bd,const StressPotential & sp,const Role r) const277 std::string Cazacu2004OrthotropicStressCriterion::computeNormalDerivative( 278 const std::string& id, 279 const BehaviourDescription& bd, 280 const StressPotential& sp, 281 const Role r) const { 282 const auto an = StressCriterion::getVariableId( 283 "a", id, StressCriterion::STRESSCRITERION); 284 const auto bn = StressCriterion::getVariableId( 285 "b", id, StressCriterion::STRESSCRITERION); 286 const auto cn = StressCriterion::getVariableId("c", id, r); 287 auto c = std::string{}; 288 if ((r == STRESSCRITERION) || (r == STRESSANDFLOWCRITERION)) { 289 #if __cplusplus >= 201703L 290 c += "const auto [seq" + id + ",dseq" + id + "_ds" + id + ",d2seq" + 291 id + "_ds" + id + "ds" + id + "] = "; 292 c += "computeCazacu2004OrthotropicStressCriterionSecondDerivative(s" + 293 id + ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 294 sp.getEquivalentStressLowerBound(bd) + ");\n"; 295 #else /* __cplusplus >= 201703L */ 296 c += "stress seq" + id + ";\n"; 297 c += "Stensor dseq" + id + "_ds" + id + ";\n"; 298 c += "Stensor4 d2seq" + id + "_ds" + id + "ds" + id + ";\n"; 299 c += "std::tie(seq" + id + ",dseq" + id + "_ds" + id + ",d2seq" + id + 300 "_ds" + id + "ds" + id + ") = "; 301 c += "computeCazacu2004OrthotropicStressCriterionSecondDerivative(s" + 302 id + ", this->" + an + ", this->" + bn + ", this->" + cn + "," + 303 sp.getEquivalentStressLowerBound(bd) + ");\n"; 304 #endif /* __cplusplus >= 201703L */ 305 } 306 if (r == STRESSANDFLOWCRITERION) { 307 c += "const auto& n" + id + " = dseq" + id + "_ds" + id + ";\n"; 308 c += "const auto& dn" + id + "_ds" + id + " = "; 309 c += "d2seq" + id + "_ds" + id + "ds" + id + ";\n"; 310 } 311 if (r == FLOWCRITERION) { 312 #if __cplusplus >= 201703L 313 c += "const auto [seqf" + id + ", n" + id + ", dn" + id + "_ds" + id + 314 "] = "; 315 c += "computeCazacu2004OrthotropicStressCriterionSecondDerivative(s" + 316 id + ", this->" + an + ", this->" + bn + ",this->" + cn + "," + 317 sp.getEquivalentStressLowerBound(bd) + ");\n"; 318 #else /* __cplusplus >= 201703L */ 319 c += "stress seqf" + id + ";\n"; 320 c += "Stensor n" + id + ";\n"; 321 c += "Stensor4 dn" + id + "_ds" + id + ";\n"; 322 c += 323 "std::tie(seqf" + id + ",n" + id + ",dn" + id + "_ds" + id + ") = "; 324 c += "computeCazacu2004OrthotropicStressCriterionSecondDerivative(s" + 325 id + +", this->" + an + ", this->" + bn + ",this->" + cn + "," + 326 sp.getEquivalentStressLowerBound(bd) + ");\n"; 327 #endif /* __cplusplus >= 201703L */ 328 } 329 return c; 330 } // end of Cazacu2004OrthotropicStressCriterion::computeNormalDerivative 331 isCoupledWithPorosityEvolution() const332 bool Cazacu2004OrthotropicStressCriterion::isCoupledWithPorosityEvolution() const { 333 return false; 334 } // end of Cazacu2004OrthotropicStressCriterion::isCoupledWithPorosityEvolution 335 isNormalDeviatoric() const336 bool Cazacu2004OrthotropicStressCriterion::isNormalDeviatoric() const { 337 return true; 338 } // end of Cazacu2004OrthotropicStressCriterion::isNormalDeviatoric 339 340 StressCriterion::PorosityEffectOnFlowRule getPorosityEffectOnEquivalentPlasticStrain() const341 Cazacu2004OrthotropicStressCriterion::getPorosityEffectOnEquivalentPlasticStrain() const { 342 return StressCriterion::NO_POROSITY_EFFECT_ON_EQUIVALENT_PLASTIC_STRAIN; 343 } // end of Cazacu2004OrthotropicStressCriterion::getPorosityEffectOnEquivalentPlasticStrain() 344 345 Cazacu2004OrthotropicStressCriterion:: 346 ~Cazacu2004OrthotropicStressCriterion() = default; 347 348 } // end of namespace bbrick 349 350 } // end of namespace mfront 351