Revision 3f51786289e66cdebc74007458f43157b17ca0b8 authored by Dmitri Naumov on 01 October 2021, 11:32:23 UTC, committed by Dmitri Naumov on 01 October 2021, 18:38:49 UTC
Update the ctest config for updated time step output. The reference files are still the same, no changes required.
1 parent f63a01b
Coulomb.h
/**
* \file
* \copyright
* Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
* Distributed under a Modified BSD License.
* See accompanying file LICENSE.txt or
* http://www.opengeosys.org/project/license
*/
#pragma once
#include <Eigen/Eigen>
#include <utility>
#include "NumLib/NewtonRaphson.h"
#include "ParameterLib/Parameter.h"
#include "FractureModelBase.h"
namespace MaterialLib
{
namespace Fracture
{
namespace Coulomb
{
template <int DisplacementDim>
struct StateVariables
: public FractureModelBase<DisplacementDim>::MaterialStateVariables
{
void setInitialConditions() { w_p = w_p_prev; }
void pushBackState() override { w_p_prev = w_p; }
/// Plastic component of the displacement jump in fracture's local
/// coordinates.
Eigen::Matrix<double, DisplacementDim, 1> w_p =
Eigen::Matrix<double, DisplacementDim, 1>::Zero();
// Initial values from previous timestep
Eigen::Matrix<double, DisplacementDim, 1> w_p_prev; ///< \copydoc w_p
EIGEN_MAKE_ALIGNED_OPERATOR_NEW;
};
template <int DisplacementDim>
class Coulomb final : public FractureModelBase<DisplacementDim>
{
public:
std::unique_ptr<
typename FractureModelBase<DisplacementDim>::MaterialStateVariables>
createMaterialStateVariables() override
{
return std::make_unique<StateVariables<DisplacementDim>>();
}
public:
/// Variables specific to the material model
struct MaterialProperties
{
using P = ParameterLib::Parameter<double>;
using X = ParameterLib::SpatialPosition;
MaterialProperties(
P const& normal_stiffness_, P const& shear_stiffness_,
P const& friction_angle_, P const& dilatancy_angle_,
P const& cohesion_)
: normal_stiffness(normal_stiffness_), shear_stiffness(shear_stiffness_),
friction_angle(friction_angle_), dilatancy_angle(dilatancy_angle_),
cohesion(cohesion_)
{
}
/// Normal stiffness given in units of stress per length.
P const& normal_stiffness;
/// Shear stiffness given in units of stress per length.
P const& shear_stiffness;
/// Governs the normal stress dependence of the shear strength.
/// \note Given in degrees (not radian).
P const& friction_angle;
/// Governs the dilatancy behaviour during the plastic deformation of
/// the fault.
/// \note Given in degrees (not radian).
P const& dilatancy_angle;
/// Fracture cohesion in units of stress.
P const& cohesion;
};
public:
explicit Coulomb(
NumLib::NewtonRaphsonSolverParameters nonlinear_solver_parameters,
double const penalty_aperture_cutoff,
bool const tension_cutoff,
MaterialProperties material_properties)
: _nonlinear_solver_parameters(std::move(nonlinear_solver_parameters)),
_penalty_aperture_cutoff(penalty_aperture_cutoff),
_tension_cutoff(tension_cutoff),
_mp(std::move(material_properties))
{
}
/**
* Computation of the constitutive relation for the Mohr-Coulomb model.
*
* @param t current time
* @param x current position in space
* @param aperture0 initial fracture's aperture
* @param sigma0 initial stress
* @param w_prev fracture displacement at previous time step
* @param w fracture displacement at current time step
* @param sigma_prev stress at previous time step
* @param sigma stress at current time step
* @param Kep tangent matrix for stress and fracture displacements
* @param material_state_variables material state variables
*/
void computeConstitutiveRelation(
double const t,
ParameterLib::SpatialPosition const& x,
double const aperture0,
Eigen::Ref<Eigen::VectorXd const>
sigma0,
Eigen::Ref<Eigen::VectorXd const>
w_prev,
Eigen::Ref<Eigen::VectorXd const>
w,
Eigen::Ref<Eigen::VectorXd const>
sigma_prev,
Eigen::Ref<Eigen::VectorXd>
sigma,
Eigen::Ref<Eigen::MatrixXd>
Kep,
typename FractureModelBase<DisplacementDim>::MaterialStateVariables&
material_state_variables) override;
private:
NumLib::NewtonRaphsonSolverParameters const _nonlinear_solver_parameters;
/// Compressive normal displacements above this value will not enter the
/// computation of the normal stiffness modulus of the fracture.
/// \note Setting this to the initial aperture value allows negative
/// apertures.
double const _penalty_aperture_cutoff;
/// If set no resistance to open the fracture over the initial aperture is
/// opposed.
bool const _tension_cutoff;
MaterialProperties _mp;
};
} // namespace Coulomb
} // namespace Fracture
} // namespace MaterialLib
namespace MaterialLib
{
namespace Fracture
{
namespace Coulomb
{
extern template class Coulomb<2>;
extern template class Coulomb<3>;
} // namespace Coulomb
} // namespace Fracture
} // namespace MaterialLib
Computing file changes ...