Skip to content

Commit

Permalink
Merge pull request #5991 from gassmoeller/grain_size_reaction_model
Browse files Browse the repository at this point in the history
Create grain size evolution reaction model
  • Loading branch information
bobmyhill authored Jul 29, 2024
2 parents 6e5f0cb + c2058e2 commit 44e61b3
Show file tree
Hide file tree
Showing 8 changed files with 1,100 additions and 822 deletions.
5 changes: 5 additions & 0 deletions doc/modules/changes/20240726_gassmoeller
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
New: The grain size evolution functionality has been extracted from the grain
size material model into the new grain size evolution reaction model. This
simplifies using the grain size evolution equations in other material models.
<br>
(Rene Gassmoeller, 2024/07/26)
175 changes: 24 additions & 151 deletions include/aspect/material_model/grain_size.h
Original file line number Diff line number Diff line change
Expand Up @@ -19,16 +19,17 @@
*/


#ifndef _aspect_model_grain_size_h
#define _aspect_model_grain_size_h
#ifndef _aspect_material_model_grain_size_h
#define _aspect_material_model_grain_size_h

#include <aspect/material_model/interface.h>
#include <aspect/material_model/reaction_model/grain_size_evolution.h>

#include <aspect/material_model/utilities.h>
#include <aspect/simulator_access.h>
#include <aspect/material_model/rheology/drucker_prager.h>

#include <deal.II/matrix_free/fe_point_evaluation.h>
#include <deal.II/sundials/arkode.h>

#include <array>

Expand Down Expand Up @@ -74,8 +75,8 @@ namespace aspect
* determine the grain size evolution in dependence of the strain rate,
* temperature, phase transitions, and the creep regime.
* This material model only works if a compositional field
* named 'grain_size' is present. In the diffusion
* creep regime, the viscosity depends on this grain size. We use the grain
* with the name 'grain_size' is present.
* The diffusion creep viscosity depends on this grain size. We use the grain
* size evolution laws described in Behn et al., 2009. Implications of grain
* size evolution on the seismic structure of the oceanic upper mantle, Earth
* Planet. Sci. Letters, 282, 178–189. Other material parameters are either
Expand Down Expand Up @@ -110,6 +111,9 @@ namespace aspect
*/
bool is_compressible () const override;

/**
* Evaluate the material model at the given input points.
*/
void evaluate(const typename Interface<dim>::MaterialModelInputs &in,
typename Interface<dim>::MaterialModelOutputs &out) const override;

Expand All @@ -130,14 +134,13 @@ namespace aspect
void
parse_parameters (ParameterHandler &prm) override;


void
create_additional_named_outputs (MaterialModel::MaterialModelOutputs<dim> &out) const override;

/**
* @}
*/

void
create_additional_named_outputs (MaterialModel::MaterialModelOutputs<dim> &out) const override;

/**
* Returns the enthalpy as calculated by HeFESTo.
*/
Expand All @@ -160,7 +163,7 @@ namespace aspect
std::array<std::pair<double, unsigned int>,2>
enthalpy_derivative (const typename Interface<dim>::MaterialModelInputs &in) const;

protected:
private:
double reference_rho;
double reference_T;
double eta;
Expand All @@ -177,125 +180,6 @@ namespace aspect
*/
double k_value;

/**
* The index of the compositional field that represents the grain size.
*/
unsigned int grain_size_index;

/**
* Parameters controlling the grain size evolution.
*/
std::vector<double> grain_growth_activation_energy;
std::vector<double> grain_growth_activation_volume;
std::vector<double> grain_growth_rate_constant;
std::vector<double> grain_growth_exponent;
double minimum_grain_size;
std::vector<double> reciprocal_required_strain;
std::vector<double> recrystallized_grain_size;

/**
* Parameters controlling the dynamic grain recrystallization.
*/
std::vector<double> grain_boundary_energy;
std::vector<double> boundary_area_change_work_fraction;
std::vector<double> geometric_constant;

/**
* A struct that contains information about which
* formulation of grain size evolution should be used.
*/
struct Formulation
{
/**
* This enum lists available options that
* determine the equations being used for grain size evolution.
*
* We currently support three approaches:
*
* 'paleowattmeter':
* Austin, N. J., & Evans, B. (2007). Paleowattmeters: A scaling
* relation for dynamically recrystallized grain size. Geology, 354), 343-346.).
*
* 'paleopiezometer':
* Hall, C. E., Parmentier, E. M. (2003). Influence of grain size
* evolution on convective instability. Geochemistry, Geophysics,
* Geosystems, 4(3).
*
* 'pinned_grain_damage':
* Mulyukova, E., & Bercovici, D. (2018). Collapse of passive margins
* by lithospheric damage and plunging grain size. Earth and Planetary
* Science Letters, 484, 341-352.
*/
enum Kind
{
paleowattmeter,
paleopiezometer,
pinned_grain_damage
};

/**
* This function translates an input string into the
* available enum options.
*/
static
Kind
parse(const std::string &input)
{
if (input == "paleowattmeter")
return Formulation::paleowattmeter;
else if (input == "paleopiezometer")
return Formulation::paleopiezometer;
else if (input == "pinned grain damage")
return Formulation::pinned_grain_damage;
else
AssertThrow(false, ExcNotImplemented());

return Formulation::Kind();
}
};

