ASPECT
introspection.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2022 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 
22 #ifndef _aspect_introspection_h
23 #define _aspect_introspection_h
24 
25 #include <deal.II/base/index_set.h>
26 #include <deal.II/base/quadrature.h>
27 #include <deal.II/fe/component_mask.h>
28 #include <deal.II/fe/fe_values_extractors.h>
29 #include <deal.II/fe/fe.h>
30 
32 #include <aspect/parameters.h>
33 
34 
35 namespace aspect
36 {
37  using namespace dealii;
38 
43  template <int dim>
44  std::vector<VariableDeclaration<dim>>
45  construct_default_variables (const Parameters<dim> &parameters);
46 
47 
48 
66  template <int dim>
67  struct Introspection: public FEVariableCollection<dim>
68  {
69  public:
73  Introspection (const std::vector<VariableDeclaration<dim>> &variables,
74  const Parameters<dim> &parameters);
75 
79  ~Introspection ();
80 
81 
82 
95  const unsigned int n_components;
96 
100  const unsigned int n_compositional_fields;
101 
107 
113 
119  {
120  unsigned int velocities[dim];
121  unsigned int pressure;
122  unsigned int temperature;
123  std::vector<unsigned int> compositional_fields;
124  };
130 
136  const unsigned int n_blocks;
137 
143  {
144  unsigned int velocities;
145  unsigned int pressure;
146  unsigned int temperature;
147  std::vector<unsigned int> compositional_fields;
148  };
154 
159  struct Extractors
160  {
161  Extractors (const ComponentIndices &component_indices);
162 
164  const FEValuesExtractors::Scalar pressure;
165  const FEValuesExtractors::Scalar temperature;
166  const std::vector<FEValuesExtractors::Scalar> compositional_fields;
167  };
173 
185  {
186  unsigned int velocities;
187  unsigned int pressure;
188  unsigned int temperature;
189  unsigned int compositional_fields;
190  };
196 
205  {
206  unsigned int velocities;
207  unsigned int temperature;
208  unsigned int compositional_fields;
209  };
215 
235  struct Quadratures
236  {
237  Quadrature<dim> velocities;
238  Quadrature<dim> pressure;
239  Quadrature<dim> temperature;
240  Quadrature<dim> compositional_fields;
241  Quadrature<dim> system;
242  };
248 
258  {
259  Quadrature<dim-1> velocities;
260  Quadrature<dim-1> pressure;
261  Quadrature<dim-1> temperature;
262  Quadrature<dim-1> compositional_fields;
263  Quadrature<dim-1> system;
264  };
270 
277  {
279 
280  ComponentMask velocities;
281  ComponentMask pressure;
282  ComponentMask temperature;
283  std::vector<ComponentMask> compositional_fields;
284  };
291 
305  std::vector<types::global_dof_index> system_dofs_per_block;
306 
312  struct IndexSets
313  {
321 
329  std::vector<IndexSet> system_partitioning;
330 
339  std::vector<IndexSet> system_relevant_partitioning;
340 
346  std::vector<IndexSet> stokes_partitioning;
347 
355 
361 
367  };
374 
380 
386  std::vector<typename Parameters<dim>::AdvectionFieldMethod::Kind> compositional_field_methods;
387 
400  unsigned int
401  compositional_index_for_name (const std::string &name) const;
402 
409  std::string
410  name_for_compositional_index (const unsigned int index) const;
411 
415  const std::vector<std::string> &
416  get_composition_names () const;
417 
422  const std::vector<typename Parameters<dim>::CompositionalFieldDescription> &
423  get_composition_descriptions () const;
424 
425 
434  bool
435  composition_type_exists (const typename Parameters<dim>::CompositionalFieldDescription::Type &type) const;
436 
437 
447  unsigned int
448  find_composition_type (const typename Parameters<dim>::CompositionalFieldDescription::Type &type) const;
449 
450 
459  bool
460  compositional_name_exists (const std::string &name) const;
461 
469  bool
470  is_stokes_component (const unsigned int component_index) const;
471 
472  private:
477  std::vector<std::string> composition_names;
478 
484  std::vector<typename Parameters<dim>::CompositionalFieldDescription> composition_descriptions;
485 
486  };
487 }
488 
489 
490 #endif
std::vector< IndexSet > system_relevant_partitioning
const FEValuesExtractors::Vector velocities
std::vector< VariableDeclaration< dim > > construct_default_variables(const Parameters< dim > &parameters)
Parameters< dim >::AdvectionFieldMethod::Kind temperature_method
std::vector< unsigned int > compositional_fields
std::vector< unsigned int > compositional_fields
const unsigned int n_blocks
const bool use_discontinuous_composition_discretization
::TrilinosWrappers::MPI::Vector Vector
Definition: global.h:223
const BaseElements base_elements
const PolynomialDegree polynomial_degree
const FEValuesExtractors::Scalar temperature
const unsigned int n_components
Definition: introspection.h:95
const Extractors extractors
std::vector< IndexSet > system_partitioning
std::vector< typename Parameters< dim >::AdvectionFieldMethod::Kind > compositional_field_methods
Quadrature< dim > compositional_fields
const unsigned int n_compositional_fields
const ComponentIndices component_indices
const bool use_discontinuous_temperature_discretization
const FaceQuadratures face_quadratures
const Quadratures quadratures
const BlockIndices block_indices
const std::vector< FEValuesExtractors::Scalar > compositional_fields
const FEValuesExtractors::Scalar pressure
const ComponentMasks component_masks
std::vector< types::global_dof_index > system_dofs_per_block
std::vector< typename Parameters< dim >::CompositionalFieldDescription > composition_descriptions
Definition: compat.h:88
std::vector< std::string > composition_names
std::vector< ComponentMask > compositional_fields
std::vector< IndexSet > stokes_partitioning