41 template <
typename DispGradType,
typename BulkType,
typename ShearType>
43 const ShearType& DeltaG)
const
45 constexpr
auto I = Identity<dim>();
46 auto K =
K0 + get<0>(DeltaK);
47 auto G =
G0 + get<0>(DeltaG);
48 auto lambda = K - (2.0 / dim) * G;
49 auto epsilon = 0.5 * (
transpose(du_dX) + du_dX);
50 return lambda *
tr(epsilon) * I + 2.0 * G * epsilon;
85 template <
typename DispGradType,
typename BulkType,
typename ShearType>
87 const ShearType& DeltaG)
const
90 constexpr
auto I = Identity<dim>();
91 auto K =
K0 + get<0>(DeltaK);
92 auto G =
G0 + get<0>(DeltaG);
93 auto lambda = K - (2.0 / dim) * G;
96 auto J = J_minus_1 + 1;
97 return (lambda *
log1p(J_minus_1) * I + G * B_minus_I) / J;
117 template <
typename... T>
119 using type = decltype((T{} + ...));
124 static constexpr
int dim = 3;
125 static constexpr
double tol = 1e-10;
138 template <
typename DisplacementGradient,
typename YieldStrength,
typename SaturationStrength,
typename StrainConstant>
139 auto operator()(
State& state,
const DisplacementGradient du_dX,
const YieldStrength sigma_y,
140 const SaturationStrength sigma_sat,
const StrainConstant strain_constant)
const
152 constexpr
auto I = Identity<dim>();
153 const double K =
E / (3.0 * (1.0 - 2.0 *
nu));
154 const double G = 0.5 *
E / (1.0 +
nu);
159 auto p = K *
tr(el_strain);
160 auto s = 2.0 * G *
dev(el_strain) * one;
165 auto residual = [eqps_old, G, *
this](
auto delta_eqps,
auto trial_mises,
auto y0,
auto ysat,
auto e0) {
166 auto Y = this->
flow_strength(eqps_old + delta_eqps, y0, ysat, e0);
167 return trial_mises - 3.0 * G * delta_eqps - Y;
176 double lower_bound = 0.0;
180 auto [delta_eqps, status] =
181 solve_scalar_equation(residual, 0.0, lower_bound, upper_bound, opts, q, sigma_y, sigma_sat, strain_constant);
183 auto Np = 1.5 * s / q;
185 s = s - 2.0 * G * delta_eqps * Np;
194 template <
typename PlasticStrain,
typename YieldStrength,
typename SaturationStrength,
typename StrainConstant>
195 auto flow_strength(
const PlasticStrain accumulated_plastic_strain,
const YieldStrength sigma_y,
196 const SaturationStrength sigma_sat,
const StrainConstant strain_constant)
const
199 return sigma_sat - (sigma_sat - sigma_y) *
exp(-accumulated_plastic_strain / strain_constant);
#define SERAC_HOST_DEVICE
Macro that evaluates to __host__ __device__ when compiling with nvcc and does nothing on a host compi...
Implementation of the quadrature-function-based functional enabling rapid development of FEM formulat...
SolidMechanics helper data types.
constexpr SERAC_HOST_DEVICE auto tr(const isotropic_tensor< T, m, m > &I)
calculate the trace of an isotropic tensor
constexpr SERAC_HOST_DEVICE auto detApIm1(const tensor< T, 2, 2 > &A)
computes det(A + I) - 1, where precision is not lost when the entries A_{ij} << 1
constexpr SERAC_HOST_DEVICE auto sym(const isotropic_tensor< T, m, m > &I)
return the symmetric part of an isotropic tensor
SERAC_HOST_DEVICE auto sqrt(dual< gradient_type > x)
implementation of square root for dual numbers
constexpr SERAC_HOST_DEVICE auto dev(const tensor< T, n, n > &A)
Calculates the deviator of a matrix (rank-2 tensor)
constexpr SERAC_HOST_DEVICE auto get_value(const T &arg)
return the "value" part from a given type. For non-dual types, this is just the identity function
auto solve_scalar_equation(function &&f, double x0, double lower_bound, double upper_bound, ScalarSolverOptions options, ParamTypes... params)
Solves a nonlinear scalar-valued equation and gives derivatives of solution to parameters.
constexpr SERAC_HOST_DEVICE auto norm(const isotropic_tensor< T, m, m > &I)
compute the Frobenius norm (sqrt(tr(dot(transpose(I), I)))) of an isotropic tensor
constexpr SERAC_HOST_DEVICE auto transpose(const isotropic_tensor< T, m, m > &I)
return the transpose of an isotropic tensor
SERAC_HOST_DEVICE auto log1p(dual< gradient_type > a)
implementation of the natural logarithm of one plus the argument function for dual numbers
SERAC_HOST_DEVICE auto exp(dual< gradient_type > a)
implementation of exponential function for dual numbers
The material and load types for the solid functional physics module.
see Nothing for a complete description of this class and when to use it
Settings for solve_scalar_equation.
double xtol
absolute tolerance on Newton correction
variables required to characterize the hysteresis response
double accumulated_plastic_strain
uniaxial equivalent plastic strain
tensor< double, dim, dim > plastic_strain
plastic strain
J2 material with Voce hardening, with hardening parameters exposed as differentiable parameters.
auto flow_strength(const PlasticStrain accumulated_plastic_strain, const YieldStrength sigma_y, const SaturationStrength sigma_sat, const StrainConstant strain_constant) const
Computes flow strength from Voce's hardening law.
double density
mass density
double nu
Poisson's ratio.
static constexpr double tol
relative tolerance on residual for accepting return map solution
static constexpr int dim
dimension of space
auto operator()(State &state, const DisplacementGradient du_dX, const YieldStrength sigma_y, const SaturationStrength sigma_sat, const StrainConstant strain_constant) const
calculate the Cauchy stress, given the displacement gradient and previous material state
Linear isotropic elasticity material model.
double K0
base bulk modulus
SERAC_HOST_DEVICE auto operator()(State &, const DispGradType &du_dX, const BulkType &DeltaK, const ShearType &DeltaG) const
stress calculation for a linear isotropic material model
double G0
base shear modulus
double density
mass density
static constexpr int numParameters()
The number of parameters in the model.
Neo-Hookean material model.
double G0
base shear modulus
SERAC_HOST_DEVICE auto operator()(State &, const DispGradType &du_dX, const BulkType &DeltaK, const ShearType &DeltaG) const
stress calculation for a NeoHookean material model
double K0
base bulk modulus
double density
mass density
static constexpr int numParameters()
The number of parameters in the model.
Infers type resulting from algebraic expressions of a group of variables.
decltype((T{}+...)) type
type of the sum of the parameters