Commit 31ef9294 authored by Antoine Jacquey's avatar Antoine Jacquey

Fixed poromechanics.

parent 654ba26d
...@@ -49,8 +49,8 @@ protected: ...@@ -49,8 +49,8 @@ protected:
const MaterialProperty<Real> & _K; const MaterialProperty<Real> & _K;
const MaterialProperty<RankFourTensor> & _tangent_modulus; const MaterialProperty<RankFourTensor> & _tangent_modulus;
const MaterialProperty<RankTwoTensor> & _strain_increment; const MaterialProperty<RankTwoTensor> & _strain_increment;
const MaterialProperty<RankTwoTensor> & _inelastic_strain; const MaterialProperty<RankTwoTensor> & _viscous_strain_incr;
const MaterialProperty<RankTwoTensor> & _inelastic_strain_old; const MaterialProperty<RankTwoTensor> & _plastic_strain_incr;
MaterialProperty<Real> & _biot; MaterialProperty<Real> & _biot;
MaterialProperty<Real> & _C_d; MaterialProperty<Real> & _C_d;
MaterialProperty<Real> & _C_biot; MaterialProperty<Real> & _C_biot;
......
...@@ -120,7 +120,7 @@ Real ...@@ -120,7 +120,7 @@ Real
LynxSolidMomentum::computeQpResidual() LynxSolidMomentum::computeQpResidual()
{ {
RealVectorValue stress_row = _stress[_qp].row(_component); RealVectorValue stress_row = _stress[_qp].row(_component);
stress_row(_component) += _biot[_qp] * _pf[_qp]; stress_row(_component) -= _biot[_qp] * _pf[_qp];
RealVectorValue grav_term = -_rho_b[_qp] * _gravity[_qp]; RealVectorValue grav_term = -_rho_b[_qp] * _gravity[_qp];
Real residual = stress_row * _grad_test[_i][_qp] + grav_term(_component) * _test[_i][_qp]; Real residual = stress_row * _grad_test[_i][_qp] + grav_term(_component) * _test[_i][_qp];
......
...@@ -271,6 +271,9 @@ LynxDeformationBase::computeStrainIncrement() ...@@ -271,6 +271,9 @@ LynxDeformationBase::computeStrainIncrement()
mooseError("Unknown strain model. Specify 'small' or 'finite'!"); mooseError("Unknown strain model. Specify 'small' or 'finite'!");
} }
// Thermal strain correction
_strain_increment[_qp] -= _thermal_strain_incr[_qp];
if (_vol_locking_correction) if (_vol_locking_correction)
vol_strain_incr += _strain_increment[_qp].trace() * _JxW[_qp] * _coord[_qp]; vol_strain_incr += _strain_increment[_qp].trace() * _JxW[_qp] * _coord[_qp];
} }
...@@ -351,7 +354,7 @@ LynxDeformationBase::initializeQpDeformation() ...@@ -351,7 +354,7 @@ LynxDeformationBase::initializeQpDeformation()
{ {
// Initialze elastic strain // Initialze elastic strain
_elastic_strain[_qp] = _elastic_strain[_qp] =
spinRotation(_elastic_strain_old[_qp]) + _strain_increment[_qp] - _thermal_strain_incr[_qp]; spinRotation(_elastic_strain_old[_qp]) + _strain_increment[_qp];
// Initialize inelastic increment // Initialize inelastic increment
_viscous_strain_incr[_qp].zero(); _viscous_strain_incr[_qp].zero();
......
...@@ -36,8 +36,8 @@ LynxHydroBase::LynxHydroBase(const InputParameters & parameters) ...@@ -36,8 +36,8 @@ LynxHydroBase::LynxHydroBase(const InputParameters & parameters)
_K(getDefaultMaterialProperty<Real>("bulk_modulus")), _K(getDefaultMaterialProperty<Real>("bulk_modulus")),
_tangent_modulus(getDefaultMaterialProperty<RankFourTensor>("tangent_modulus")), _tangent_modulus(getDefaultMaterialProperty<RankFourTensor>("tangent_modulus")),
_strain_increment(getDefaultMaterialProperty<RankTwoTensor>("strain_increment")), _strain_increment(getDefaultMaterialProperty<RankTwoTensor>("strain_increment")),
_inelastic_strain(getDefaultMaterialProperty<RankTwoTensor>("inelastic_strain")), _viscous_strain_incr(getDefaultMaterialProperty<RankTwoTensor>("viscous_strain_increment")),
_inelastic_strain_old(getMaterialPropertyOld<RankTwoTensor>("inelastic_strain")), _plastic_strain_incr(getDefaultMaterialProperty<RankTwoTensor>("plastic_strain_increment")),
_biot(declareProperty<Real>("biot_coefficient")), _biot(declareProperty<Real>("biot_coefficient")),
_C_d(declareProperty<Real>("bulk_compressibility")), _C_d(declareProperty<Real>("bulk_compressibility")),
_C_biot(declareProperty<Real>("biot_compressibility")), _C_biot(declareProperty<Real>("biot_compressibility")),
...@@ -97,7 +97,7 @@ LynxHydroBase::computeQpPoroMech() ...@@ -97,7 +97,7 @@ LynxHydroBase::computeQpPoroMech()
{ {
Real K_cto = _tangent_modulus[_qp].sum3x3() / 9.0; Real K_cto = _tangent_modulus[_qp].sum3x3() / 9.0;
RankTwoTensor e_tot = _strain_increment[_qp] / _dt; RankTwoTensor e_tot = _strain_increment[_qp] / _dt;
RankTwoTensor e_in = (_inelastic_strain[_qp] - _inelastic_strain_old[_qp]) / _dt; RankTwoTensor e_in = (_viscous_strain_incr[_qp] + _plastic_strain_incr[_qp]) / _dt;
_poro_mech[_qp] = _biot[_qp] * e_tot.trace() + (1.0 - _biot[_qp]) * e_in.trace(); _poro_mech[_qp] = _biot[_qp] * e_tot.trace() + (1.0 - _biot[_qp]) * e_in.trace();
_poro_mech_jac[_qp] = _biot[_qp] + (1.0 - _biot[_qp]) * (1.0 - K_cto / _K[_qp]); _poro_mech_jac[_qp] = _biot[_qp] + (1.0 - _biot[_qp]) * (1.0 - K_cto / _K[_qp]);
......
...@@ -27,7 +27,7 @@ InputParameters ...@@ -27,7 +27,7 @@ InputParameters
validParams<LynxHydroConstant>() validParams<LynxHydroConstant>()
{ {
InputParameters params = validParams<LynxHydroBase>(); InputParameters params = validParams<LynxHydroBase>();
params.addClassDescription("Constant thermal properties."); params.addClassDescription("Constant hydraulic properties.");
params.addRequiredParam<std::vector<Real>>("permeability", "The permeability of the matrix."); params.addRequiredParam<std::vector<Real>>("permeability", "The permeability of the matrix.");
params.addRequiredParam<std::vector<Real>>("fluid_viscosity", "The viscosity of the fluid."); params.addRequiredParam<std::vector<Real>>("fluid_viscosity", "The viscosity of the fluid.");
params.addParam<std::vector<Real>>("fluid_modulus", "The bulk modulus of the fluid phase."); params.addParam<std::vector<Real>>("fluid_modulus", "The bulk modulus of the fluid phase.");
......
...@@ -42,21 +42,21 @@ ...@@ -42,21 +42,21 @@
variable = disp_x variable = disp_x
component = 0 component = 0
displacements = 'disp_x disp_y disp_z' displacements = 'disp_x disp_y disp_z'
fluid_pressure = pf fluid_pressure = 'pf'
[../] [../]
[./mech_y] [./mech_y]
type = LynxSolidMomentum type = LynxSolidMomentum
variable = disp_y variable = disp_y
component = 1 component = 1
displacements = 'disp_x disp_y disp_z' displacements = 'disp_x disp_y disp_z'
fluid_pressure = pf fluid_pressure = 'pf'
[../] [../]
[./mech_z] [./mech_z]
type = LynxSolidMomentum type = LynxSolidMomentum
variable = disp_z variable = disp_z
component = 2 component = 2
displacements = 'disp_x disp_y disp_z' displacements = 'disp_x disp_y disp_z'
fluid_pressure = pf fluid_pressure = 'pf'
[../] [../]
[] []
...@@ -147,7 +147,7 @@ ...@@ -147,7 +147,7 @@
[../] [../]
[./hydro_mat] [./hydro_mat]
type = LynxHydroConstant type = LynxHydroConstant
porosity = porosity porosity = 'porosity'
fluid_modulus = 8 fluid_modulus = 8
solid_modulus = 2.5 solid_modulus = 2.5
permeability = 1.5e-03 permeability = 1.5e-03
...@@ -276,4 +276,4 @@ ...@@ -276,4 +276,4 @@
interval = 3 interval = 3
type = CSV type = CSV
[../] [../]
[] []
\ No newline at end of file
#[Tests] [Tests]
# [./mandel] [./mandel]
# type = 'Exodiff' type = 'Exodiff'
# input = 'mandel.i' input = 'mandel.i'
# exodiff = 'mandel_out.e' exodiff = 'mandel_out.e'
# [../] [../]
#[] []
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment