MaterialLib::Solids::Phasefield Namespace Reference

Classes
class	AT_DegradationDerivative
class	COHESIVE_DegradationDerivative
class	DegradationDerivative

Enumerations
enum class	PhaseFieldModel { AT1 , AT2 , COHESIVE }
enum class	SofteningCurve { Linear , Exponential }
enum class	EnergySplitModel {   Isotropic , VolDev , EffectiveStress , OrthoVolDev ,   OrthoMasonry }

Functions
double	heaviside (double const v)
double	macaulayTensile (double const v)
double	macaulayCompressive (double v)
template<int DisplacementDim>
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > >	calculateVolDevDegradedStress (double const degradation, double const bulk_modulus, double const mu, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps)
template<int DisplacementDim>
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > >	calculateIsotropicDegradedStress (double const degradation, double const bulk_modulus, double const mu, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps)
template<int DisplacementDim>
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > >	calculateIsotropicDegradedStressWithRankineEnergy (double const degradation, double const bulk_modulus, double const mu, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps)
template<int DisplacementDim>
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > >	calculateOrthoVolDevDegradedStress (double const degradation, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > const &C_ortho)
template<int DisplacementDim>
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > >	calculateOrthoMasonryDegradedStress (double const degradation, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > const &C_ortho)
template<int DisplacementDim>
PhaseFieldModel	convertStringToPhaseFieldModel (std::string const &phasefield_model)
template<int DisplacementDim>
SofteningCurve	convertStringToSofteningCurve (std::optional< std::string > const &softening_curve)
template<int DisplacementDim>
EnergySplitModel	convertStringToEnergySplitModel (std::string const &energy_split_model)
template<int DisplacementDim>
std::unique_ptr< DegradationDerivative >	creatDegradationDerivative (PhaseFieldModel const &phasefield_model, double const lch, SofteningCurve const &softening_curve)
template<typename T_DNDX , typename T_N , typename T_W , typename T_D , typename T_LOCAL_JAC , typename T_LOCAL_RHS >
void	calculateCrackLocalJacobianAndResidual (T_DNDX &dNdx, T_N &N, T_W &w, T_D &d, T_LOCAL_JAC &local_Jac, T_LOCAL_RHS local_rhs, double const gc, double const ls, PhaseFieldModel &phasefield_model)
template<typename T_VECTOR , typename T_MATRIX , int DisplacementDim>
void	calculateStress (T_VECTOR &sigma, T_VECTOR &sigma_tensile, T_VECTOR &sigma_compressive, T_VECTOR &eps_tensile, T_MATRIX &D, T_MATRIX &C_tensile, T_MATRIX &C_compressive, double &strain_energy_tensile, double &elastic_energy, double const degradation, T_VECTOR const &eps, EnergySplitModel const &energy_split_model, double const t, ParameterLib::SpatialPosition const &x, MaterialLib::Solids::MechanicsBase< DisplacementDim > const &solid_material)

Enumeration Type Documentation

EnergySplitModel

enum class MaterialLib::Solids::Phasefield::EnergySplitModel

strong

Enumerator
Isotropic
VolDev
EffectiveStress
OrthoVolDev
OrthoMasonry

Definition at line 85 of file PhaseFieldBase.h.

{
    Isotropic,
    VolDev,
    EffectiveStress,
    OrthoVolDev,
    OrthoMasonry
};

PhaseFieldModel

enum class MaterialLib::Solids::Phasefield::PhaseFieldModel

strong

Enumerator
AT1
AT2
COHESIVE

Definition at line 28 of file PhaseFieldBase.h.

{
    AT1,
    AT2,
    COHESIVE
};

SofteningCurve

enum class MaterialLib::Solids::Phasefield::SofteningCurve

strong

Enumerator
Linear
Exponential

Definition at line 57 of file PhaseFieldBase.h.

{
    Linear,
    Exponential
};

Function Documentation

calculateCrackLocalJacobianAndResidual()

template<typename T_DNDX , typename T_N , typename T_W , typename T_D , typename T_LOCAL_JAC , typename T_LOCAL_RHS >

