ASPECT
introspection.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 
22 #ifndef _aspect_introspection_h
23 #define _aspect_introspection_h
24 
25 #include <deal.II/base/index_set.h>
26 #include <deal.II/fe/component_mask.h>
27 #include <deal.II/fe/fe_values_extractors.h>
28 #include <deal.II/fe/fe.h>
29 
31 #include <aspect/parameters.h>
32 
33 
34 namespace aspect
35 {
36  using namespace dealii;
37 
42  template <int dim>
43  std::vector<VariableDeclaration<dim> >
44  construct_default_variables (const Parameters<dim> &parameters);
45 
46 
47 
65  template <int dim>
66  struct Introspection: public FEVariableCollection<dim>
67  {
68  public:
72  Introspection (const std::vector<VariableDeclaration<dim> > &variables,
73  const Parameters<dim> &parameters);
74 
78  ~Introspection ();
79 
80 
81 
94  const unsigned int n_components;
95 
99  const unsigned int n_compositional_fields;
100 
106 
112 
118  {
119  unsigned int velocities[dim];
120  unsigned int pressure;
121  unsigned int temperature;
122  std::vector<unsigned int> compositional_fields;
123  };
129 
135  const unsigned int n_blocks;
136 
142  {
143  unsigned int velocities;
144  unsigned int pressure;
145  unsigned int temperature;
146  std::vector<unsigned int> compositional_fields;
147  };
153 
158  struct Extractors
159  {
160  Extractors (const ComponentIndices &component_indices);
161 
165  const std::vector<FEValuesExtractors::Scalar> compositional_fields;
166  };
172 
184  {
185  unsigned int velocities;
186  unsigned int pressure;
187  unsigned int temperature;
188  unsigned int compositional_fields;
189  };
195 
204  {
205  unsigned int velocities;
206  unsigned int temperature;
207  unsigned int compositional_fields;
208  };
214 
221  {
223 
227  std::vector<ComponentMask> compositional_fields;
228  };
235 
249  std::vector<types::global_dof_index> system_dofs_per_block;
250 
256  struct IndexSets
257  {
265 
273  std::vector<IndexSet> system_partitioning;
274 
283  std::vector<IndexSet> system_relevant_partitioning;
284 
290  std::vector<IndexSet> stokes_partitioning;
291 
299 
305 
311  };
318 
324  std::vector<typename Parameters<dim>::AdvectionFieldMethod::Kind> compositional_field_methods;
325 
338  unsigned int
339  compositional_index_for_name (const std::string &name) const;
340 
347  std::string
348  name_for_compositional_index (const unsigned int index) const;
349 
353  const std::vector<std::string> &
354  get_composition_names () const;
355 
364  bool
365  compositional_name_exists (const std::string &name) const;
366 
374  bool
375  is_stokes_component (const unsigned int component_index) const;
376 
377  private:
382  std::vector<std::string> composition_names;
383 
384  };
385 }
386 
387 
388 #endif
std::vector< IndexSet > system_relevant_partitioning
const FEValuesExtractors::Vector velocities
std::vector< VariableDeclaration< dim > > construct_default_variables(const Parameters< dim > &parameters)
std::vector< unsigned int > compositional_fields
std::vector< unsigned int > compositional_fields
const unsigned int n_blocks
const bool use_discontinuous_composition_discretization
const BaseElements base_elements
const PolynomialDegree polynomial_degree
const FEValuesExtractors::Scalar temperature
const unsigned int n_components
Definition: introspection.h:94
const Extractors extractors
std::vector< IndexSet > system_partitioning
std::vector< typename Parameters< dim >::AdvectionFieldMethod::Kind > compositional_field_methods
const unsigned int n_compositional_fields
Definition: introspection.h:99
const ComponentIndices component_indices
Definition: compat.h:53
const bool use_discontinuous_temperature_discretization
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< std::string > composition_names
std::vector< ComponentMask > compositional_fields
std::vector< IndexSet > stokes_partitioning