ASPECT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
utilities.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2024 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_utilities_h
22 #define _aspect_material_model_utilities_h
23 
24 #include <aspect/global.h>
25 #include <deal.II/base/point.h>
26 #include <deal.II/base/symmetric_tensor.h>
27 #include <deal.II/fe/component_mask.h>
28 #include <deal.II/base/signaling_nan.h>
29 #include <deal.II/base/parameter_handler.h>
30 
31 #include <mpi.h>
32 
33 
34 namespace aspect
35 {
36  template <int dim> class SimulatorAccess;
37  namespace Utilities
38  {
53  using namespace ::Utilities;
54 
55  template <int dim>
56  class StructuredDataLookup;
57  }
58  namespace MaterialModel
59  {
60  template <int dim> class MaterialModelOutputs;
61  template <int dim> struct EquationOfStateOutputs;
62 
68  namespace MaterialUtilities
69  {
70  namespace Lookup
71  {
79  {
80  public:
81 
82  double
83  specific_heat(const double temperature,
84  const double pressure) const;
85 
86  double
87  density(const double temperature,
88  const double pressure) const;
89 
90  double
91  thermal_expansivity(const double temperature,
92  const double pressure) const;
93 
94  double
95  seismic_Vp(const double temperature,
96  const double pressure) const;
97 
98  double
99  seismic_Vs(const double temperature,
100  const double pressure) const;
101 
102  double
103  enthalpy(const double temperature,
104  const double pressure) const;
105 
111  double
112  dHdT (const double temperature,
113  const double pressure) const;
114 
120  double
121  dHdp (const double temperature,
122  const double pressure) const;
123 
136  std::array<std::pair<double, unsigned int>,2>
137  enthalpy_derivatives(const std::vector<double> &temperatures,
138  const std::vector<double> &pressures,
139  const unsigned int n_substeps = 1) const;
140 
141  double
142  dRhodp (const double temperature,
143  const double pressure) const;
144 
149  unsigned int
150  dominant_phase (const double temperature,
151  const double pressure) const;
152 
157  bool
158  has_dominant_phase() const;
159 
164  std::vector<std::string>
165  phase_volume_column_names() const;
166 
171  double
172  phase_volume_fraction(const int phase_id,
173  const double temperature,
174  const double pressure) const;
175 
180  std::array<double,2>
181  get_pT_steps() const;
182 
183 
190  const std::vector<std::string> &
191  get_dominant_phase_names() const;
192 
193  protected:
201  double
202  value (const double temperature,
203  const double pressure,
204  const Table<2, double> &values,
205  const bool interpol) const;
206 
212  unsigned int
213  value (const double temperature,
214  const double pressure,
215  const Table<2, unsigned int> &values) const;
216 
220  double get_nT(const double temperature) const;
221 
225  double get_np(const double pressure) const;
226 
227  ::Table<2,double> density_values;
228  ::Table<2,double> thermal_expansivity_values;
229  ::Table<2,double> specific_heat_values;
230  ::Table<2,double> vp_values;
231  ::Table<2,double> vs_values;
232  ::Table<2,double> enthalpy_values;
233  ::Table<2,unsigned int> dominant_phase_indices;
234 
241  std::vector<std::string> phase_column_names;
242  std::vector<::Table<2,double>> phase_volume_fractions;
243 
244  double delta_press;
245  double min_press;
246  double max_press;
247  double delta_temp;
248  double min_temp;
249  double max_temp;
250  unsigned int n_temperature;
251  unsigned int n_pressure;
252  unsigned int n_phases;
253  unsigned int n_columns;
256  std::vector<std::string> dominant_phase_names;
257  };
258 
264  {
265  public:
266  HeFESToReader(const std::string &material_filename,
267  const std::string &derivatives_filename,
268  const bool interpol,
269  const MPI_Comm comm);
270  };
271 
277  {
278  public:
279  PerplexReader(const std::string &filename,
280  const bool interpol,
281  const MPI_Comm comm);
282  };
283 
289  {
290  public:
291 
295  void
296  initialize(const MPI_Comm comm,
297  const std::string &data_directory,
298  const std::string &material_file_name);
299 
303  double
304  specific_heat(const double entropy,
305  const double pressure) const;
306 
310  double
311  density(const double entropy,
312  const double pressure) const;
313 
317  double
318  thermal_expansivity(const double entropy,
319  const double pressure) const;
320 
324  double
325  temperature(const double entropy,
326  const double pressure) const;
327 
331  double
332  seismic_vp(const double entropy,
333  const double pressure) const;
334 
338  double
339  seismic_vs(const double entropy,
340  const double pressure) const;
341 
345  Tensor<1, 2>
346  density_gradient(const double entropy,
347  const double pressure) const;
348 
349  private:
353  std::unique_ptr<Utilities::StructuredDataLookup<2>> material_lookup;
354  };
355  }
356 
372  std::vector<double>
373  compute_only_composition_fractions(const std::vector<double> &compositional_fields,
374  const std::vector<unsigned int> &indices_to_use);
375 
394  std::vector<double>
395  compute_composition_fractions(const std::vector<double> &compositional_fields,
396  const ComponentMask &field_mask = ComponentMask());
397 
405  std::vector<double>
406  compute_volumes_from_masses(const std::vector<double> &masses,
407  const std::vector<double> &densities,
408  const bool return_as_fraction);
409 
419  {
424  };
425 
426 
427 
434  parse_compositional_averaging_operation (const std::string &parameter_name,
435  const ParameterHandler &prm);
436 
437 
438 
458  double average_value (const std::vector<double> &volume_fractions,
459  const std::vector<double> &parameter_values,
460  const CompositionalAveragingOperation &average_type);
461 
462 
463 
478  template <int dim>
479  void
481  const std::vector<double> &mass_fractions,
482  const std::vector<double> &volume_fractions,
483  const unsigned int i,
485 
486 
487 
491  namespace PhaseUtilities
492  {
499  {
502  };
503  }
504 
521  double phase_average_value (const std::vector<double> &phase_function_values,
522  const std::vector<unsigned int> &n_phase_transitions_per_composition,
523  const std::vector<double> &parameter_values,
524  const unsigned int composition_index,
526 
527 
528 
534  template <int dim>
536  {
541  PhaseFunctionInputs(const double temperature,
542  const double pressure,
543  const double depth,
544  const double pressure_depth_derivative,
545  const unsigned int phase_transition_index);
546 
547  double temperature;
548  double pressure;
549  double depth;
551 
564  };
565 
574  template <int dim>
576  {
577  public:
578 
583  void initialize();
584 
591  double compute_value (const PhaseFunctionInputs<dim> &in) const;
592 
598  double compute_derivative () const;
599 
603  unsigned int n_phase_transitions () const;
604 
608  unsigned int n_phases () const;
609 
613  unsigned int n_phases_over_all_chemical_compositions () const;
614 
618  const std::vector<unsigned int> &
619  n_phase_transitions_for_each_chemical_composition () const;
620 
624  const std::vector<unsigned int> &
625  n_phases_for_each_chemical_composition () const;
626 
633  const std::vector<unsigned int> &
634  n_phase_transitions_for_each_composition () const;
635 
642  const std::vector<unsigned int> &
643  n_phases_for_each_composition () const;
644 
651  static
652  void
653  declare_parameters (ParameterHandler &prm);
654 
661  void
662  parse_parameters (ParameterHandler &prm);
663 
664 
665  private:
669  std::string data_directory;
670 
674  std::vector<std::string> material_file_names;
675 
679  std::vector<double> minimum_temperature;
680 
684  std::vector<double> maximum_temperature;
685 
689  std::vector<double> interval_temperature;
690 
694  std::vector<double> minimum_pressure;
695 
699  std::vector<double> maximum_pressure;
700 
704  std::vector<double> interval_pressure;
705 
710  std::vector<std::unique_ptr<Utilities::StructuredDataLookup<2>>> material_lookup;
711 
718  std::vector<unsigned int> transition_indicators;
719 
723  std::unique_ptr<std::vector<unsigned int>> n_phase_transitions_per_composition;
724 
728  std::vector<unsigned int> n_phases_per_composition;
729 
734 
738  std::vector<unsigned int> n_phases_per_chemical_composition;
739 
743  unsigned int n_phases_total;
744 
749  };
750 
758  template <int dim>
760  {
761  public:
767  double compute_value (const PhaseFunctionInputs<dim> &in) const;
768 
773  double compute_derivative (const PhaseFunctionInputs<dim> &in) const;
774 
778  unsigned int n_phase_transitions () const;
779 
783  unsigned int n_phases () const;
784 
788  unsigned int n_phases_over_all_chemical_compositions () const;
789 
794  double get_transition_slope (const unsigned int phase_transition_index) const;
795 
799  double get_transition_depth (const unsigned int phase_transition_index) const;
800 
804  const std::vector<unsigned int> &
805  n_phase_transitions_for_each_chemical_composition () const;
806 
810  const std::vector<unsigned int> &
811  n_phases_for_each_chemical_composition () const;
812 
819  const std::vector<unsigned int> &
820  n_phase_transitions_for_each_composition () const;
821 
828  const std::vector<unsigned int> &
829  n_phases_for_each_composition () const;
830 
837  static
838  void
839  declare_parameters (ParameterHandler &prm);
840 
847  void
848  parse_parameters (ParameterHandler &prm);
849 
850 
851  private:
856  std::vector<double> transition_depths;
857  std::vector<double> transition_pressures;
858  std::vector<double> transition_temperatures;
859  std::vector<double> transition_widths;
860  std::vector<double> transition_pressure_widths;
861  std::vector<double> transition_slopes;
864 
872 
876  std::unique_ptr<std::vector<unsigned int>> n_phase_transitions_per_composition;
877 
881  std::vector<unsigned int> n_phases_per_composition;
882 
887 
891  std::vector<unsigned int> n_phases_per_chemical_composition;
892 
896  unsigned int n_phases_total;
897 
902  };
903  }
904  }
905 }
906 
907 
908 #endif
void fill_averaged_equation_of_state_outputs(const EquationOfStateOutputs< dim > &eos_outputs, const std::vector< double > &mass_fractions, const std::vector< double > &volume_fractions, const unsigned int i, MaterialModelOutputs< dim > &out)
std::vector< unsigned int > n_phases_per_composition
Definition: utilities.h:881
std::unique_ptr< Utilities::StructuredDataLookup< 2 > > material_lookup
Definition: utilities.h:353
std::vector<::Table< 2, double > > phase_volume_fractions
Definition: utilities.h:242
void declare_parameters(ParameterHandler &prm)
std::unique_ptr< std::vector< unsigned int > > n_phase_transitions_per_composition
Definition: utilities.h:723
std::unique_ptr< std::vector< unsigned int > > n_phase_transitions_per_composition
Definition: utilities.h:876
std::vector< unsigned int > n_phase_transitions_per_chemical_composition
Definition: utilities.h:733
double average_value(const std::vector< double > &volume_fractions, const std::vector< double > &parameter_values, const CompositionalAveragingOperation &average_type)
std::vector< std::unique_ptr< Utilities::StructuredDataLookup< 2 > > > material_lookup
Definition: utilities.h:710
double phase_average_value(const std::vector< double > &phase_function_values, const std::vector< unsigned int > &n_phase_transitions_per_composition, const std::vector< double > &parameter_values, const unsigned int composition_index, const PhaseUtilities::PhaseAveragingOperation operation=PhaseUtilities::arithmetic)
std::vector< unsigned int > n_phase_transitions_per_chemical_composition
Definition: utilities.h:886
CompositionalAveragingOperation parse_compositional_averaging_operation(const std::string &parameter_name, const ParameterHandler &prm)
std::vector< unsigned int > n_phases_per_chemical_composition
Definition: utilities.h:891
std::vector< double > compute_volumes_from_masses(const std::vector< double > &masses, const std::vector< double > &densities, const bool return_as_fraction)
std::vector< double > compute_only_composition_fractions(const std::vector< double > &compositional_fields, const std::vector< unsigned int > &indices_to_use)
std::vector< double > compute_composition_fractions(const std::vector< double > &compositional_fields, const ComponentMask &field_mask=ComponentMask())