void MaterialLib::Solids::Phasefield::calculateCrackLocalJacobianAndResidual	(	T_DNDX &	dNdx,
		T_N &	N,
		T_W &	w,
		T_D &	d,
		T_LOCAL_JAC &	local_Jac,
		T_LOCAL_RHS	local_rhs,
		double const	gc,
		double const	ls,
		PhaseFieldModel &	phasefield_model )

Definition at line 292 of file PhaseFieldBase.h.

{
    switch (phasefield_model)
    {
        case PhaseFieldModel::AT1:
        {
            auto const local_Jac_AT1 =
                (gc * 0.75 * dNdx.transpose() * dNdx * ls * w).eval();
            local_Jac.noalias() += local_Jac_AT1;
 
            local_rhs.noalias() -=
                gc * (-0.375 * N.transpose() / ls) * w + local_Jac_AT1 * d;
            break;
        }
        case PhaseFieldModel::AT2:
        {
            auto const local_Jac_AT2 =
                (gc * (N.transpose() * N / ls + dNdx.transpose() * dNdx * ls) *
                 w)
                    .eval();
            local_Jac.noalias() += local_Jac_AT2;
 
            local_rhs.noalias() -=
                local_Jac_AT2 * d - gc * (N.transpose() / ls) * w;
            break;
        }
        case PhaseFieldModel::COHESIVE:
        {
            auto const local_Jac_COHESIVE =
                (2.0 / std::numbers::pi * gc *
                 (-N.transpose() * N / ls + dNdx.transpose() * dNdx * ls) * w)
                    .eval();
 
            local_Jac.noalias() += local_Jac_COHESIVE;
 
            local_rhs.noalias() -= local_Jac_COHESIVE * d;
            break;
        }
        default:
        {
            OGS_FATAL("Invalid phase field model specified.");
        }
    }
};

References AT1, AT2, COHESIVE, and OGS_FATAL.

calculateIsotropicDegradedStress()

template<int DisplacementDim>

std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Phasefield::calculateIsotropicDegradedStress	(	double const	degradation,
		double const	bulk_modulus,
		double const	mu,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps )

Definition at line 138 of file LinearElasticIsotropicPhaseField.h.

{
    static int const KelvinVectorSize =
        MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim);
    using KelvinVector =
        MathLib::KelvinVector::KelvinVectorType<DisplacementDim>;
    using KelvinMatrix =
        MathLib::KelvinVector::KelvinMatrixType<DisplacementDim>;
    using Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>;
    // calculation of deviatoric parts
    auto const& P_dev = Invariants::deviatoric_projection;
    KelvinVector const epsd_curr = P_dev * eps;
 
    // Hydrostatic part for the stress and the tangent.
    double const eps_curr_trace = Invariants::trace(eps);
 
    KelvinMatrix C_tensile = KelvinMatrix::Zero();
    KelvinMatrix C_compressive = KelvinMatrix::Zero();
 
    double const strain_energy_tensile =
        bulk_modulus / 2 * eps_curr_trace * eps_curr_trace +
        mu * epsd_curr.transpose() * epsd_curr;
    double const elastic_energy = degradation * strain_energy_tensile;
    KelvinVector const sigma_tensile =
        bulk_modulus * eps_curr_trace * Invariants::identity2 +
        2 * mu * epsd_curr;
    C_tensile.template topLeftCorner<3, 3>().setConstant(bulk_modulus);
    C_tensile.noalias() += 2 * mu * P_dev * KelvinMatrix::Identity();
    KelvinVector const sigma_real = degradation * sigma_tensile;
 
    KelvinMatrix const D = degradation * C_tensile + C_compressive;
 
    return std::make_tuple(sigma_real, sigma_tensile, D, strain_energy_tensile,
                           elastic_energy, C_tensile, C_compressive);
}

References MathLib::KelvinVector::kelvin_vector_dimensions().

Referenced by calculateStress().

calculateIsotropicDegradedStressWithRankineEnergy()

template<int DisplacementDim>

std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Phasefield::calculateIsotropicDegradedStressWithRankineEnergy	(	double const	degradation,
		double const	bulk_modulus,
		double const	mu,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps )

