Skip to content

Commit

Permalink
Add cutoff pressures to tian parameterization
Browse files Browse the repository at this point in the history
  • Loading branch information
danieldouglas92 committed May 30, 2024
1 parent 065d22b commit 91ccfca
Show file tree
Hide file tree
Showing 6 changed files with 25 additions and 18 deletions.
5 changes: 5 additions & 0 deletions include/aspect/material_model/reactive_fluid_transport.h
Original file line number Diff line number Diff line change
Expand Up @@ -179,6 +179,11 @@ namespace aspect
std::vector<double> csat_sediment_poly_coeffs {-0.150662, 0.301807, 1.01867};
std::vector<double> Td_sediment_poly_coeffs {2.83277, -24.7593, 85.9090, 524.898};

// The polynomials breakdown above certain pressures, 10 GPa for peridotite, 26 GPa for gabbro, 16 GPa for MORB,
// and 50 GPa for sediment. These cutoff pressures were determined by extending the pressure range in Tian et al. (2019)
// and observing where the maximum allowed water contents jump towards infinite values.
const std::vector<double> pressure_cutoffs {10, 26, 16, 50};

std::vector<std::vector<double>> devolatilization_enthalpy_changes {LR_peridotite_poly_coeffs, LR_gabbro_poly_coeffs, \
LR_MORB_poly_coeffs, LR_sediment_poly_coeffs
};
Expand Down
10 changes: 6 additions & 4 deletions source/material_model/reactive_fluid_transport.cc
Original file line number Diff line number Diff line change
Expand Up @@ -56,10 +56,6 @@ namespace aspect
tian_equilibrium_bound_water_content (const MaterialModel::MaterialModelInputs<dim> &in,
unsigned int q) const
{
// Pressure, which must be in GPa for the parametrization, or GPa^-1
const double pressure = in.pressure[q]<=0 ? 1e-12 : in.pressure[q]/1.e9;
const double inverse_pressure = std::pow(pressure, -1);

// Create arrays that will store the values of the polynomials at the current pressure
std::vector<double> LR_values(4);
std::vector<double> csat_values(4);
Expand All @@ -71,6 +67,12 @@ namespace aspect
// Td polynomials are defined in equations 15, B3, B11, and B19.
for (unsigned int i = 0; i<devolatilization_enthalpy_changes.size(); ++i)
{
// Pressure, which must be in GPa for the parametrization, or GPa^-1. The polynomials for each lithology
// breaks down above certain pressures, make sure that we cap the pressure just before this break down.
// Introduce minimum pressure to avoid a division by 0.
const double minimum_pressure = 1e-12;
const double pressure = std::min(std::max(minimum_pressure, in.pressure[q]/1.e9), pressure_cutoffs[i]);
const double inverse_pressure = 1.0/pressure;
for (unsigned int j = 0; j<devolatilization_enthalpy_changes[i].size(); ++j)
{
LR_values[i] += devolatilization_enthalpy_changes[i][j] * std::pow(inverse_pressure, devolatilization_enthalpy_changes[i].size() - 1 - j);
Expand Down
6 changes: 3 additions & 3 deletions tests/tian_MORB/screen-output
Original file line number Diff line number Diff line change
Expand Up @@ -31,8 +31,8 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 0, 8.51473e-17, 0, 4.64648e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64648e-10
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 0, 8.51473e-17, 0, 4.64647e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64647e-10


Postprocessing:
Expand Down Expand Up @@ -76,7 +76,7 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.5218e-12, 1.52682e-10, 5.02896e-12, 0, 0, 1.07956e-16, 0, 1.87227e-09
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.5218e-12, 1.52682e-10, 5.02898e-12, 0, 0, 1.07956e-16, 0, 1.87227e-09
Relative nonlinear residual (total system) after nonlinear iteration 3: 1.87227e-09


Expand Down
8 changes: 4 additions & 4 deletions tests/tian_gabbro/screen-output
Original file line number Diff line number Diff line change
Expand Up @@ -31,8 +31,8 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 8.51473e-17, 0, 0, 4.64648e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64648e-10
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 8.51473e-17, 0, 0, 4.64647e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64647e-10


Postprocessing:
Expand Down Expand Up @@ -76,8 +76,8 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.771e-12, 1.64474e-10, 4.67296e-12, 0, 9.75727e-17, 0, 0, 2.23762e-09
Relative nonlinear residual (total system) after nonlinear iteration 3: 2.23762e-09
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.77098e-12, 1.64474e-10, 4.67294e-12, 0, 9.75727e-17, 0, 0, 2.23761e-09
Relative nonlinear residual (total system) after nonlinear iteration 3: 2.23761e-09


Postprocessing:
Expand Down
8 changes: 4 additions & 4 deletions tests/tian_peridotite/screen-output
Original file line number Diff line number Diff line change
Expand Up @@ -31,8 +31,8 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 8.51473e-17, 0, 0, 0, 4.64648e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64648e-10
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 8.51473e-17, 0, 0, 0, 4.64647e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64647e-10


Postprocessing:
Expand Down Expand Up @@ -63,7 +63,7 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 7+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.181e-05, 0.000574823, 0.000145666, 1.02675e-16, 0, 0, 0, 0.000630417
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 2.181e-05, 0.000574823, 0.000145666, 9.16123e-17, 0, 0, 0, 0.000630417
Relative nonlinear residual (total system) after nonlinear iteration 2: 0.000630417

Solving temperature system... 4 iterations.
Expand All @@ -76,7 +76,7 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 4+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 4.66691e-08, 2.35667e-06, 3.02875e-07, 1.25585e-16, 0, 0, 0, 2.19117e-06
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 4.66691e-08, 2.35667e-06, 3.02875e-07, 1.22071e-16, 0, 0, 0, 2.19117e-06
Relative nonlinear residual (total system) after nonlinear iteration 3: 2.35667e-06


Expand Down
6 changes: 3 additions & 3 deletions tests/tian_sediment/screen-output
Original file line number Diff line number Diff line change
Expand Up @@ -31,8 +31,8 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 0, 0, 8.51473e-17, 4.64648e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64648e-10
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.54103e-16, 0, 1.39664e-16, 0, 0, 0, 8.51473e-17, 4.64647e-10
Relative nonlinear residual (total system) after nonlinear iteration 2: 4.64647e-10


Postprocessing:
Expand Down Expand Up @@ -76,7 +76,7 @@ Number of degrees of freedom: 989 (162+73+162+25+81+81+81+81+81+81+81)
Rebuilding Stokes preconditioner...
Solving Stokes system... 0+0 iterations.
Solving fluid velocity system... 1 iterations.
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.77101e-12, 1.37697e-10, 2.93736e-12, 0, 0, 0, 1.04501e-16, 1.53067e-09
Relative nonlinear residuals (temperature, compositional fields, Stokes system): 1.77098e-12, 1.37697e-10, 2.93735e-12, 0, 0, 0, 1.04501e-16, 1.53067e-09
Relative nonlinear residual (total system) after nonlinear iteration 3: 1.53067e-09


Expand Down

0 comments on commit 91ccfca

Please sign in to comment.