ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2018 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_interface_h
22 #define _aspect_material_model_interface_h
23 
24 #include <aspect/global.h>
25 #include <aspect/plugins.h>
27 
28 #include <deal.II/base/point.h>
34 #include <deal.II/fe/mapping.h>
35 #include <deal.II/fe/fe_values.h>
39 
40 namespace aspect
41 {
42  template <int dim>
43  struct Introspection;
44 
45 
53  namespace MaterialModel
54  {
55  using namespace dealii;
56 
61  namespace NonlinearDependence
62  {
85  {
87 
88  none = 1,
90  pressure = 4,
93 
95  };
96 
97 
102  const Dependence d2)
103  {
104  return Dependence((int)d1 | (int)d2);
105  }
106 
108  const Dependence d2)
109  {
110  d1 = (d1 | d2);
111  return d1;
112  }
113 
119  {
125 
131 
137 
143 
149 
155  ModelDependence ();
156  };
157 
166  bool
167  identifies_single_variable(const Dependence dependence);
168  }
169 
170 
171  template <int dim> class AdditionalMaterialInputs;
172 
173 
180  template <int dim>
182  {
193  MaterialModelInputs(const unsigned int n_points,
194  const unsigned int n_comp);
195 
205  const Introspection<dim> &introspection,
206  const bool use_strain_rate = true);
207 
208 
222  const typename DoFHandler<dim>::active_cell_iterator &cell,
223  const Introspection<dim> &introspection,
224  const LinearAlgebra::BlockVector &solution_vector,
225  const bool use_strain_rates = true);
226 
231 
236  void reinit(const FEValuesBase<dim,dim> &fe_values,
237  const typename DoFHandler<dim>::active_cell_iterator &cell,
238  const Introspection<dim> &introspection,
239  const LinearAlgebra::BlockVector &solution_vector,
240  const bool use_strain_rates = true);
241 
242 
247  std::vector<Point<dim> > position;
248 
252  std::vector<double> temperature;
253 
257  std::vector<double> pressure;
258 
263  std::vector<Tensor<1,dim> > pressure_gradient;
264 
272  std::vector<Tensor<1,dim> > velocity;
273 
279  std::vector<std::vector<double> > composition;
280 
293  std::vector<SymmetricTensor<2,dim> > strain_rate;
294 
308 
325 
331  std::vector<std::shared_ptr<AdditionalMaterialInputs<dim> > > additional_inputs;
332 
339  template <class AdditionalInputType>
340  AdditionalInputType *get_additional_input();
341 
345  template <class AdditionalInputType>
346  const AdditionalInputType *get_additional_input() const;
347 
348  private:
355  MaterialModelInputs &operator=(const MaterialModelInputs &material);
356  };
357 
358 
359  template <int dim> class AdditionalMaterialOutputs;
360 
361 
368  template <int dim>
370  {
381  MaterialModelOutputs (const unsigned int n_points,
382  const unsigned int n_comp);
383 
387  std::vector<double> viscosities;
388 
392  std::vector<double> densities;
393 
398  std::vector<double> thermal_expansion_coefficients;
399 
403  std::vector<double> specific_heat;
404 
408  std::vector<double> thermal_conductivities;
409 
414  std::vector<double> compressibilities;
415 
421  std::vector<double> entropy_derivative_pressure;
422 
429  std::vector<double> entropy_derivative_temperature;
430 
467  std::vector<std::vector<double> > reaction_terms;
468 
474  std::vector<std::shared_ptr<AdditionalMaterialOutputs<dim> > > additional_outputs;
475 
483  template <class AdditionalOutputType>
484  AdditionalOutputType *get_additional_output();
485 
489  template <class AdditionalOutputType>
490  const AdditionalOutputType *get_additional_output() const;
491  };
492 
493 
508  namespace MaterialAveraging
509  {
545  {
553  };
554 
555 
562  std::string get_averaging_operation_names ();
563 
569  AveragingOperation parse_averaging_operation_name (const std::string &s);
570 
576  template <int dim>
577  void average (const AveragingOperation operation,
578  const typename DoFHandler<dim>::active_cell_iterator &cell,
579  const Quadrature<dim> &quadrature_formula,
580  const Mapping<dim> &mapping,
581  MaterialModelOutputs<dim> &values_out);
582 
588  void average_property (const AveragingOperation operation,
589  const FullMatrix<double> &projection_matrix,
590  const FullMatrix<double> &expansion_matrix,
591  std::vector<double> &values_out);
592  }
593 
594 
613  template <int dim>
615  {
616  public:
618  {}
619 
624  virtual void
625  fill (const LinearAlgebra::BlockVector &solution,
626  const FEValuesBase<dim> &fe_values,
627  const Introspection<dim> &introspection) = 0;
628  };
629 
630 
657  template <int dim>
659  {
660  public:
662  {}
663 
664  virtual void average (const MaterialAveraging::AveragingOperation /*operation*/,
665  const FullMatrix<double> &/*projection_matrix*/,
666  const FullMatrix<double> &/*expansion_matrix*/)
667  {}
668  };
669 
670 
692  template <int dim>
694  {
695  public:
704  NamedAdditionalMaterialOutputs(const std::vector<std::string> &output_names);
705 
710 
715  const std::vector<std::string> &get_names() const;
716 
721  virtual std::vector<double> get_nth_output(const unsigned int idx) const = 0;
722 
723  virtual void average (const MaterialAveraging::AveragingOperation /*operation*/,
724  const FullMatrix<double> &/*projection_matrix*/,
725  const FullMatrix<double> &/*expansion_matrix*/)
726  {}
727 
728  private:
729  const std::vector<std::string> names;
730  };
731 
732 
738  template <int dim>
740  {
741  public:
742  SeismicAdditionalOutputs(const unsigned int n_points);
743 
744  virtual std::vector<double> get_nth_output(const unsigned int idx) const;
745 
751  std::vector<double> vs;
752 
758  std::vector<double> vp;
759  };
760 
761 
782  template <int dim>
784  {
785  public:
786  ReactionRateOutputs (const unsigned int n_points,
787  const unsigned int n_comp);
788 
789  virtual std::vector<double> get_nth_output(const unsigned int idx) const;
790 
798  std::vector<std::vector<double> > reaction_rates;
799  };
800 
801 
820  template <int dim>
822  {
823  public:
824  PrescribedFieldOutputs (const unsigned int n_points,
825  const unsigned int n_comp);
826 
827  virtual std::vector<double> get_nth_output(const unsigned int idx) const;
828 
836  std::vector<std::vector<double> > prescribed_field_outputs;
837  };
838 
839 
846  template <int dim>
848  {
849  public:
850  AdditionalMaterialOutputsStokesRHS(const unsigned int n_points)
851  : rhs_u(n_points), rhs_p(n_points), rhs_melt_pc(n_points)
852  {}
853 
855  {}
856 
857  virtual void average (const MaterialAveraging::AveragingOperation /*operation*/,
858  const FullMatrix<double> &/*projection_matrix*/,
859  const FullMatrix<double> &/*expansion_matrix*/)
860  {
861  // TODO: not implemented
862  }
863 
869  std::vector<Tensor<1,dim> > rhs_u;
870 
875  std::vector<double> rhs_p;
876 
881  std::vector<double> rhs_melt_pc;
882  };
883 
890  template <int dim>
892  {
893  public:
894  ElasticOutputs(const unsigned int n_points)
895  : elastic_force(n_points, numbers::signaling_nan<Tensor<2,dim> >() )
896  {}
897 
898  virtual ~ElasticOutputs()
899  {}
900 
901  virtual void average (const MaterialAveraging::AveragingOperation operation,
902  const FullMatrix<double> &/*projection_matrix*/,
903  const FullMatrix<double> &/*expansion_matrix*/)
904  {
905  AssertThrow(operation == MaterialAveraging::AveragingOperation::none,ExcNotImplemented());
906  return;
907  }
908 
914  std::vector<Tensor<2,dim> > elastic_force;
915  };
916 
917 
918 
934  template <int dim>
935  class Interface
936  {
937  public:
952 
957  virtual ~Interface();
958 
964  virtual
965  void
966  initialize ();
967 
972  virtual void update ();
973 
986  get_model_dependence () const;
987 
997  virtual bool is_compressible () const = 0;
1024  virtual double reference_viscosity () const = 0;
1034  virtual
1035  void evaluate (const MaterialModel::MaterialModelInputs<dim> &in,
1047  static
1048  void
1050 
1057  virtual
1058  void
1059  parse_parameters (ParameterHandler &prm);
1069  virtual
1070  void
1071  create_additional_named_outputs (MaterialModelOutputs &outputs) const;
1072 
1073 
1081  virtual
1082  void
1083  fill_additional_material_model_inputs(MaterialModel::MaterialModelInputs<dim> &input,
1084  const LinearAlgebra::BlockVector &solution,
1085  const FEValuesBase<dim> &fe_values,
1086  const Introspection<dim> &introspection) const;
1087 
1088  protected:
1105  };
1106 
1122  template <int dim>
1123  void
1124  register_material_model (const std::string &name,
1125  const std::string &description,
1126  void (*declare_parameters_function) (ParameterHandler &),
1127  Interface<dim> *(*factory_function) ());
1128 
1139  template <int dim>
1140  Interface<dim> *
1141  create_material_model (const std::string &model_name);
1142 
1143 
1155  template <int dim>
1156  Interface<dim> *
1158 
1159 
1166  template <int dim>
1167  std::string
1169 
1170 
1176  template <int dim>
1177  void
1179 
1180 
1181 
1191  template <int dim>
1192  void
1193  write_plugin_graph (std::ostream &output_stream);
1194 
1195 
1196 
1197 // --------------------- template function definitions ----------------------------------
1198 
1199  template <int dim>
1200  template <class AdditionalInputType>
1202  {
1203  for (unsigned int i=0; i<additional_inputs.size(); ++i)
1204  {
1205  AdditionalInputType *result = dynamic_cast<AdditionalInputType *> (additional_inputs[i].get());
1206  if (result)
1207  return result;
1208  }
1209  return nullptr;
1210  }
1211 
1212 
1213  template <int dim>
1214  template <class AdditionalInputType>
1215  const AdditionalInputType *MaterialModelInputs<dim>::get_additional_input() const
1216  {
1217  for (unsigned int i=0; i<additional_inputs.size(); ++i)
1218  {
1219  const AdditionalInputType *result = dynamic_cast<const AdditionalInputType *> (additional_inputs[i].get());
1220  if (result)
1221  return result;
1222  }
1223  return nullptr;
1224  }
1225 
1226 
1227  template <int dim>
1228  template <class AdditionalOutputType>
1230  {
1231  for (unsigned int i=0; i<additional_outputs.size(); ++i)
1232  {
1233  AdditionalOutputType *result = dynamic_cast<AdditionalOutputType *> (additional_outputs[i].get());
1234  if (result)
1235  return result;
1236  }
1237  return nullptr;
1238  }
1239 
1240 
1241  template <int dim>
1242  template <class AdditionalOutputType>
1243  const AdditionalOutputType *MaterialModelOutputs<dim>::get_additional_output() const
1244  {
1245  for (unsigned int i=0; i<additional_outputs.size(); ++i)
1246  {
1247  const AdditionalOutputType *result = dynamic_cast<const AdditionalOutputType *> (additional_outputs[i].get());
1248  if (result)
1249  return result;
1250  }
1251  return nullptr;
1252  }
1253 
1254 
1262 #define ASPECT_REGISTER_MATERIAL_MODEL(classname,name,description) \
1263  template class classname<2>; \
1264  template class classname<3>; \
1265  namespace ASPECT_REGISTER_MATERIAL_MODEL_ ## classname \
1266  { \
1267  aspect::internal::Plugins::RegisterHelper<aspect::MaterialModel::Interface<2>,classname<2> > \
1268  dummy_ ## classname ## _2d (&aspect::MaterialModel::register_material_model<2>, \
1269  name, description); \
1270  aspect::internal::Plugins::RegisterHelper<aspect::MaterialModel::Interface<3>,classname<3> > \
1271  dummy_ ## classname ## _3d (&aspect::MaterialModel::register_material_model<3>, \
1272  name, description); \
1273  }
1274  }
1275 }
1276 
1277 
1278 #endif
void write_plugin_graph(std::ostream &output_stream)
void average_property(const AveragingOperation operation, const FullMatrix< double > &projection_matrix, const FullMatrix< double > &expansion_matrix, std::vector< double > &values_out)
AdditionalInputType * get_additional_input()
Definition: interface.h:1201
std::vector< SymmetricTensor< 2, dim > > strain_rate
Definition: interface.h:293
std::vector< double > thermal_expansion_coefficients
Definition: interface.h:398
virtual void average(const MaterialAveraging::AveragingOperation, const FullMatrix< double > &, const FullMatrix< double > &)
Definition: interface.h:664
std::vector< std::vector< double > > prescribed_field_outputs
Definition: interface.h:836
MaterialModel::MaterialModelOutputs< dim > MaterialModelOutputs
Definition: interface.h:951
std::vector< std::vector< double > > composition
Definition: interface.h:279
std::vector< std::shared_ptr< AdditionalMaterialInputs< dim > > > additional_inputs
Definition: interface.h:331
std::vector< Point< dim > > position
Definition: interface.h:247
#define AssertThrow(cond, exc)
ElasticOutputs(const unsigned int n_points)
Definition: interface.h:894
AdditionalOutputType * get_additional_output()
Definition: interface.h:1229
std::vector< std::vector< double > > reaction_rates
Definition: interface.h:798
void declare_parameters(ParameterHandler &prm)
bool identifies_single_variable(const Dependence dependence)
void average(const AveragingOperation operation, const typename DoFHandler< dim >::active_cell_iterator &cell, const Quadrature< dim > &quadrature_formula, const Mapping< dim > &mapping, MaterialModelOutputs< dim > &values_out)
const std::vector< std::string > names
Definition: interface.h:729
std::vector< std::vector< double > > reaction_terms
Definition: interface.h:467
MaterialModel::MaterialModelInputs< dim > MaterialModelInputs
Definition: interface.h:944
AveragingOperation parse_averaging_operation_name(const std::string &s)
DoFHandler< dim >::active_cell_iterator current_cell
Definition: interface.h:324
virtual void average(const MaterialAveraging::AveragingOperation, const FullMatrix< double > &, const FullMatrix< double > &)
Definition: interface.h:723
AdditionalMaterialOutputsStokesRHS(const unsigned int n_points)
Definition: interface.h:850
Dependence operator|=(Dependence &d1, const Dependence d2)
Definition: interface.h:107
std::vector< Tensor< 1, dim > > pressure_gradient
Definition: interface.h:263
Interface< dim > * create_material_model(const std::string &model_name)
std::vector< Tensor< 2, dim > > elastic_force
Definition: interface.h:914
Definition: compat.h:38
std::vector< double > entropy_derivative_pressure
Definition: interface.h:421
std::vector< double > thermal_conductivities
Definition: interface.h:408
virtual void average(const MaterialAveraging::AveragingOperation, const FullMatrix< double > &, const FullMatrix< double > &)
Definition: interface.h:857
std::vector< double > entropy_derivative_temperature
Definition: interface.h:429
T signaling_nan()
NonlinearDependence::ModelDependence model_dependence
Definition: interface.h:1104
std::string get_valid_model_names_pattern()
static::ExceptionBase & ExcNotImplemented()
std::vector< Tensor< 1, dim > > velocity
Definition: interface.h:272
Dependence operator|(const Dependence d1, const Dependence d2)
Definition: interface.h:101
void register_material_model(const std::string &name, const std::string &description, void(*declare_parameters_function)(ParameterHandler &), Interface< dim > *(*factory_function)())
const DoFHandler< dim >::active_cell_iterator *cell DEAL_II_DEPRECATED
Definition: interface.h:307
std::vector< std::shared_ptr< AdditionalMaterialOutputs< dim > > > additional_outputs
Definition: interface.h:474
virtual void average(const MaterialAveraging::AveragingOperation operation, const FullMatrix< double > &, const FullMatrix< double > &)
Definition: interface.h:901