Definition at line 191 of file LinearElasticIsotropicPhaseField.h.

{
    static int const KelvinVectorSize =
        MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim);
    using KelvinVector =
        MathLib::KelvinVector::KelvinVectorType<DisplacementDim>;
    using KelvinMatrix =
        MathLib::KelvinVector::KelvinMatrixType<DisplacementDim>;
    using Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>;
    // calculation of deviatoric parts
    auto const& P_dev = Invariants::deviatoric_projection;
    KelvinVector const epsd_curr = P_dev * eps;
 
    // Hydrostatic part for the stress and the tangent.
    double const eps_curr_trace = Invariants::trace(eps);
 
    KelvinMatrix C_tensile = KelvinMatrix::Zero();
    KelvinMatrix C_compressive = KelvinMatrix::Zero();
 
    KelvinVector const sigma_tensile =
        bulk_modulus * eps_curr_trace * Invariants::identity2 +
        2 * mu * epsd_curr;
    // The maximum principal stress
    auto sigma_tensor =
        MathLib::KelvinVector::kelvinVectorToTensor(sigma_tensile);
    Eigen::SelfAdjointEigenSolver<decltype(sigma_tensor)>
        selfAdjoint_eigen_solver(sigma_tensor);
    double const max_principle_stress =
        selfAdjoint_eigen_solver
            .eigenvalues()[2];  // selfAdjoint_eigen_solver function gives an
                                // ascending sequence.
    //
    double const E0 = 9.0 * mu * bulk_modulus / (3.0 * bulk_modulus + mu);
    auto hs = [&](double const v) { return heaviside(v); };
    //
    double const strain_energy_tensile = 0.5 * max_principle_stress *
                                         max_principle_stress *
                                         hs(max_principle_stress) / E0;
    double const elastic_energy =
        degradation * (bulk_modulus / 2 * eps_curr_trace * eps_curr_trace +
                       mu * epsd_curr.transpose() * epsd_curr);
    //
    C_tensile.template topLeftCorner<3, 3>().setConstant(bulk_modulus);
    C_tensile.noalias() += 2 * mu * P_dev * KelvinMatrix::Identity();
    KelvinVector const sigma_real = degradation * sigma_tensile;
 
    KelvinMatrix const D = degradation * C_tensile + C_compressive;
    return std::make_tuple(sigma_real, sigma_tensile, D, strain_energy_tensile,
                           elastic_energy, C_tensile, C_compressive);
}

References heaviside(), MathLib::KelvinVector::kelvin_vector_dimensions(), and MathLib::KelvinVector::kelvinVectorToTensor().

calculateOrthoMasonryDegradedStress()

template<int DisplacementDim>

std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Phasefield::calculateOrthoMasonryDegradedStress	(	double const	degradation,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps,
		MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > const &	C_ortho )

Definition at line 127 of file LinearElasticOrthotropicPhaseField.h.

