v0.13.0
SmallStrainPlasticityMaterialModels.hpp
/** \file SmallStrainPlasticityMaterialModels.hpp
* \ingroup small_strain_plasticity
* \brief Small Strain plasticity material models
* \example SmallStrainPlasticityMaterialModels.hpp
*/
/*
* This file is part of MoFEM.
* MoFEM is free software: you can redistribute it and/or modify it under
* the terms of the GNU Lesser General Public License as published by the
* Free Software Foundation, either version 3 of the License, or (at your
* option) any later version.
*
* MoFEM is distributed in the hope that it will be useful, but WITHOUT
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
* License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
#ifndef __SMALL_STRAIN_PLASTICITY_MATERIAL_MODELS_HPP
#define __SMALL_STRAIN_PLASTICITY_MATERIAL_MODELS_HPP
#ifndef WITH_ADOL_C
#error "MoFEM need to be compiled with ADOL-C"
#endif
/** \brief J2 plasticity (Kinematic Isotropic (Linear) Hardening)
* \ingroup small_strain_plasticity
*/
internalVariables.resize(7,false);
}
double mu; //< Lamé parameter
double lambda; //< Lamé parameters
double H; //< Isotropic hardening
double K; //< Kinematic hardening
double phi; //< Combined isotropic/kinematic hardening
double sIgma_y; //< Yield stress
/** \brief Available energy
\f[
\psi =
\frac{1}{2} \lambda \textrm{tr}[\varepsilon]^2 + \mu \varepsilon : \varepsilon
+
\sigma_y\alpha
+
\frac{1}{2} \phi H \alpha^2
+
\frac{1}{2} (1-\phi)K \beta^2
\f]
Isotropic hardening variable \f$\alpha\f$ is in first index of
internal variables vector. Kinematic hardening variable is in
the remaining indices of internal variables vector.
*/
virtual PetscErrorCode freeHemholtzEnergy() {
PetscFunctionBegin;
//elastic part
adouble tR = 0;
for(int dd = 0;dd<3;dd++) {
}
a_w = 0.5*lambda*tR*tR;
adouble t;
int dd = 0;
for(;dd<3;dd++) {
a_w += mu*t*t;
}
for(;dd<6;dd++) {
a_w += 0.5*mu*t*t;
}
//plastic part
//isotropic
//kinematic
for(unsigned int dd = 1;dd<a_internalVariables.size();dd++) {
}
PetscFunctionReturn(0);
}
/** \brief Auxiliary function
\f[
\eta=\textrm{dev}[\sigma]-\overline{\beta}
\f]
\f[
f = \frac{3}{2} \eta:\eta
\f]
This is \f$3J_2\f$.
*/
adouble t,f;
PetscErrorCode evalF() {
PetscFunctionBegin;
adouble tR = 0;
for(int dd = 0;dd<3;dd++) {
tR += a_sTress[dd];
}
// Fix the constant to mach uniaxial test
const double c = 3./2.;
tR /= 3.;
f = 0;
for(int dd = 0;dd<3;dd++) {
f += c*t*t;
}
for(int dd = 3;dd<6;dd++) {
f += c*2.*t*t; // s:s off diagonal terms need to go twice
}
PetscFunctionReturn(0);
}
/** \brief Evaluate yield function
\f[
y = \sqrt{f} - \overline{\alpha}
\f]
where \f$f\f$ is defined in \ref evalF.
*/
virtual PetscErrorCode yieldFunction() {
PetscFunctionBegin;
ierr = evalF(); CHKERRQ(ierr);
a_y = sqrt(f) - a_internalFluxes[0];
PetscFunctionReturn(0);
}
/** \brief Associated flow potential
See \ref yieldFunction.
*/
virtual PetscErrorCode flowPotential() {
PetscFunctionBegin;
ierr = evalF(); CHKERRQ(ierr);
a_h = sqrt(f) - a_internalFluxes[0];
PetscFunctionReturn(0);
}
};
/** \brief Small strain plasticity with paraboloidal yield criterion (Isotropic Hardening)
* \ingroup small_strain_plasticity
*/
internalVariables.resize(2,false);
}
double mu; //< Lamé parameter
double lambda; //< Lamé parameters
double nup; //< Plastic Poisson’s ratio
double Ht, Hc; //< Isotropic hardening for tension and compression
double sIgma_yt, sIgma_yc; //< Yield stress in tension and compression
double nt, nc; //< Hardening parameters for tension and compression
/** \brief Available energy
\f[
\psi =
\frac{1}{2} \lambda \textrm{tr}[\varepsilon]^2 + \mu \varepsilon : \varepsilon
+
\sigma_{yt}\alpha_0
+
\frac{1}{2} H_t \alpha_0^2
+
\sigma_{yc}\alpha_1
+
\frac{1}{2} H_c \alpha_1^2
\f]
//Energy with linear hardening
// */
// virtual PetscErrorCode freeHemholtzEnergy() {
// PetscFunctionBegin;
// //elastic part
// adouble tR = 0;
// for(int dd = 0;dd<3;dd++) {
// tR += a_sTrain[dd]-a_plasticStrain[dd];
// }
// a_w = 0.5*lambda*tR*tR;
// adouble t;
// int dd = 0;
// for(;dd<3;dd++) {
// t = a_sTrain[dd]-a_plasticStrain[dd];
// a_w += mu*t*t;
// }
// for(;dd<6;dd++) {
// t = a_sTrain[dd]-a_plasticStrain[dd];
// a_w += 0.5*mu*t*t;
// }
// //plastic part
// //isotropic
// a_w += a_internalVariables[0]*sIgma_yt + 0.5*Ht*a_internalVariables[0]*a_internalVariables[0];
// a_w += a_internalVariables[1]*sIgma_yc + 0.5*Hc*a_internalVariables[1]*a_internalVariables[1];
// PetscFunctionReturn(0);
// }
//Energy with exponential hardening
// */
virtual PetscErrorCode freeHemholtzEnergy() {
PetscFunctionBegin;
//elastic part
adouble tR = 0;
for(int dd = 0;dd<3;dd++) {
}
a_w = 0.5*lambda*tR*tR;
adouble t;
int dd = 0;
for(;dd<3;dd++) {
a_w += mu*t*t;
}
for(;dd<6;dd++) {
a_w += 0.5*mu*t*t;
}
//plastic part
//isotropic
// a_w += a_internalVariables[0]*sIgma_yt + 0.5*Ht*a_internalVariables[0]*a_internalVariables[0];
// a_w += a_internalVariables[1]*sIgma_yc + 0.5*Hc*a_internalVariables[1]*a_internalVariables[1];
PetscFunctionReturn(0);
}
/** \brief Auxiliary function
\f[
I_1 = \textrm{tr} (\boldsymbol{\sigma})
\f]
\f[
\eta=\textrm{dev}[\boldsymbol{\sigma}]
\f]
\f[
J_2 = \frac{1}{2} \eta:\eta
\f]
*/
adouble I1,J2;
PetscErrorCode evalF() {
PetscFunctionBegin;
adouble tR = 0;
for(int dd = 0;dd<3;dd++) {
tR += a_sTress[dd];
}
I1=tR;
tR /= 3.;
J2 = 0;
adouble t;
for(int dd = 0;dd<3;dd++) {
t = (a_sTress[dd]-tR);
J2 += t*t;
}
for(int dd = 3;dd<6;dd++) {
J2 += 2.*a_sTress[dd]*a_sTress[dd]; // s:s off diagonal terms need to go twice
}
J2=0.5*J2;
PetscFunctionReturn(0);
}
/** \brief Evaluate yield function
\f[
y = 6J_2 + 2I_1\left(\overline{\alpha_1}-\overline{\alpha_0}\right) - 2\overline{\alpha_0} \,\overline{\alpha_1}
\f]
where
\f[
\overline{\alpha_0}=\frac{\partial \psi}{\partial \alpha_0}=\sigma_{yt} + H_t \alpha_0
\f]
\f[
\overline{\alpha_1}=\frac{\partial \psi}{\partial \alpha_1}=\sigma_{yc} + H_c \alpha_1
\f]
*/
virtual PetscErrorCode yieldFunction() {
PetscFunctionBegin;
ierr = evalF(); CHKERRQ(ierr);
PetscFunctionReturn(0);
}
/** \brief Flow potential
\f[
\Psi = 6J_2 + 2\alpha I_1 \left(\overline{\alpha_1}-\overline{\alpha_0}\right) - 2\overline{\alpha_0} \,\overline{\alpha_1}
\f]
\f[
\alpha= \frac{1-2\nu_p}{1+\nu_p}
\f]
*/
virtual PetscErrorCode flowPotential() {
PetscFunctionBegin;
ierr = evalF(); CHKERRQ(ierr);
double alpha= (1-2*nup)/(1+nup); //relation between alpha and plastic Poission's ratio
PetscFunctionReturn(0);
}
};
#endif //__SMALL_STRAIN_PLASTICITY_MATERIAL_MODELS_HPP
constexpr double alpha
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
constexpr double t
plate stiffness
Definition: plate.cpp:76
J2 plasticity (Kinematic Isotropic (Linear) Hardening)
virtual PetscErrorCode freeHemholtzEnergy()
Available energy.
virtual PetscErrorCode flowPotential()
Associated flow potential.
virtual PetscErrorCode yieldFunction()
Evaluate yield function.
Small strain plasticity with paraboloidal yield criterion (Isotropic Hardening)
virtual PetscErrorCode yieldFunction()
Evaluate yield function.
virtual PetscErrorCode flowPotential()
Flow potential.
virtual PetscErrorCode freeHemholtzEnergy()
Available energy.
Small strain finite element implementation.