ASPECT
lateral_averaging.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_lateral_averaging_h
23 #define _aspect_lateral_averaging_h
24 
26 
27 #include <deal.II/fe/fe_values.h>
28 
29 namespace aspect
30 {
31  using namespace dealii;
32 
33  namespace internal
34  {
43  template <int dim>
45  {
46  public:
51  virtual
52  bool
53  need_material_properties() const;
54 
59  virtual
60  void
61  create_additional_material_model_outputs (const unsigned int n_points,
63 
68  virtual
69  void
70  setup(const unsigned int q_points);
71 
79  virtual
80  void
81  operator()(const MaterialModel::MaterialModelInputs<dim> &in,
83  const FEValues<dim> &fe_values,
84  const LinearAlgebra::BlockVector &solution,
85  std::vector<double> &output) = 0;
86 
90  virtual
91  ~FunctorBase();
92  };
93  }
94 
109  template <int dim>
110  class LateralAveraging : public SimulatorAccess<dim>
111  {
112  public:
127  std::vector<std::vector<double> >
128  get_averages(const unsigned int n_slices,
129  const std::vector<std::string> &property_names) const;
130 
140  void
141  get_temperature_averages(std::vector<double> &values) const;
142 
155  void
156  get_composition_averages(const unsigned int composition_index,
157  std::vector<double> &values) const;
158 
166  void
167  get_viscosity_averages(std::vector<double> &values) const;
168 
176  void
177  get_velocity_magnitude_averages(std::vector<double> &values) const;
178 
186  void
187  get_sinking_velocity_averages(std::vector<double> &values) const;
188 
196  void
197  get_Vs_averages(std::vector<double> &values) const;
198 
206  void
207  get_Vp_averages(std::vector<double> &values) const;
208 
217  void
218  get_vertical_heat_flux_averages(std::vector<double> &values) const;
219 
231  void
232  get_vertical_mass_flux_averages(std::vector<double> &values) const;
233 
234  private:
250  std::vector<std::vector<double> >
251  compute_lateral_averages(const unsigned int n_slices,
252  std::vector<std::unique_ptr<internal::FunctorBase<dim> > > &functors) const;
253  };
254 }
255 
256 
257 #endif
Definition: compat.h:53