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 , Spectral , EffectiveStress ,   OrthoVolDev , OrthoMasonry }

Functions
double	evaluateHTensSpectral (int const i, int const j, Eigen::Matrix< double, 3, 1 > const &principal_strain)
	H terms in the Spectral decomposition:
double	evaluateHCompSpectral (int const i, int const j, Eigen::Matrix< double, 3, 1 > const &principal_strain)
template<>
MathLib::KelvinVector::KelvinMatrixType< 3 >	aOdotB< 3 > (MathLib::KelvinVector::KelvinVectorType< 3 > const &A, MathLib::KelvinVector::KelvinVectorType< 3 > const &B)
template<>
MathLib::KelvinVector::KelvinMatrixType< 2 >	aOdotB< 2 > (MathLib::KelvinVector::KelvinVectorType< 2 > const &A, MathLib::KelvinVector::KelvinVectorType< 2 > const &B)
template<int DisplacementDim>
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >	aOdotB (MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &A, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &B)
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 > >	calculateSpectralDegradedStress (double const degradation, double const lambda, 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
Spectral
EffectiveStress
OrthoVolDev
OrthoMasonry

Definition at line 79 of file PhaseFieldBase.h.

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

PhaseFieldModel

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

strong

Enumerator
AT1
AT2
COHESIVE

Definition at line 22 of file PhaseFieldBase.h.

{
    AT1,
    AT2,
    COHESIVE
};

SofteningCurve

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

strong

Enumerator
Linear
Exponential

Definition at line 51 of file PhaseFieldBase.h.

{
    Linear,
    Exponential
};

Function Documentation

aOdotB()

template<int DisplacementDim>

MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > MaterialLib::Solids::Phasefield::aOdotB	(	MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	A,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	B )

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.

Referenced by calculateSpectralDegradedStress().

aOdotB< 2 >()

template<>

MathLib::KelvinVector::KelvinMatrixType< 2 > MaterialLib::Solids::Phasefield::aOdotB< 2 >	(	MathLib::KelvinVector::KelvinVectorType< 2 > const &	A,
		MathLib::KelvinVector::KelvinVectorType< 2 > const &	B )

Definition at line 27 of file LinearElasticIsotropicPhaseField.cpp.

{
    MathLib::KelvinVector::KelvinMatrixType<2> result;
 
    result(0, 0) = A(0) * B(0);
    result(0, 1) = result(1, 0) = A(3) * B(3) / 2.;
    result(0, 2) = result(2, 0) = 0;
    result(0, 3) = result(3, 0) = (A(0) * B(3) + A(3) * B(0)) / 2;
 
    result(1, 1) = A(1) * B(1);
    result(1, 2) = result(2, 1) = 0;
    result(1, 3) = result(3, 1) = (A(1) * B(3) + A(3) * B(1)) / 2.;
 
    result(2, 2) = A(2) * B(2);
    result(2, 3) = result(3, 2) = 0;
 
    result(3, 3) = (A(0) * B(1) + A(1) * B(0) + A(3) * B(3)) / 2.;
 
    return result;
}

References heaviside(), and macaulayCompressive().

aOdotB< 3 >()

template<>

MathLib::KelvinVector::KelvinMatrixType< 3 > MaterialLib::Solids::Phasefield::aOdotB< 3 >	(	MathLib::KelvinVector::KelvinVectorType< 3 > const &	A,
		MathLib::KelvinVector::KelvinVectorType< 3 > const &	B )

Double-minor symmetrized tensor product

aOdotB computes the symmetric 4th-order tensor product

\[(A \odot B)_{ijkl} = \tfrac{1}{4} (A_{ik}B_{jl} + A_{il}B_{jk} + A_{jk}B_{il} + A_{jl}B_{ik}), \]

which is symmetric in both index pairs \((i,j)\) and \((k,l)\) (“double-minor symmetry”).

Reference Mehrabadi, M.M., & Cowin, S.C. (1990). Eigentensors of linear anisotropic elastic materials. Q. J. Mech. Appl. Math., 43(1), 15-41.

Definition at line 27 of file LinearElasticIsotropicPhaseField.cpp.

{
    MathLib::KelvinVector::KelvinMatrixType<3> result;
 
    result(0, 0) = A(0) * B(0);
    result(0, 1) = result(1, 0) = A(3) * B(3) / 2.;
    result(0, 2) = result(2, 0) = A(5) * B(5) / 2.;
    result(0, 3) = result(3, 0) = (A(0) * B(3) + A(3) * B(0)) / 2;
    result(0, 4) = result(4, 0) =
        (A(3) * B(5) + A(5) * B(3)) / (2. * std::sqrt(2.));
    result(0, 5) = result(5, 0) = (A(0) * B(5) + A(5) * B(0)) / 2;
 
    result(1, 1) = A(1) * B(1);
    result(1, 2) = result(2, 1) = A(4) * B(4) / 2.;
    result(1, 3) = result(3, 1) = (A(1) * B(3) + A(3) * B(1)) / 2.;
    result(1, 4) = result(4, 1) = (A(1) * B(4) + A(4) * B(1)) / 2.;
    result(1, 5) = result(5, 1) =
        (A(3) * B(4) + A(4) * B(3)) / (2. * std::sqrt(2.));
 
    result(2, 2) = A(2) * B(2);
    result(2, 3) = result(3, 2) =
        (A(4) * B(5) + A(5) * B(4)) / (2. * std::sqrt(2.));
    result(2, 4) = result(4, 2) = (A(2) * B(4) + A(4) * B(2)) / 2.;
    result(2, 5) = result(5, 2) = (A(2) * B(5) + A(5) * B(2)) / 2.;
 
    result(3, 3) = (A(0) * B(1) + A(1) * B(0) + A(3) * B(3)) / 2.;
    result(3, 4) = result(4, 3) =
        (A(3) * B(4) + A(4) * B(3)) / 4. +
        (A(1) * B(5) + A(5) * B(1)) / (2. * std::sqrt(2.));
    result(3, 5) = result(5, 3) =
        (A(3) * B(5) + A(5) * B(3)) / 4. +
        (A(0) * B(4) + A(4) * B(0)) / (2. * std::sqrt(2.));
 
    result(4, 4) = (A(1) * B(2) + A(2) * B(1) + A(4) * B(4)) / 2.;
    result(4, 5) = result(5, 4) =
        (A(4) * B(5) + A(5) * B(4)) / 4. +
        (A(2) * B(3) + A(3) * B(2)) / (2. * std::sqrt(2.));
 
    result(5, 5) = (A(0) * B(2) + A(2) * B(0) + A(5) * B(5)) / 2.;
    return result;
}

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 291 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 267 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 320 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 124 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 35 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().

calculateSpectralDegradedStress()

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::calculateSpectralDegradedStress	(	double const	degradation,
		double const	lambda,
		double const	mu,
		MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &	eps )

Definition at line 146 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>;
 
    KelvinMatrix C_tensile = KelvinMatrix::Zero();
    KelvinMatrix C_compressive = KelvinMatrix::Zero();
 
    // non-const for Eigen solver.
    auto eps_tensor = MathLib::KelvinVector::kelvinVectorToTensor(eps);
 
    Eigen::EigenSolver<decltype(eps_tensor)> eigen_solver(eps_tensor);
    Eigen::Matrix<double, 3, 1> const principal_strain =
        eigen_solver.eigenvalues().real();
    double const sum_strain = principal_strain.sum();
 
    std::array<KelvinVector, 3> M_kelvin;
 
    for (int i = 0; i < 3; i++)
    {
        M_kelvin[i] = MathLib::KelvinVector::tensorToKelvin<DisplacementDim>(
            eigen_solver.eigenvectors().real().col(i).normalized() *
            eigen_solver.eigenvectors().real().col(i).normalized().transpose());
    }
 
    auto strain_energy_computation = [&](auto&& macaulay)
    {
        auto macaulay_squared = [&macaulay](double x)
        { return boost::math::pow<2>(macaulay(x)); };
        return lambda / 2 * macaulay_squared(sum_strain) +
               mu * principal_strain.unaryExpr(macaulay_squared).sum();
    };
 
    auto stress_computation = [&](auto&& macaulay)
    {
        KelvinVector stress =
            lambda * macaulay(sum_strain) * Invariants::identity2;
        for (int i = 0; i < 3; i++)
        {
            stress += 2 * mu * macaulay(principal_strain[i]) * M_kelvin[i];
        }
        return stress;
    };
 
    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(mt);
 
    double const strain_energy_compressive = strain_energy_computation(mc);
 
    KelvinVector const sigma_tensile = stress_computation(mt);
 
    KelvinVector const sigma_compressive = stress_computation(mc);
 
    C_tensile.template topLeftCorner<3, 3>().setConstant(lambda *
                                                         hs(sum_strain));
 
    for (int i = 0; i < 3; i++)
    {
        C_tensile.noalias() += 2 * mu * hs(principal_strain[i]) * M_kelvin[i] *
                               M_kelvin[i].transpose();
        for (int j = 0; j < 3; j++)
        {
            C_tensile.noalias() +=
                mu * evaluateHTensSpectral(i, j, principal_strain) *
                aOdotB<DisplacementDim>(M_kelvin[i], M_kelvin[j]);
        }
    }
 
    C_compressive.template topLeftCorner<3, 3>().setConstant(
        lambda * (1 - hs(sum_strain)));
    KelvinMatrix C_temp = KelvinMatrix::Zero();
    for (int i = 0; i < 3; i++)
    {
        C_compressive.noalias() += 2 * mu * (1 - hs(principal_strain[i])) *
                                   M_kelvin[i] * M_kelvin[i].transpose();
        C_temp.noalias() = M_kelvin[i] * M_kelvin[i].transpose();
        for (int j = 0; j < 3; j++)
        {
            C_compressive.noalias() +=
                mu * evaluateHCompSpectral(i, j, principal_strain) *
                aOdotB<DisplacementDim>(M_kelvin[i], M_kelvin[j]);
        }
    }
 
    double const elastic_energy =
        degradation * strain_energy_tensile + strain_energy_compressive;
 
    KelvinVector 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 aOdotB(), evaluateHCompSpectral(), evaluateHTensSpectral(), heaviside(), MathLib::KelvinVector::kelvin_vector_dimensions(), MathLib::KelvinVector::kelvinVectorToTensor(), macaulayCompressive(), macaulayTensile(), and MathLib::KelvinVector::tensorToKelvin().

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 341 of file PhaseFieldBase.h.

{
    auto linear_elastic_mp =
        static_cast<MaterialLib::Solids::LinearElasticIsotropic<
            DisplacementDim> const&>(solid_material)
            .getMaterialProperties();
 
    auto const lambda = linear_elastic_mp.lambda(t, x);
    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::Spectral:
        {
            std::tie(sigma, sigma_tensile, D, strain_energy_tensile,
                     elastic_energy, C_tensile, C_compressive) =
                MaterialLib::Solids::Phasefield::
                    calculateSpectralDegradedStress<DisplacementDim>(
                        degradation, lambda, 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(), calculateSpectralDegradedStress(), calculateVolDevDegradedStress(), EffectiveStress, Isotropic, OGS_FATAL, OrthoMasonry, OrthoVolDev, Spectral, and VolDev.

Referenced by ProcessLib::HMPhaseField::IntegrationPointData< BMatricesType, ShapeMatricesType, DisplacementDim >::updateConstitutiveRelation(), and ProcessLib::PhaseField::IntegrationPointData< BMatricesType, ShapeMatricesType, 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 69 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().

convertStringToEnergySplitModel()

template<int DisplacementDim>

EnergySplitModel MaterialLib::Solids::Phasefield::convertStringToEnergySplitModel

(

std::string const &

energy_split_model

)

Definition at line 90 of file PhaseFieldBase.h.

{
    if (energy_split_model == "Isotropic")
    {
        return EnergySplitModel::Isotropic;
    }
    if (energy_split_model == "VolumetricDeviatoric")
    {
        return EnergySplitModel::VolDev;
    }
    if (energy_split_model == "Spectral")
    {
        return EnergySplitModel::Spectral;
    }
    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, Spectral, 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 30 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 58 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 275 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.

evaluateHCompSpectral()

double MaterialLib::Solids::Phasefield::evaluateHCompSpectral	(	int const	i,
		int const	j,
		Eigen::Matrix< double, 3, 1 > const &	principal_strain )

Definition at line 27 of file LinearElasticIsotropicPhaseField.cpp.

{
    if (i == j)
    {
        return 0.0;
    }
    double num_zero = 1.e-14;
 
    if (std::fabs(principal_strain[i] - principal_strain[j]) < num_zero)
    {
        return 2 * (1 - heaviside(principal_strain[i]));
    }
    return 2 *
           (macaulayCompressive(principal_strain[i]) -
            macaulayCompressive(principal_strain[j])) /
           (principal_strain[i] - principal_strain[j]);
}

Referenced by calculateSpectralDegradedStress().

evaluateHTensSpectral()

double MaterialLib::Solids::Phasefield::evaluateHTensSpectral	(	int const	i,
		int const	j,
		Eigen::Matrix< double, 3, 1 > const &	principal_strain )

H terms in the Spectral decomposition:

Definition at line 8 of file LinearElasticIsotropicPhaseField.cpp.

{
    if (i == j)
    {
        return 0.0;
    }
    if (std::fabs(principal_strain[i] - principal_strain[j]) <
        std::numeric_limits<double>::epsilon())
    {
        return 2 * heaviside(principal_strain[i]);
    }
    return 2 *
           (macaulayTensile(principal_strain[i]) -
            macaulayTensile(principal_strain[j])) /
           (principal_strain[i] - principal_strain[j]);
}

References heaviside(), and macaulayTensile().

Referenced by calculateSpectralDegradedStress().

heaviside()

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

inline

heaviside function returns 1.0 if the argument is positive and 0.0 if negative

Definition at line 31 of file LinearElasticIsotropicPhaseField.h.

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

Referenced by aOdotB< 2 >(), calculateIsotropicDegradedStressWithRankineEnergy(), calculateOrthoMasonryDegradedStress(), calculateSpectralDegradedStress(), calculateVolDevDegradedStress(), evaluateHTensSpectral(), macaulayCompressive(), and macaulayTensile().

macaulayCompressive()

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

inline

Definition at line 42 of file LinearElasticIsotropicPhaseField.h.

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

References heaviside().

Referenced by aOdotB< 2 >(), calculateOrthoMasonryDegradedStress(), calculateSpectralDegradedStress(), 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 38 of file LinearElasticIsotropicPhaseField.h.

{
    return v * heaviside(v);
}

References heaviside().

Referenced by calculateOrthoMasonryDegradedStress(), calculateSpectralDegradedStress(), calculateVolDevDegradedStress(), and evaluateHTensSpectral().