{
    using KelvinVector =
        MathLib::KelvinVector::KelvinVectorType<DisplacementDim>;
    using KelvinMatrix =
        MathLib::KelvinVector::KelvinMatrixType<DisplacementDim>;
 
    KelvinMatrix const C = C_ortho;
 
    // calculation of square root of elasticity tensor
    Eigen::SelfAdjointEigenSolver<KelvinMatrix> es(C);
    KelvinMatrix const sqrt_C = es.operatorSqrt();
    Eigen::SelfAdjointEigenSolver<KelvinMatrix> es_inverse(C.inverse());
    KelvinMatrix const sqrt_C_inverse = es_inverse.operatorSqrt();
 
    // strain in a transformed space
    KelvinVector const epst = sqrt_C * eps;
 
    // strain tensor in 3D
    Eigen::Matrix3d const eps_3D =
        MathLib::KelvinVector::kelvinVectorToTensor(epst);
 
    Eigen::SelfAdjointEigenSolver<Eigen::Matrix3d> const es_eps_3D(eps_3D);
    Eigen::Vector3d const eigen_values_eps_3D = es_eps_3D.eigenvalues().real();
    Eigen::Matrix3d const eigen_vectors_eps_3D = es_eps_3D.eigenvectors();
    Eigen::Matrix3d const E1_eigenp =
        eigen_vectors_eps_3D.col(0) * eigen_vectors_eps_3D.col(0).transpose();
 
    Eigen::Matrix3d const E3_eigenp =
        eigen_vectors_eps_3D.col(2) * eigen_vectors_eps_3D.col(2).transpose();
 
    KelvinVector const E1_vec =
        MathLib::KelvinVector::tensorToKelvin<DisplacementDim>(E1_eigenp);
 
    KelvinVector const E3_vec =
        MathLib::KelvinVector::tensorToKelvin<DisplacementDim>(E3_eigenp);
 
    KelvinMatrix const E1oE1 = E1_vec * E1_vec.transpose();
    KelvinMatrix const E3oE3 = E3_vec * E3_vec.transpose();
 
    KelvinMatrix I_p = KelvinMatrix::Zero();
    KelvinMatrix I_n = KelvinMatrix::Zero();
 
    KelvinMatrix const I_S = KelvinMatrix::Identity();
    if (DisplacementDim == 2)
    {
        if (std::abs(eigen_values_eps_3D(2) - eigen_values_eps_3D(0)) >
            std::numeric_limits<double>::epsilon())
        {
            I_p = (macaulayTensile(eigen_values_eps_3D(0)) -
                   macaulayTensile(eigen_values_eps_3D(2))) /
                      (eigen_values_eps_3D(0) - eigen_values_eps_3D(2)) *
                      (I_S - (E1oE1 + E3oE3)) +
                  (heaviside(eigen_values_eps_3D(0)) * E1oE1 +
                   heaviside(eigen_values_eps_3D(2)) * E3oE3);
            I_n = (macaulayCompressive(eigen_values_eps_3D(0)) -
                   macaulayCompressive(eigen_values_eps_3D(2))) /
                      (eigen_values_eps_3D(0) - eigen_values_eps_3D(2)) *
                      (I_S - (E1oE1 + E3oE3)) +
                  (heaviside(-eigen_values_eps_3D(0)) * E1oE1 +
                   heaviside(-eigen_values_eps_3D(2)) * E3oE3);
        }
        else
        {
            I_p = heaviside(eigen_values_eps_3D(0)) * I_S;
            I_n = heaviside(-eigen_values_eps_3D(0)) * I_S;
        }
    }
 
    // strain tensile / compressive
    KelvinVector const eps_tensile = sqrt_C_inverse * (I_p * sqrt_C) * eps;
    KelvinVector const eps_compressive = sqrt_C_inverse * (I_n * sqrt_C) * eps;
 
    // C tensile / compressive
    KelvinMatrix const C_tensile = sqrt_C * I_p * sqrt_C;
    KelvinMatrix const C_compressive = sqrt_C * I_n * sqrt_C;
 
    // stress tensile / compressive
    KelvinVector const sigma_tensile = C * eps_tensile;
    KelvinVector const sigma_compressive = C * eps_compressive;
 
    // decomposition of strain energy density
    double const strain_energy_tensile =
        0.5 * sigma_tensile.adjoint() * eps_tensile;
    double const strain_energy_compressive =
        0.5 * sigma_compressive.adjoint() * eps_compressive;
 
    double const elastic_energy =
        degradation * strain_energy_tensile + strain_energy_compressive;
 
    KelvinVector const sigma_real =
        degradation * sigma_tensile + sigma_compressive;
 
    KelvinMatrix const D = degradation * C_tensile + C_compressive;
 
    return std::make_tuple(eps_tensile, sigma_real, sigma_tensile, D,
                           strain_energy_tensile, elastic_energy, C_tensile,
                           C_compressive);
}

References heaviside(), MathLib::KelvinVector::kelvinVectorToTensor(), macaulayCompressive(), macaulayTensile(), and MathLib::KelvinVector::tensorToKelvin().

Referenced by calculateStress().

calculateOrthoVolDevDegradedStress()

template<int DisplacementDim>

std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Phasefield::calculateOrthoVolDevDegradedStress	(	double const	degradation,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps,
		MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > const &	C_ortho )

Definition at line 38 of file LinearElasticOrthotropicPhaseField.h.

