Commit 54e433b8 authored by Antoine Jacquey's avatar Antoine Jacquey

Added adiabatic heating

parent 5f040348
......@@ -56,6 +56,8 @@ protected:
std::vector<VariableName> _temperature;
bool _has_inelastic_heat_var;
std::vector<VariableName> _inelastic_heat;
bool _has_pressure_var;
std::vector<VariableName> _pressure;
std::vector<std::string> _aux_variables;
std::vector<std::string> _aux_kernels;
......
......@@ -21,6 +21,7 @@
using LynxADAdvectionBase<compute_stage>::_cr_stabilization; \
using LynxADAdvectionBase<compute_stage>::_value_old; \
using LynxADAdvectionBase<compute_stage>::_value_older; \
using LynxADAdvectionBase<compute_stage>::_vel; \
using LynxADAdvectionBase<compute_stage>::_vel_old; \
using LynxADAdvectionBase<compute_stage>::_vel_older; \
using LynxADAdvectionBase<compute_stage>::_entropy_old; \
......
......@@ -31,7 +31,7 @@ protected:
virtual void computeEntropyResidual();
const Real _coeff_Hs;
const ADVariableGradient & _grad_pressure;
const ADMaterialProperty(Real) & _thermal_diff;
const ADMaterialProperty(Real) & _rhoC;
const bool _has_inelastic_heat_mat;
......@@ -39,6 +39,7 @@ protected:
const ADMaterialProperty(Real) * _inelastic_heat_mat;
const bool _coupled_inelastic_heat;
const ADVariableValue & _inelastic_heat;
const ADMaterialProperty(Real) & _thermal_exp;
usingAdvectionBaseMembers;
};
\ No newline at end of file
......@@ -30,13 +30,16 @@ protected:
virtual ADReal computeQpResidual() override;
const Real _coeff_Hs;
unsigned int _nvel;
std::vector<const ADVariableValue *> _vel;
const ADVariableGradient & _grad_pressure;
const ADMaterialProperty(Real) & _rhoC_b;
const bool _has_inelastic_heat_mat;
const ADMaterialProperty(Real) * _radiogenic_heat;
const ADMaterialProperty(Real) * _inelastic_heat_mat;
const bool _coupled_inelastic_heat;
const ADVariableValue & _inelastic_heat;
const ADMaterialProperty(Real) & _thermal_exp;
usingKernelMembers;
};
\ No newline at end of file
......@@ -39,6 +39,8 @@ validParams<LynxADAdvectionAction>()
params.addParam<std::vector<VariableName>>(
"inelastic_heat",
"The auxiliary variable holding the inelastic heat value for running in a subApp.");
params.addParam<std::vector<VariableName>>("pressure",
"The pressure variable");
MooseEnum element_length_type_options("min=0 max=1 average=2", "min");
params.addParam<MooseEnum>(
"element_length_type", element_length_type_options, "The diameter of a single cell.");
......@@ -71,6 +73,9 @@ LynxADAdvectionAction::LynxADAdvectionAction(InputParameters params)
_has_inelastic_heat_var = isParamValid("inelastic_heat");
if (_has_inelastic_heat_var)
_inelastic_heat = getParam<std::vector<VariableName>>("inelastic_heat");
_has_pressure_var = isParamValid("pressure");
if (_has_pressure_var)
_pressure = getParam<std::vector<VariableName>>("pressure");
}
void
......@@ -449,6 +454,8 @@ LynxADAdvectionAction::createKernelActions()
params.set<std::vector<VariableName>>("velocities") = _velocities;
if (_has_inelastic_heat_var)
params.set<std::vector<VariableName>>("inelastic_heat") = _inelastic_heat;
if (_has_pressure_var)
params.set<std::vector<VariableName>>("pressure") = _pressure;
params.set<MooseEnum>("element_length_type") = _element_length_type;
params.set<Real>("beta_stabilization") = _beta_stabilization;
params.set<Real>("cr_stabilization") = _cr_stabilization;
......
......@@ -21,19 +21,23 @@ defineADValidParams(LynxADAdvectionTemperature, LynxADAdvectionBase,
"coeff_shear_heating", 0.0, "The coefficient in front of the shear heating generation.");
params.addCoupledVar(
"inelastic_heat",
"The auxiliary variable holding the inelastic heat value for running in a subApp."););
"The auxiliary variable holding the inelastic heat value for running in a subApp.");
params.addCoupledVar("pressure",
"The pressure variable"););
template <ComputeStage compute_stage>
LynxADAdvectionTemperature<compute_stage>::LynxADAdvectionTemperature(const InputParameters & parameters)
: LynxADAdvectionBase<compute_stage>(parameters),
_coeff_Hs(getParam<Real>("coeff_shear_heating")),
_grad_pressure(isCoupled("pressure") ? adCoupledGradient("pressure") : adZeroGradient()),
_thermal_diff(getADMaterialProperty<Real>("thermal_diffusivity")),
_rhoC(getADMaterialProperty<Real>("bulk_specific_heat")),
_has_inelastic_heat_mat(hasMaterialProperty<Real>("inelastic_heat")),
_radiogenic_heat(_has_inelastic_heat_mat ? &getADMaterialProperty<Real>("radiogenic_heat_production") : nullptr),
_inelastic_heat_mat(_has_inelastic_heat_mat ? &getADMaterialProperty<Real>("inelastic_heat") : nullptr),
_coupled_inelastic_heat(isCoupled("inelastic_heat")),
_inelastic_heat(_coupled_inelastic_heat ? adCoupledValue("inelastic_heat") : adZeroValue())
_inelastic_heat(_coupled_inelastic_heat ? adCoupledValue("inelastic_heat") : adZeroValue()),
_thermal_exp(getADMaterialProperty<Real>("thermal_expansion_coefficient"))
{
}
......@@ -109,7 +113,10 @@ LynxADAdvectionTemperature<compute_stage>::computeEntropyResidual()
else
Hs *= 0.0;
ADReal heat_sources = Hr + Hs;
ADRealVectorValue vel((*_vel[0])[_qp], (*_vel[1])[_qp], (*_vel[2])[_qp]);
ADReal Ha = _thermal_exp[_qp] * _u[_qp] * vel * _grad_pressure[_qp];
ADReal heat_sources = Hr + Hs + Ha;
if (_rhoC[_qp] != 0.0)
heat_sources /= _rhoC[_qp];
......
......@@ -25,19 +25,31 @@ defineADValidParams(
"The coefficient in front of the shear heating generation.");
params.addCoupledVar(
"inelastic_heat",
"The auxiliary variable holding the inelastic heat value for running in a subApp."););
"The auxiliary variable holding the inelastic heat value for running in a subApp.");
params.addCoupledVar("velocities",
"The string of velocities suitable for the problem statement");
params.addCoupledVar("pressure",
"The pressure variable"););
template <ComputeStage compute_stage>
LynxADHeatSources<compute_stage>::LynxADHeatSources(const InputParameters & parameters)
: ADKernel<compute_stage>(parameters),
_coeff_Hs(getParam<Real>("coeff_shear_heating")),
_nvel(coupledComponents("velocities")),
_vel(3),
_grad_pressure(isCoupled("pressure") ? adCoupledGradient("pressure") : adZeroGradient()),
_rhoC_b(getADMaterialProperty<Real>("bulk_specific_heat")),
_has_inelastic_heat_mat(hasMaterialProperty<Real>("inelastic_heat")),
_radiogenic_heat(_has_inelastic_heat_mat ? &getADMaterialProperty<Real>("radiogenic_heat_production") : nullptr),
_inelastic_heat_mat(_has_inelastic_heat_mat ? &getADMaterialProperty<Real>("inelastic_heat") : nullptr),
_coupled_inelastic_heat(isCoupled("inelastic_heat")),
_inelastic_heat(_coupled_inelastic_heat ? adCoupledValue("inelastic_heat") : adZeroValue())
_inelastic_heat(_coupled_inelastic_heat ? adCoupledValue("inelastic_heat") : adZeroValue()),
_thermal_exp(getADMaterialProperty<Real>("thermal_expansion_coefficient"))
{
for (unsigned i = 0; i < _nvel; ++i)
_vel[i] = &adCoupledValue("velocities", i);
for (unsigned i = _nvel; i < 3; ++i)
_vel[i] = &adZeroValue();
}
template <ComputeStage compute_stage>
......@@ -52,8 +64,10 @@ LynxADHeatSources<compute_stage>::computeQpResidual()
Hs *= (*_inelastic_heat_mat)[_qp];
else
Hs *= 0.0;
ADRealVectorValue vel((*_vel[0])[_qp], (*_vel[1])[_qp], (*_vel[2])[_qp]);
ADReal Ha = _thermal_exp[_qp] * _u[_qp] * vel * _grad_pressure[_qp];
ADReal heat_sources = Hr + Hs;
ADReal heat_sources = Hr + Hs + Ha;
if (_rhoC_b[_qp] != 0.0)
heat_sources /= _rhoC_b[_qp];
......
......@@ -167,10 +167,10 @@ LynxADElasticDeformation<compute_stage>::computeQpThermalSources()
inelastic_strain_incr += (*_viscous_strain_incr)[_qp];
if (_has_plastic)
inelastic_strain_incr += (*_plastic_strain_incr)[_qp];
if (_coupled_temp && !_coupled_temp_aux)
inelastic_strain_incr.addIa((*_thermal_exp)[_qp] * _temp_dot[_qp] * _dt / 3.0);
else if (_coupled_temp_aux && !_coupled_temp)
inelastic_strain_incr.addIa((*_thermal_exp)[_qp] * _temp_dot_aux[_qp] * _dt / 3.0);
// if (_coupled_temp && !_coupled_temp_aux)
// inelastic_strain_incr.addIa((*_thermal_exp)[_qp] * _temp_dot[_qp] * _dt / 3.0);
// else if (_coupled_temp_aux && !_coupled_temp)
// inelastic_strain_incr.addIa((*_thermal_exp)[_qp] * _temp_dot_aux[_qp] * _dt / 3.0);
_inelastic_heat[_qp] = _stress[_qp].doubleContraction(inelastic_strain_incr) / _dt;
}
\ No newline at end of file
......@@ -897,7 +897,7 @@ LynxDeformationBase::rootNewtonSafe(iterative_viscous & viscous_model, const Rea
void
LynxDeformationBase::computeQpThermalSources()
{
RankTwoTensor inelastic_strain_incr = _viscous_strain_incr[_qp] + _plastic_strain_incr[_qp] + _thermal_strain_incr[_qp];
RankTwoTensor inelastic_strain_incr = _viscous_strain_incr[_qp] + _plastic_strain_incr[_qp];
_inelastic_heat[_qp] = _stress[_qp].doubleContraction(inelastic_strain_incr) / _dt;
// if (_fe_problem.currentlyComputingJacobian())
......
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