ASPECT
solution_evaluator.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2024 - 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_solution_evaluator_h
22 #define _aspect_solution_evaluator_h
23 
24 #include <aspect/global.h>
25 
26 
28 
29 #include <deal.II/base/array_view.h>
30 #include <deal.II/matrix_free/fe_point_evaluation.h>
31 #include <deal.II/non_matching/mapping_info.h>
32 #include <deal.II/dofs/dof_handler.h>
33 
34 namespace aspect
35 {
36  using namespace dealii;
37 
38  namespace internal
39  {
44  template <int dim>
46  {
47  public:
51  DynamicFEPointEvaluation(const unsigned int first_component, const unsigned int n_components)
52  : first_component (first_component),
53  n_components (n_components)
54  {}
55 
59  virtual ~DynamicFEPointEvaluation() = default;
60 
67  virtual void evaluate(const ArrayView<double> &solution_values,
68  const EvaluationFlags::EvaluationFlags flags) = 0;
69 
73  virtual
75  get_gradient(const unsigned int evaluation_point) const = 0;
76 
81  virtual
82  void
83  get_gradient(const unsigned int evaluation_point,
84  const ArrayView<Tensor<1,dim>> &gradients) const = 0;
85 
89  virtual
91  get_value(const unsigned int evaluation_point) const = 0;
92 
97  virtual
98  void
99  get_value(const unsigned int evaluation_point,
100  const ArrayView<double> &solution) const = 0;
101 
105  virtual
106  unsigned int
107  get_first_component() const final
108  {
109  return first_component;
110  }
111 
115  virtual
116  unsigned int
117  get_n_components() const final
118  {
119  return n_components;
120  }
121 
122  private:
126  unsigned int first_component;
127 
131  unsigned int n_components;
132  };
133  }
134 
140  template <int dim>
142  {
143  public:
149  SolutionEvaluator(const SimulatorAccess<dim> &simulator,
150  const UpdateFlags update_flags);
151 
156  void
157  reinit(const typename DoFHandler<dim>::active_cell_iterator &cell,
158  const ArrayView<Point<dim>> &positions);
159 
169  void
170  evaluate(const ArrayView<double> &solution_values,
171  const std::vector<EvaluationFlags::EvaluationFlags> &evaluation_flags);
172 
185  void get_solution(const unsigned int evaluation_point,
186  const ArrayView<double> &solution,
187  const std::vector<EvaluationFlags::EvaluationFlags> &evaluation_flags) const;
188 
201  void get_gradients(const unsigned int evaluation_point,
202  const ArrayView<Tensor<1, dim>> &gradients,
203  const std::vector<EvaluationFlags::EvaluationFlags> &evaluation_flags) const;
204 
209  FEPointEvaluation<dim, dim> &
210  get_velocity_or_fluid_velocity_evaluator(const bool use_fluid_velocity);
211 
215  NonMatching::MappingInfo<dim> &
216  get_mapping_info();
217 
221  unsigned int
222  n_components() const;
223 
224  private:
228  NonMatching::MappingInfo<dim> mapping_info;
229 
235  FEPointEvaluation<dim, dim> velocity;
236  std::unique_ptr<FEPointEvaluation<1, dim>> pressure;
237  FEPointEvaluation<1, dim> temperature;
238 
243  std::vector<std::unique_ptr<internal::DynamicFEPointEvaluation<dim>>> compositions;
244 
251  std::unique_ptr<FEPointEvaluation<dim, dim>> fluid_velocity;
252  std::unique_ptr<FEPointEvaluation<1, dim>> compaction_pressure;
253  std::unique_ptr<FEPointEvaluation<1, dim>> fluid_pressure;
254 
261  std::array<unsigned int, 3> melt_component_indices;
262 
268  };
269 
273  template <int dim>
274  std::unique_ptr<SolutionEvaluator<dim>>
275  construct_solution_evaluator(const SimulatorAccess<dim> &simulator_access,
276  const UpdateFlags update_flags);
277 }
278 
279 #endif
DynamicFEPointEvaluation(const unsigned int first_component, const unsigned int n_components)
virtual unsigned int get_n_components() const final
std::unique_ptr< FEPointEvaluation< 1, dim > > compaction_pressure
FEPointEvaluation< dim, dim > velocity
FEPointEvaluation< 1, dim > temperature
std::array< unsigned int, 3 > melt_component_indices
std::unique_ptr< FEPointEvaluation< dim, dim > > fluid_velocity
Definition: compat.h:59
NonMatching::MappingInfo< dim > mapping_info
std::unique_ptr< FEPointEvaluation< 1, dim > > pressure
virtual unsigned int get_first_component() const final
boost::container::small_vector< T, N > small_vector
Definition: global.h:245
Definition: compat.h:42
std::unique_ptr< FEPointEvaluation< 1, dim > > fluid_pressure
const SimulatorAccess< dim > & simulator_access
std::vector< std::unique_ptr< internal::DynamicFEPointEvaluation< dim > > > compositions
std::unique_ptr< SolutionEvaluator< dim > > construct_solution_evaluator(const SimulatorAccess< dim > &simulator_access, const UpdateFlags update_flags)