{
    static constexpr int KelvinVectorSize =
        MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim);
    using KelvinVector =
        MathLib::KelvinVector::KelvinVectorType<DisplacementDim>;
    using KelvinMatrix =
        MathLib::KelvinVector::KelvinMatrixType<DisplacementDim>;
    using Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>;
    // calculation of deviatoric parts
    auto const& P_dev = Invariants::deviatoric_projection;
 
    KelvinMatrix const C = C_ortho;
    // calculation of square root of elasticity tensor
    Eigen::SelfAdjointEigenSolver<KelvinMatrix> es(C);
    KelvinMatrix const sqrt_C = es.operatorSqrt();
    Eigen::SelfAdjointEigenSolver<KelvinMatrix> es_inverse(C.inverse());
    KelvinMatrix const sqrt_C_inverse = es_inverse.operatorSqrt();
 
    // strain in a transformed space
    KelvinVector const epst = sqrt_C * eps;
    double const epst_curr_trace = Invariants::trace(epst);
 
    // projection tensors in transformed space
    KelvinMatrix teps_p = KelvinMatrix::Zero();
    KelvinMatrix teps_n = KelvinMatrix::Zero();
    if (epst_curr_trace >= 0) /* QQQ */
    {
        teps_p.template topLeftCorner<3, 3>().setConstant(1. / 3);
    }
    else
    {
        teps_n.template topLeftCorner<3, 3>().setConstant(1. / 3);
    }
 
    teps_p.noalias() += P_dev * KelvinMatrix::Identity();
 
    // strain tensile / compressive
    KelvinVector const eps_tensile = sqrt_C_inverse * (teps_p * sqrt_C) * eps;
    KelvinVector const eps_compressive =
        sqrt_C_inverse * (teps_n * sqrt_C) * eps;
 
    // projection tensors in original space
    KelvinMatrix const der_eps_p = sqrt_C_inverse * (teps_p * sqrt_C);
    KelvinMatrix const der_eps_n = sqrt_C_inverse * (teps_n * sqrt_C);
    // C tensile / compressive
    KelvinMatrix const C_tensile = der_eps_p.transpose() * C * der_eps_p;
    KelvinMatrix const C_compressive = der_eps_n.transpose() * C * der_eps_n;
 
    // stress tensile / compressive
    KelvinVector const sigma_tensile = C_tensile * eps;
    KelvinVector const sigma_compressive = C_compressive * eps;
 
    // decomposition of strain energy density
    double const strain_energy_tensile =
        0.5 * sigma_tensile.adjoint() * eps_tensile;
    double const strain_energy_compressive =
        0.5 * sigma_compressive.adjoint() * eps_compressive;
 
    double const elastic_energy =
        degradation * strain_energy_tensile + strain_energy_compressive;
 
    KelvinVector const sigma_real =
        degradation * sigma_tensile + sigma_compressive;
    KelvinMatrix const D = degradation * C_tensile + C_compressive;
 
    return std::make_tuple(eps_tensile, sigma_real, sigma_tensile,
                           sigma_compressive, D, strain_energy_tensile,
                           elastic_energy, C_tensile, C_compressive);
}

References MathLib::KelvinVector::kelvin_vector_dimensions().

Referenced by calculateStress().

calculateStress()

template<typename T_VECTOR , typename T_MATRIX , int DisplacementDim>

void MaterialLib::Solids::Phasefield::calculateStress	(	T_VECTOR &	sigma,
		T_VECTOR &	sigma_tensile,
		T_VECTOR &	sigma_compressive,
		T_VECTOR &	eps_tensile,
		T_MATRIX &	D,
		T_MATRIX &	C_tensile,
		T_MATRIX &	C_compressive,
		double &	strain_energy_tensile,
		double &	elastic_energy,
		double const	degradation,
		T_VECTOR const &	eps,
		EnergySplitModel const &	energy_split_model,
		double const	t,
		ParameterLib::SpatialPosition const &	x,
		MaterialLib::Solids::MechanicsBase< DisplacementDim > const &	solid_material )

Definition at line 342 of file PhaseFieldBase.h.

