31 static constexpr
int dim = 3;
52 double order_parameter,
double transition_temperature,
double N_b_squared)
61 SLIC_ERROR_ROOT_IF(
density <= 0.0,
"Density must be positive in the LCE material model.");
62 SLIC_ERROR_ROOT_IF(
shear_modulus_ <= 0.0,
"Shear modulus must be positive in the LCE material model.");
63 SLIC_ERROR_ROOT_IF(
bulk_modulus_ <= 0.0,
"Bulk modulus must be positive in the LCE material model.");
64 SLIC_ERROR_ROOT_IF(
order_constant_ <= 0.0,
"Order constant must be positive in the LCE material model.");
66 "The transition temperature must be positive in the LCE material model.");
82 template <
typename DispGradType,
typename OrderParamType,
typename GammaAngleType>
84 OrderParamType temperature_tuple, GammaAngleType gamma_tuple)
const
89 auto temperature = get<0>(temperature_tuple);
90 auto gamma = get<0>(gamma_tuple);
92 auto I = Identity<dim>();
103 auto F = displacement_grad + I;
105 auto F_hat =
dot(F,
inv(F_old));
110 auto mu0 = (
N_b_squared_ / 3.0) * ((1 - q0) * I) + (3 * q0 *
outer(normal, normal));
145 template <
typename S,
typename T,
typename U>
156 constexpr
auto I = Identity<dim>();
157 auto n_dyad =
outer(normal, normal);
162 auto Q_old = 0.5 * ((1.0 - q_old) * I + 3.0 * q_old * n_dyad);
163 auto Q = 0.5 * ((1.0 - q) * I + 3.0 * q * n_dyad);
167 auto R_hat = F_hat *
inv(U_hat);
172 auto mu_hat = mu_old + alpha * (Q - Q_old);
214 static constexpr
int dim = 3;
226 double beta_parameter)
233 SLIC_ERROR_ROOT_IF(
density <= 0.0,
"Density must be positive in the LCE material model.");
234 SLIC_ERROR_ROOT_IF(
young_modulus_ <= 0.0,
"Bulk modulus must be positive in the LCE material model.");
235 SLIC_ERROR_ROOT_IF(
poisson_ratio_ <= 0.0,
"Poisson ratio must be positive in the LCE material model.");
237 "Initial order parameter must be positive in the LCE material model.");
238 SLIC_ERROR_ROOT_IF(
beta_parameter_ <= 0.0,
"The beta parameter must be positive in the LCE material model.");
257 template <
typename DispGradType,
typename OrderParamType,
typename GammaAngleType,
typename EtaAngleType>
259 OrderParamType inst_order_param_tuple, GammaAngleType gamma_tuple,
260 EtaAngleType eta_tuple)
const
265 auto gamma = get<0>(gamma_tuple);
266 auto eta = get<0>(eta_tuple);
270 auto St = get<0>(inst_order_param_tuple);
277 auto I = Identity<dim>();
278 auto F = displacement_grad + I;
280 auto E = 0.5 * (C - I);
284 auto S_stress_1 = lambda *
tr(E) * I;
285 auto S_stress_2 = 2 * mu * E;
289 auto stress = 1.0 / J * F * (S_stress_1 + S_stress_2 + S_stress_3) *
transpose(F);
311 template <
typename DispGradType,
typename orderParamType,
typename GammaAngleType,
typename EtaAngleType>
313 orderParamType inst_order_param_tuple, GammaAngleType gamma_tuple,
314 EtaAngleType eta_tuple)
const
319 auto gamma = get<0>(gamma_tuple);
320 auto eta = get<0>(eta_tuple);
324 auto St = get<0>(inst_order_param_tuple);
331 auto I = Identity<dim>();
332 auto F = displacement_grad + I;
334 auto E = 0.5 * (C - I);
337 auto strain_energy_1 = 0.5 * lambda *
tr(E) *
tr(E);
338 auto strain_energy_2 = mu *
inner(E, E);
341 auto strain_energy = strain_energy_1 + strain_energy_2 + strain_energy_3;
343 return strain_energy;
#define SERAC_HOST_DEVICE
Macro that evaluates to __host__ __device__ when compiling with nvcc and does nothing on a host compi...
This file contains the declaration of a dual number class.
Accelerator functionality.
SERAC_HOST_DEVICE auto sin(dual< gradient_type > a)
implementation of sine for dual numbers
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 dot(const isotropic_tensor< S, m, m > &I, const tensor< T, m, n... > &A)
dot product between an isotropic and (nonisotropic) tensor
SERAC_HOST_DEVICE auto cos(dual< gradient_type > a)
implementation of cosine for dual numbers
constexpr SERAC_HOST_DEVICE auto inner(const tensor< S, m, n > &A, const tensor< T, m, n > &B)
this function contracts over all indices of the two tensor arguments
auto matrix_sqrt(const tensor< T, dim, dim > &A)
compute the matrix square root of a square, real-valued, symmetric matrix i.e. given A,...
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
constexpr SERAC_HOST_DEVICE auto det(const isotropic_tensor< T, m, m > &I)
compute the determinant of an isotropic tensor
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 inv(const isotropic_tensor< T, m, m > &I)
return the inverse 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 log(dual< gradient_type > a)
implementation of the natural logarithm function for dual numbers
constexpr SERAC_HOST_DEVICE auto outer(double A, tensor< T, n > B)
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
internal variables for the liquid crystal elastomer model
tensor< double, dim, dim > distribution_tensor
mu from the last timestep
tensor< double, dim, dim > deformation_gradient
F from the last timestep.
double temperature
temperature at the last timestep
Brighenti's liquid crystal elastomer model.
double bulk_modulus_
bulk modulus in stress-free configuration
static constexpr int dim
this model is only intended to be used in 3D
double shear_modulus_
shear modulus in stress-free configuration
auto calculateDistributionTensor(const State &state, const tensor< S, dim, dim > &F_hat, const T theta, const tensor< U, dim > &normal) const
Compute the distribution tensor using Brighenti's model.
double transition_temperature_
Transition temperature.
double N_b_squared_
Kuhn segment parameters.
double initial_order_parameter_
initial value of order parameter
SERAC_HOST_DEVICE auto operator()(State &state, const tensor< DispGradType, dim, dim > &displacement_grad, OrderParamType temperature_tuple, GammaAngleType gamma_tuple) const
Material response.
double order_constant_
Order constant.
double density
mass density
LiquidCrystElastomerBrighenti(double rho, double shear_modulus, double bulk_modulus, double order_constant, double order_parameter, double transition_temperature, double N_b_squared)
Constructor.
Bertoldi's liquid crystal elastomer model Paper: Li, S., Librandi, G., Yao, Y., Richard,...
auto calculateStrainEnergy(const State &, const tensor< DispGradType, dim, dim > &displacement_grad, orderParamType inst_order_param_tuple, GammaAngleType gamma_tuple, EtaAngleType eta_tuple) const
Compute the strain energy.
SERAC_HOST_DEVICE auto operator()(State &, const tensor< DispGradType, dim, dim > &displacement_grad, OrderParamType inst_order_param_tuple, GammaAngleType gamma_tuple, EtaAngleType eta_tuple) const
Material response.
double initial_order_parameter_
initial value of order parameter
double poisson_ratio_
poisson's ratio
static constexpr int dim
this model is only intended to be used in 3D
double density
mass density
double beta_parameter_
Degree of coupling between elastic and nematic energies.
double young_modulus_
Young's modulus in stress-free configuration.
LiquidCrystalElastomerBertoldi(double rho, double young_modulus, double poisson_ratio, double initial_order_parameter, double beta_parameter)
Constructor.
Implementation of the tensor class used by Functional.
Implements a std::tuple-like object that works in CUDA kernels.