Skip to content

Commit e544695

Browse files
gassmoellergassmoeller
committed
Add shear heating plugin
1 parent 6b74606 commit e544695

File tree

4 files changed

+190
-226
lines changed

4 files changed

+190
-226
lines changed

‎cookbooks/anisotropic_viscosity/av_material.cc‎

Lines changed: 0 additions & 113 deletions
Original file line numberDiff line numberDiff line change
@@ -48,8 +48,6 @@
4848
#include <deal.II/base/geometry_info.h>
4949

5050
#include <aspect/material_model/simple.h>
51-
#include <aspect/heating_model/interface.h>
52-
#include <aspect/heating_model/shear_heating.h>
5351
#include <aspect/gravity_model/interface.h>
5452
#include <aspect/simulator/assemblers/stokes.h>
5553
#include <aspect/simulator_signals.h>
@@ -60,103 +58,6 @@
6058

6159
namespace aspect
6260
{
63-
namespace HeatingModel
64-
{
65-
template <int dim>
66-
class ShearHeatingAnisotropicViscosity : public Interface<dim>, public ::aspect::SimulatorAccess<dim>
67-
{
68-
public:
69-
/**
70-
* Compute the heating model outputs for this class.
71-
*/
72-
void
73-
evaluate (const MaterialModel::MaterialModelInputs<dim> &material_model_inputs,
74-
const MaterialModel::MaterialModelOutputs<dim> &material_model_outputs,
75-
HeatingModel::HeatingModelOutputs &heating_model_outputs) const override;
76-
77-
/**
78-
* Allow the heating model to attach additional material model outputs.
79-
*/
80-
void
81-
create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &material_model_outputs) const override;
82-
};
83-
84-
85-
86-
template <int dim>
87-
void
88-
ShearHeatingAnisotropicViscosity<dim>::
89-
evaluate (const MaterialModel::MaterialModelInputs<dim> &material_model_inputs,
90-
const MaterialModel::MaterialModelOutputs<dim> &material_model_outputs,
91-
HeatingModel::HeatingModelOutputs &heating_model_outputs) const
92-
{
93-
Assert(heating_model_outputs.heating_source_terms.size() == material_model_inputs.position.size(),
94-
ExcMessage ("Heating outputs need to have the same number of entries as the material model inputs."));
95-
96-
// Some material models provide dislocation viscosities and boundary area work fractions
97-
// as additional material outputs. If they are attached, use them.
98-
const std::shared_ptr<const ShearHeatingOutputs<dim>> shear_heating_out
99-
= material_model_outputs.template get_additional_output_object<ShearHeatingOutputs<dim>>();
100-
101-
const std::shared_ptr<const MaterialModel::AnisotropicViscosity<dim>> anisotropic_viscosity
102-
= material_model_outputs.template get_additional_output_object<MaterialModel::AnisotropicViscosity<dim>>();
103-
104-
for (unsigned int q=0; q<heating_model_outputs.heating_source_terms.size(); ++q)
105-
{
106-
// If there is an anisotropic viscosity, use it to compute the correct stress
107-
const SymmetricTensor<2,dim> &directed_strain_rate = ((anisotropic_viscosity != nullptr)
108-
?
109-
anisotropic_viscosity->stress_strain_directors[q]
110-
* material_model_inputs.strain_rate[q]
111-
:
112-
material_model_inputs.strain_rate[q]);
113-
114-
const SymmetricTensor<2,dim> stress =
115-
2 * material_model_outputs.viscosities[q] *
116-
(this->get_material_model().is_compressible()
117-
?
118-
directed_strain_rate - 1./3. * trace(directed_strain_rate) * unit_symmetric_tensor<dim>()
119-
:
120-
directed_strain_rate);
121-
122-
const SymmetricTensor<2,dim> deviatoric_strain_rate =
123-
(this->get_material_model().is_compressible()
124-
?
125-
material_model_inputs.strain_rate[q]
126-
- 1./3. * trace(material_model_inputs.strain_rate[q]) * unit_symmetric_tensor<dim>()
127-
:
128-
material_model_inputs.strain_rate[q]);
129-
130-
heating_model_outputs.heating_source_terms[q] = stress * deviatoric_strain_rate;
131-
132-
// If shear heating work fractions are provided, reduce the
133-
// overall heating by this amount (which is assumed to be converted into other forms of energy)
134-
if (shear_heating_out != nullptr)
135-
heating_model_outputs.heating_source_terms[q] *= shear_heating_out->shear_heating_work_fractions[q];
136-
137-
heating_model_outputs.lhs_latent_heat_terms[q] = 0.0;
138-
}
139-
}
140-
141-
142-
143-
template <int dim>
144-
void
145-
ShearHeatingAnisotropicViscosity<dim>::
146-
create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &material_model_outputs) const
147-
{
148-
const unsigned int n_points = material_model_outputs.viscosities.size();
149-
150-
if (material_model_outputs.template has_additional_output_object<MaterialModel::AnisotropicViscosity<dim>>() == false)
151-
{
152-
material_model_outputs.additional_outputs.push_back(
153-
std::make_unique<MaterialModel::AnisotropicViscosity<dim>> (n_points));
154-
}
155-
156-
this->get_material_model().create_additional_named_outputs(material_model_outputs);
157-
}
158-
}
159-
16061
namespace MaterialModel
16162
{
16263
// The AV material model calculates an anisotropic viscosity tensor from director vectors and the normal and shear
@@ -481,20 +382,6 @@ namespace aspect
481382
// explicit instantiations
482383
namespace aspect
483384
{
484-
namespace HeatingModel
485-
{
486-
ASPECT_REGISTER_HEATING_MODEL(ShearHeatingAnisotropicViscosity,
487-
"anisotropic shear heating",
488-
"Implementation of a standard model for shear heating. "
489-
"Adds the term: "
490-
"$ 2 \\eta \\left( \\varepsilon - \\frac{1}{3} \\text{tr} "
491-
"\\varepsilon \\mathbf 1 \\right) : \\left( \\varepsilon - \\frac{1}{3} "
492-
"\\text{tr} \\varepsilon \\mathbf 1 \\right)$ to the "
493-
"right-hand side of the temperature equation.")
494-
}
495-
496-
497-
498385
namespace MaterialModel
499386
{
500387
ASPECT_REGISTER_MATERIAL_MODEL(AV,

‎include/aspect/heating_model/shear_heating_anisotropic_viscosity.h‎

Lines changed: 61 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,61 @@
1+
/*
2+
Copyright (C) 2025 by the authors of the ASPECT code.
3+
4+
This file is part of ASPECT.
5+
6+
ASPECT is free software; you can redistribute it and/or modify
7+
it under the terms of the GNU General Public License as published by
8+
the Free Software Foundation; either version 2, or (at your option)
9+
any later version.
10+
11+
ASPECT is distributed in the hope that it will be useful,
12+
but WITHOUT ANY WARRANTY; without even the implied warranty of
13+
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14+
GNU General Public License for more details.
15+
16+
You should have received a copy of the GNU General Public License
17+
along with ASPECT; see the file LICENSE. If not see
18+
<http://www.gnu.org/licenses/>.
19+
*/
20+
21+
#ifndef _aspect_heating_model_shear_heating_anisotropic_viscosity_h
22+
#define _aspect_heating_model_shear_heating_anisotropic_viscosity_h
23+
24+
#include <aspect/heating_model/interface.h>
25+
26+
#include <aspect/simulator_access.h>
27+
28+
namespace aspect
29+
{
30+
namespace HeatingModel
31+
{
32+
/**
33+
* A class that implements a standard model for shear heating extended for an
34+
* anisotropic viscosity tensor. If the material model provides a stress-
35+
* strain director tensor, then the strain-rate is multiplied with this
36+
* tensor to compute the stress that is used when computing the shear heating.
37+
*
38+
* @ingroup HeatingModels
39+
*/
40+
template <int dim>
41+
class ShearHeatingAnisotropicViscosity : public Interface<dim>, public ::aspect::SimulatorAccess<dim>
42+
{
43+
public:
44+
/**
45+
* Compute the heating model outputs for this class.
46+
*/
47+
void
48+
evaluate (const MaterialModel::MaterialModelInputs<dim> &material_model_inputs,
49+
const MaterialModel::MaterialModelOutputs<dim> &material_model_outputs,
50+
HeatingModel::HeatingModelOutputs &heating_model_outputs) const override;
51+
52+
/**
53+
* Allow the heating model to attach additional material model outputs.
54+
*/
55+
void
56+
create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &material_model_outputs) const override;
57+
};
58+
}
59+
}
60+
61+
#endif

‎source/heating_model/shear_heating_anisotropic_viscosity.cc‎

Lines changed: 123 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,123 @@
1+
/*
2+
Copyright (C) 2015 - 2023 by the authors of the ASPECT code.
3+
4+
This file is part of ASPECT.
5+
6+
ASPECT is free software; you can redistribute it and/or modify
7+
it under the terms of the GNU General Public License as published by
8+
the Free Software Foundation; either version 2, or (at your option)
9+
any later version.
10+
11+
ASPECT is distributed in the hope that it will be useful,
12+
but WITHOUT ANY WARRANTY; without even the implied warranty of
13+
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14+
GNU General Public License for more details.
15+
16+
You should have received a copy of the GNU General Public License
17+
along with ASPECT; see the file LICENSE. If not see
18+
<http://www.gnu.org/licenses/>.
19+
*/
20+
21+
#include <aspect/heating_model/shear_heating_anisotropic_viscosity.h>
22+
23+
#include <aspect/material_model/interface.h>
24+
#include <aspect/heating_model/shear_heating.h>
25+
#include <aspect/simulator/assemblers/stokes_anisotropic_viscosity.h>
26+
27+
#include <deal.II/base/symmetric_tensor.h>
28+
#include <deal.II/base/signaling_nan.h>
29+
30+
namespace aspect
31+
{
32+
namespace HeatingModel
33+
{
34+
template <int dim>
35+
void
36+
ShearHeatingAnisotropicViscosity<dim>::
37+
evaluate (const MaterialModel::MaterialModelInputs<dim> &material_model_inputs,
38+
const MaterialModel::MaterialModelOutputs<dim> &material_model_outputs,
39+
HeatingModel::HeatingModelOutputs &heating_model_outputs) const
40+
{
41+
Assert(heating_model_outputs.heating_source_terms.size() == material_model_inputs.position.size(),
42+
ExcMessage ("Heating outputs need to have the same number of entries as the material model inputs."));
43+
44+
// Some material models provide dislocation viscosities and boundary area work fractions
45+
// as additional material outputs. If they are attached, use them.
46+
const std::shared_ptr<const ShearHeatingOutputs<dim>> shear_heating_out =
47+
material_model_outputs.template get_additional_output_object<ShearHeatingOutputs<dim>>();
48+
49+
const std::shared_ptr<const MaterialModel::AnisotropicViscosity<dim>> anisotropic_viscosity =
50+
material_model_outputs.template get_additional_output_object<MaterialModel::AnisotropicViscosity<dim>>();
51+
52+
for (unsigned int q=0; q<heating_model_outputs.heating_source_terms.size(); ++q)
53+
{
54+
// If there is an anisotropic viscosity, use it to compute the correct stress
55+
const SymmetricTensor<2,dim> &directed_strain_rate = ((anisotropic_viscosity != nullptr)
56+
?
57+
anisotropic_viscosity->stress_strain_directors[q]
58+
* material_model_inputs.strain_rate[q]
59+
:
60+
material_model_inputs.strain_rate[q]);
61+
62+
const SymmetricTensor<2,dim> stress =
63+
2 * material_model_outputs.viscosities[q] *
64+
(this->get_material_model().is_compressible()
65+
?
66+
directed_strain_rate - 1./3. * trace(directed_strain_rate) * unit_symmetric_tensor<dim>()
67+
:
68+
directed_strain_rate);
69+
70+
const SymmetricTensor<2,dim> deviatoric_strain_rate =
71+
(this->get_material_model().is_compressible()
72+
?
73+
material_model_inputs.strain_rate[q]
74+
- 1./3. * trace(material_model_inputs.strain_rate[q]) * unit_symmetric_tensor<dim>()
75+
:
76+
material_model_inputs.strain_rate[q]);
77+
78+
heating_model_outputs.heating_source_terms[q] = stress * deviatoric_strain_rate;
79+
80+
// If shear heating work fractions are provided, reduce the
81+
// overall heating by this amount (which is assumed to be converted into other forms of energy)
82+
if (shear_heating_out != nullptr)
83+
heating_model_outputs.heating_source_terms[q] *= shear_heating_out->shear_heating_work_fractions[q];
84+
85+
heating_model_outputs.lhs_latent_heat_terms[q] = 0.0;
86+
}
87+
}
88+
89+
90+
91+
template <int dim>
92+
void
93+
ShearHeatingAnisotropicViscosity<dim>::
94+
create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &material_model_outputs) const
95+
{
96+
const unsigned int n_points = material_model_outputs.viscosities.size();
97+
98+
if (material_model_outputs.template has_additional_output_object<MaterialModel::AnisotropicViscosity<dim>>() == false)
99+
{
100+
material_model_outputs.additional_outputs.push_back(
101+
std::make_unique<MaterialModel::AnisotropicViscosity<dim>> (n_points));
102+
}
103+
104+
this->get_material_model().create_additional_named_outputs(material_model_outputs);
105+
}
106+
}
107+
}
108+
109+
110+
111+
// explicit instantiations
112+
namespace aspect
113+
{
114+
namespace HeatingModel
115+
{
116+
ASPECT_REGISTER_HEATING_MODEL(ShearHeatingAnisotropicViscosity,
117+
"anisotropic shear heating",
118+
"Implementation of a standard model for shear heating extended for an "
119+
"anisotropic viscosity tensor. If the material model provides a stress-"
120+
"strain director tensor, then the strain-rate is multiplied with this "
121+
"tensor to compute the stress that is used when computing the shear heating.")
122+
}
123+
}

0 commit comments

Comments
 (0)