{
    auto linear_elastic_mp =
        static_cast<MaterialLib::Solids::LinearElasticIsotropic<
            DisplacementDim> const&>(solid_material)
            .getMaterialProperties();
 
    auto const bulk_modulus = linear_elastic_mp.bulk_modulus(t, x);
    auto const mu = linear_elastic_mp.mu(t, x);
    switch (energy_split_model)
    {
        case EnergySplitModel::Isotropic:
        {
            std::tie(sigma, sigma_tensile, D, strain_energy_tensile,
                     elastic_energy, C_tensile, C_compressive) =
                MaterialLib::Solids::Phasefield::
                    calculateIsotropicDegradedStress<DisplacementDim>(
                        degradation, bulk_modulus, mu, eps);
            break;
        }
        case EnergySplitModel::VolDev:
        {
            std::tie(sigma, sigma_tensile, D, strain_energy_tensile,
                     elastic_energy, C_tensile, C_compressive) =
                MaterialLib::Solids::Phasefield::calculateVolDevDegradedStress<
                    DisplacementDim>(degradation, bulk_modulus, mu, eps);
            break;
        }
        case EnergySplitModel::EffectiveStress:
        {
            std::tie(sigma, sigma_tensile, D, strain_energy_tensile,
                     elastic_energy, C_tensile, C_compressive) =
                MaterialLib::Solids::Phasefield::
                    calculateIsotropicDegradedStressWithRankineEnergy<
                        DisplacementDim>(degradation, bulk_modulus, mu, eps);
            break;
        }
        case EnergySplitModel::OrthoVolDev:
        {
            double temperature = 0.;
            MathLib::KelvinVector::KelvinMatrixType<DisplacementDim> C_ortho =
                static_cast<MaterialLib::Solids::LinearElasticOrthotropic<
                    DisplacementDim> const&>(solid_material)
                    .getElasticTensor(t, x, temperature);
 
            std::tie(eps_tensile, sigma, sigma_tensile, sigma_compressive, D,
                     strain_energy_tensile, elastic_energy, C_tensile,
                     C_compressive) = MaterialLib::Solids::Phasefield::
                calculateOrthoVolDevDegradedStress<DisplacementDim>(
                    degradation, eps, C_ortho);
            break;
        }
        case EnergySplitModel::OrthoMasonry:
        {
            double temperature = 0.;
            MathLib::KelvinVector::KelvinMatrixType<DisplacementDim> C_ortho =
                static_cast<MaterialLib::Solids::LinearElasticOrthotropic<
                    DisplacementDim> const&>(solid_material)
                    .getElasticTensor(t, x, temperature);
 
            std::tie(eps_tensile, sigma, sigma_tensile, D,
                     strain_energy_tensile, elastic_energy, C_tensile,
                     C_compressive) = MaterialLib::Solids::Phasefield::
                calculateOrthoMasonryDegradedStress<DisplacementDim>(
                    degradation, eps, C_ortho);
            break;
        }
        default:
            OGS_FATAL("Invalid energy split model specified.");
    }
};

References calculateIsotropicDegradedStress(), calculateOrthoMasonryDegradedStress(), calculateOrthoVolDevDegradedStress(), calculateVolDevDegradedStress(), EffectiveStress, Isotropic, OGS_FATAL, OrthoMasonry, OrthoVolDev, and VolDev.

Referenced by ProcessLib::HMPhaseField::IntegrationPointData< BMatricesType, ShapeMatrixType, DisplacementDim >::updateConstitutiveRelation(), and ProcessLib::PhaseField::IntegrationPointData< BMatricesType, ShapeMatrixType, DisplacementDim >::updateConstitutiveRelation().

calculateVolDevDegradedStress()

template<int DisplacementDim>

std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > > MaterialLib::Solids::Phasefield::calculateVolDevDegradedStress	(	double const	degradation,
		double const	bulk_modulus,
		double const	mu,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps )

Definition at line 61 of file LinearElasticIsotropicPhaseField.h.