/**
* A variable that records the formulation of how to evolve grain size.
* See the type of this variable for a description of available options.
*/
typename Formulation::Kind grain_size_evolution_formulation;

/**
* This function returns the fraction of shear heating energy partitioned
* into grain damage using the implementation by Mulyukova and Bercovici (2018)
* Collapse of passive margins by lithospheric damage
* and plunging grain size. Earth and Planetary Science Letters, 484, 341-352.
*/
double compute_partitioning_fraction (const double temperature) const;

/**
* Parameters controlling the partitioning of energy
* into grain damage in the pinned state.
*/
double grain_size_reduction_work_fraction_exponent;
double minimum_grain_size_reduction_work_fraction;
double maximum_grain_size_reduction_work_fraction;
double temperature_minimum_partitioning_power;
double temperature_maximum_partitioning_power;

/**
* Functions and parameters controlling conversion from interface roughness to grain size,
* used in pinned state formulation of grain damage. This conversion depends on the
* proportion of the two mineral phases.
*
* A detailed description of this approach can be found in Appendix H.1, in Equation (8) in
* the main manuscript, and in equation (F.28) of Bercovici, David, and Yanick Ricard (2012).
* Mechanisms for the generation of plate tectonics by two-phase grain-damage
* and pinning. Physics of the Earth and Planetary Interiors 202 (2012): 27-55.
*/
double phase_distribution;

/**
* The factor used to convert roughness into the equivalent mean grain
* size for a given volume fraction of a mineral in the two-phase damage model.
*/
double roughness_to_grain_size;

/**
* Parameters controlling the viscosity.
*/
Expand Down Expand Up @@ -324,8 +208,12 @@ namespace aspect
double min_thermal_expansivity;
double max_thermal_expansivity;
unsigned int max_latent_heat_substeps;
double min_grain_size;
double minimum_grain_size;

/**
* Calculate the diffusion viscosity in dependence of temperature,
* pressure, grain size, and phase.
*/
double diffusion_viscosity (const double temperature,
const double adiabatic_temperature,
const double adiabatic_pressure,
Expand Down Expand Up @@ -388,25 +276,6 @@ namespace aspect
const std::vector<double> &compositional_fields,
const Point<dim> &position) const;

/**
* Rate of grain size growth (Ostwald ripening) or reduction
* (due to dynamic recrystallization and phase transformations)
* in dependence on temperature, pressure, strain rate, mineral
* phase and creep regime.
* We use the grain size growth laws as for example described
* in Behn, M. D., Hirth, G., & Elsenbeck, J. R. (2009). Implications
* of grain size evolution on the seismic structure of the oceanic
* upper mantle. Earth and Planetary Science Letters, 282(1), 178-189.
*
* For the rate of grain size reduction due to dynamic crystallization
* there is the choice between the paleowattmeter (Austins and
* Evans, 2007) and the paleopiezometer (Hall and Parmentier, 2003)
* as described in the parameter use_paleowattmeter.
*/
std::vector<std::vector<double>>
grain_size_change (const typename Interface<dim>::MaterialModelInputs &in,
const std::vector<double> &adiabatic_pressure,
const std::vector<unsigned int> &phase_indices) const;

/**
* Function that returns the phase for a given
Expand All @@ -431,7 +300,7 @@ namespace aspect
* transition in the model, given the temperature, pressure, depth,
* and density gradient.
*/
MaterialUtilities::PhaseFunction<dim> phase_function;
std::shared_ptr<MaterialUtilities::PhaseFunction<dim>> phase_function;


/**
Expand Down Expand Up @@ -473,7 +342,6 @@ namespace aspect
mutable std::unique_ptr<FEPointEvaluation<1, dim>> temperature_evaluator;
mutable std::unique_ptr<FEPointEvaluation<1, dim>> pressure_evaluator;

private:
/*
* Object for computing plastic stresses, viscosities, and additional outputs,
* as well as an object for the required input parameters.
Expand All @@ -482,6 +350,11 @@ namespace aspect
bool use_adiabatic_pressure_for_yielding;
Rheology::DruckerPrager<dim> drucker_prager_plasticity;
Rheology::DruckerPragerParameters drucker_prager_parameters;

/**
* The reaction model that calculates the grain size change over time.
*/
std::unique_ptr<ReactionModel::GrainSizeEvolution<dim>> grain_size_evolution;
};

}
Expand Down
Loading

0 comments on commit 44e61b3

Please sign in to comment.