ASPECT
peierls_creep.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2020 - 2021 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_material_model_rheology_peierls_creep_h
22 #define _aspect_material_model_rheology_peierls_creep_h
23 
24 #include <aspect/global.h>
26 
27 namespace aspect
28 {
29  namespace MaterialModel
30  {
31  using namespace dealii;
32 
33  namespace Rheology
34  {
39  {
45  double prefactor;
53  double stress_cutoff;
54  };
55 
64  template <int dim>
66  {
67  public:
71  PeierlsCreep();
72 
77  compute_creep_parameters (const unsigned int composition,
78  const std::vector<double> &phase_function_values = std::vector<double>(),
79  const std::vector<unsigned int> &n_phases_per_composition = std::vector<unsigned int>()) const;
80 
84  static
85  void
86  declare_parameters (ParameterHandler &prm);
87 
95  void
96  parse_parameters (ParameterHandler &prm,
97  const std::unique_ptr<std::vector<unsigned int>> &expected_n_phases_per_composition = nullptr);
98 
102  double
103  compute_approximate_viscosity (const double strain_rate,
104  const double pressure,
105  const double temperature,
106  const unsigned int composition,
107  const std::vector<double> &phase_function_values = std::vector<double>(),
108  const std::vector<unsigned int> &n_phases_per_composition = std::vector<unsigned int>()) const;
109 
113  double
114  compute_exact_viscosity (const double strain_rate,
115  const double pressure,
116  const double temperature,
117  const unsigned int composition,
118  const std::vector<double> &phase_function_values = std::vector<double>(),
119  const std::vector<unsigned int> &n_phases_per_composition = std::vector<unsigned int>()) const;
120 
126  double
127  compute_viscosity (const double strain_rate,
128  const double pressure,
129  const double temperature,
130  const unsigned int composition,
131  const std::vector<double> &phase_function_values = std::vector<double>(),
132  const std::vector<unsigned int> &n_phases_per_composition = std::vector<unsigned int>()) const;
133 
138  std::pair<double, double>
139  compute_approximate_strain_rate_and_derivative (const double stress,
140  const double pressure,
141  const double temperature,
142  const PeierlsCreepParameters creep_parameters) const;
143 
148  std::pair<double, double>
149  compute_exact_strain_rate_and_derivative (const double stress,
150  const double pressure,
151  const double temperature,
152  const PeierlsCreepParameters creep_parameters) const;
153 
161  std::pair<double, double>
162  compute_strain_rate_and_derivative (const double stress,
163  const double pressure,
164  const double temperature,
165  const PeierlsCreepParameters creep_parameters) const;
166 
167  private:
181  {
183  exact
184  } peierls_creep_flow_law;
185 
189  std::vector<double> prefactors;
190 
194  std::vector<double> stress_exponents;
195 
199  std::vector<double> activation_energies;
200 
204  std::vector<double> activation_volumes;
205 
209  std::vector<double> peierls_stresses;
210 
214  std::vector<double> fitting_parameters;
215 
220  std::vector<double> glide_parameters_p;
221 
226  std::vector<double> glide_parameters_q;
227 
228  std::vector<double> stress_cutoffs;
229 
236 
237  };
238  }
239  }
240 }
241 #endif
void declare_parameters(ParameterHandler &prm)
Definition: compat.h:88