{
    static int const KelvinVectorSize =
        MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim);
    using KelvinVector =
        MathLib::KelvinVector::KelvinVectorType<DisplacementDim>;
    using KelvinMatrix =
        MathLib::KelvinVector::KelvinMatrixType<DisplacementDim>;
    using Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>;
    // calculation of deviatoric parts
    auto const& P_dev = Invariants::deviatoric_projection;
    KelvinVector const epsd_curr = P_dev * eps;
 
    // Hydrostatic part for the stress and the tangent.
    double const eps_curr_trace = Invariants::trace(eps);
 
    KelvinMatrix C_tensile = KelvinMatrix::Zero();
    KelvinMatrix C_compressive = KelvinMatrix::Zero();
 
    auto strain_energy_computation_vol = [&](auto&& macaulay)
    {
        auto macaulay_squared = [&macaulay](double x)
        { return boost::math::pow<2>(macaulay(x)); };
        return bulk_modulus / 2 * macaulay_squared(eps_curr_trace);
    };
 
    auto stress_computation_vol = [&](auto&& macaulay)
    { return bulk_modulus * macaulay(eps_curr_trace) * Invariants::identity2; };
 
    auto hs = [&](double const v) { return heaviside(v); };
 
    auto mt = [&](double const v) { return macaulayTensile(v); };
 
    auto mc = [&](double const v) { return macaulayCompressive(v); };
 
    double const strain_energy_tensile = strain_energy_computation_vol(mt) +
                                         mu * epsd_curr.transpose() * epsd_curr;
 
    KelvinVector const sigma_tensile =
        stress_computation_vol(mt) + 2 * mu * epsd_curr;
 
    KelvinVector const sigma_compressive = stress_computation_vol(mc);
 
    C_tensile.template topLeftCorner<3, 3>().setConstant(bulk_modulus *
                                                         hs(eps_curr_trace));
    C_tensile.noalias() += 2 * mu * P_dev * KelvinMatrix::Identity();
 
    C_compressive.template topLeftCorner<3, 3>().setConstant(
        bulk_modulus * (1 - hs(eps_curr_trace)));
 
    double const elastic_energy =
        bulk_modulus / 2 * eps_curr_trace * eps_curr_trace +
        mu * epsd_curr.transpose() * epsd_curr;
    KelvinVector const sigma_real =
        degradation * sigma_tensile + sigma_compressive;
 
    KelvinMatrix const D = degradation * C_tensile + C_compressive;
 
    return std::make_tuple(sigma_real, sigma_tensile, D, strain_energy_tensile,
                           elastic_energy, C_tensile, C_compressive);
}

References heaviside(), MathLib::KelvinVector::kelvin_vector_dimensions(), macaulayCompressive(), and macaulayTensile().

Referenced by calculateStress(), and ProcessLib::ThermoMechanicalPhaseField::IntegrationPointData< BMatricesType, ShapeMatrixType, DisplacementDim >::updateConstitutiveRelation().

convertStringToEnergySplitModel()

template<int DisplacementDim>

EnergySplitModel MaterialLib::Solids::Phasefield::convertStringToEnergySplitModel

(

std::string const &

energy_split_model

)

Definition at line 95 of file PhaseFieldBase.h.

{
    if (energy_split_model == "Isotropic")
    {
        return EnergySplitModel::Isotropic;
    }
    if (energy_split_model == "VolumetricDeviatoric")
    {
        return EnergySplitModel::VolDev;
    }
    if (energy_split_model == "EffectiveStress")
    {
        return EnergySplitModel::EffectiveStress;
    }
    if (energy_split_model == "OrthoVolDev")
    {
        return EnergySplitModel::OrthoVolDev;
    }
    if (energy_split_model == "OrthoMasonry")
    {
        return EnergySplitModel::OrthoMasonry;
    }
    OGS_FATAL(
        "energy_split_model must be 'Isotropic', 'VolumetricDeviatoric', "
        "'EffectiveStress', 'OrthoVolDev' or 'OrthoMasonry' but '{}' was given",
        energy_split_model.c_str());
};

References EffectiveStress, Isotropic, OGS_FATAL, OrthoMasonry, OrthoVolDev, and VolDev.

Referenced by ProcessLib::HMPhaseField::createHMPhaseFieldProcess(), and ProcessLib::PhaseField::createPhaseFieldProcess().

convertStringToPhaseFieldModel()

template<int DisplacementDim>

PhaseFieldModel MaterialLib::Solids::Phasefield::convertStringToPhaseFieldModel

(

std::string const &

phasefield_model

)

Definition at line 36 of file PhaseFieldBase.h.

{
    if (phasefield_model == "AT1")
    {
        return PhaseFieldModel::AT1;
    }
    if (phasefield_model == "AT2")
    {
        return PhaseFieldModel::AT2;
    }
    if (phasefield_model == "COHESIVE")
    {
        return PhaseFieldModel::COHESIVE;
    }
    OGS_FATAL(
        "phasefield_model must be 'AT1', 'AT2' or 'COHESIVE' but '{}' "
        "was given",
        phasefield_model.c_str());
};

References AT1, AT2, COHESIVE, and OGS_FATAL.

Referenced by ProcessLib::HMPhaseField::createHMPhaseFieldProcess(), and ProcessLib::PhaseField::createPhaseFieldProcess().

convertStringToSofteningCurve()

template<int DisplacementDim>

SofteningCurve MaterialLib::Solids::Phasefield::convertStringToSofteningCurve

(

std::optional< std::string > const &

softening_curve

)

Definition at line 64 of file PhaseFieldBase.h.

{
    if (softening_curve)
    {
        if (*softening_curve == "Linear")
        {
            return SofteningCurve::Linear;
        }
        if (*softening_curve == "Exponential")
        {
            return SofteningCurve::Exponential;
        }
        OGS_FATAL(
            "softening_curve must be 'Linear' or 'Exponential' but '{}' "
            "was given",
            softening_curve->c_str());
    }
    return SofteningCurve::Linear;  // default
};

References Exponential, Linear, and OGS_FATAL.

Referenced by ProcessLib::HMPhaseField::createHMPhaseFieldProcess(), and ProcessLib::PhaseField::createPhaseFieldProcess().

creatDegradationDerivative()

template<int DisplacementDim>

std::unique_ptr< DegradationDerivative > MaterialLib::Solids::Phasefield::creatDegradationDerivative	(	PhaseFieldModel const &	phasefield_model,
		double const	lch,
		SofteningCurve const &	softening_curve )

Definition at line 276 of file PhaseFieldBase.h.

{
    switch (phasefield_model)
    {
        case PhaseFieldModel::COHESIVE:
            return std::make_unique<COHESIVE_DegradationDerivative>(
                lch, softening_curve);
        default:
            return std::make_unique<AT_DegradationDerivative>();
    }
};

References COHESIVE.

heaviside()

double MaterialLib::Solids::Phasefield::heaviside ( double const v )

inline

Decompose the stiffness into tensile and compressive part. Judging by the physical observations, compression perpendicular to a crack does not cause crack propagation. Thus, the phase-field parameter is only involved into the tensile part to degrade the elastic strain energy. heaviside function returns 1.0 if the argument is positive and 0.0 if negative

Definition at line 32 of file LinearElasticIsotropicPhaseField.h.

{
    return (v < 0) ? 0.0 : 1.0;
}

Referenced by calculateIsotropicDegradedStressWithRankineEnergy(), calculateOrthoMasonryDegradedStress(), calculateVolDevDegradedStress(), macaulayCompressive(), and macaulayTensile().

macaulayCompressive()

double MaterialLib::Solids::Phasefield::macaulayCompressive ( double v )

inline

Definition at line 43 of file LinearElasticIsotropicPhaseField.h.

{
    return v * (1 - heaviside(v));
}

References heaviside().

Referenced by calculateOrthoMasonryDegradedStress(), and calculateVolDevDegradedStress().

macaulayTensile()

double MaterialLib::Solids::Phasefield::macaulayTensile ( double const v )

inline

Macaulay brackets: positive strain is tensile and negative strain for compressive

Definition at line 39 of file LinearElasticIsotropicPhaseField.h.

{
    return v * heaviside(v);
}

References heaviside().

Referenced by calculateOrthoMasonryDegradedStress(), and calculateVolDevDegradedStress().