ASPECT
simulator_access.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2021 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_simulator_access_h
23 #define _aspect_simulator_access_h
24 
25 #include <aspect/global.h>
26 #include <aspect/parameters.h>
27 #include <aspect/introspection.h>
28 
29 #include <deal.II/base/table_handler.h>
30 #include <deal.II/base/timer.h>
31 #include <deal.II/base/conditional_ostream.h>
32 #include <deal.II/distributed/tria.h>
33 #include <deal.II/dofs/dof_handler.h>
34 #include <deal.II/fe/fe.h>
35 #include <deal.II/fe/mapping_q.h>
36 #include <deal.II/lac/affine_constraints.h>
37 
38 namespace WorldBuilder
39 {
40  class World;
41 }
42 
43 namespace aspect
44 {
45  using namespace dealii;
46 
47  // forward declarations:
48  template <int dim> class Simulator;
49  template <int dim> struct SimulatorSignals;
50  template <int dim> class LateralAveraging;
51  template <int dim> struct RotationProperties;
52 
53  namespace GravityModel
54  {
55  template <int dim> class Interface;
56  }
57 
58  namespace HeatingModel
59  {
60  template <int dim> class Manager;
61  }
62 
63  namespace MaterialModel
64  {
65  template <int dim> class Interface;
66  }
67 
68  namespace InitialTemperature
69  {
70  template <int dim> class Manager;
71  template <int dim> class Interface;
72  }
73 
74  namespace BoundaryTemperature
75  {
76  template <int dim> class Manager;
77  template <int dim> class Interface;
78  }
79 
80  namespace BoundaryHeatFlux
81  {
82  template <int dim> class Interface;
83  }
84 
85  namespace BoundaryComposition
86  {
87  template <int dim> class Manager;
88  template <int dim> class Interface;
89  }
90 
91  namespace BoundaryTraction
92  {
93  template <int dim> class Interface;
94  }
95 
96  namespace BoundaryVelocity
97  {
98  template <int dim> class Manager;
99  template <int dim> class Interface;
100  }
101 
102  namespace InitialComposition
103  {
104  template <int dim> class Manager;
105  template <int dim> class Interface;
106  }
107 
108  namespace InitialTopographyModel
109  {
110  template <int dim> class Interface;
111  }
112 
113  namespace MeshRefinement
114  {
115  template <int dim> class Manager;
116  }
117 
118  namespace AdiabaticConditions
119  {
120  template <int dim> class Interface;
121  }
122 
123  namespace Postprocess
124  {
125  template <int dim> class Manager;
126  }
127 
128  template <int dim> class MeltHandler;
129  template <int dim> class VolumeOfFluidHandler;
130 
131  namespace MeshDeformation
132  {
133  template <int dim> class MeshDeformationHandler;
134  }
135 
136  template <int dim> class NewtonHandler;
137 
138  template <int dim> class StokesMatrixFreeHandler;
139 
140  namespace Particle
141  {
142  template <int dim> class World;
143  }
144 
166  template <int dim>
167  class SimulatorAccess
168  {
169  public:
176  SimulatorAccess ();
177 
182  SimulatorAccess (const Simulator<dim> &simulator_object);
183 
188  virtual ~SimulatorAccess () = default;
189 
199  virtual void initialize_simulator (const Simulator<dim> &simulator_object);
200 
211  const Introspection<dim> &
212  introspection () const;
213 
221  const Simulator<dim> &
222  get_simulator () const;
223 
228  const Parameters<dim> &
229  get_parameters () const;
230 
235  get_signals() const;
236 
240  MPI_Comm
241  get_mpi_communicator () const;
242 
249  TimerOutput &
250  get_computing_timer () const;
251 
257  const ConditionalOStream &
258  get_pcout () const;
259 
263  double get_time () const;
264 
268  double
269  get_timestep () const;
270 
274  double
275  get_old_timestep () const;
276 
280  unsigned int
281  get_timestep_number () const;
282 
286  unsigned int
287  get_nonlinear_iteration () const;
288 
293  const parallel::distributed::Triangulation<dim> &
294  get_triangulation () const;
295 
299  double
300  get_volume () const;
301 
306  const Mapping<dim> &
307  get_mapping () const;
308 
314  std::string
315  get_output_directory () const;
316 
320  bool
321  include_adiabatic_heating () const;
322 
326  bool
327  include_latent_heat () const;
328 
332  bool
333  include_melt_transport () const;
334 
338  int
339  get_stokes_velocity_degree () const;
340 
344  double
345  get_adiabatic_surface_temperature () const;
346 
350  double
351  get_surface_pressure () const;
352 
358  bool
359  convert_output_to_years () const;
360 
370  unsigned int
371  get_pre_refinement_step () const;
372 
377  unsigned int
378  n_compositional_fields () const;
379 
386  void
387  get_refinement_criteria(Vector<float> &estimated_error_per_cell) const;
388 
399  void
400  get_artificial_viscosity(Vector<float> &viscosity_per_cell,
401  const bool skip_interior_cells = false) const;
402 
407  void
408  get_artificial_viscosity_composition(Vector<float> &viscosity_per_cell,
409  const unsigned int compositional_variable) const;
428  get_current_linearization_point () const;
429 
441  get_solution () const;
442 
452  get_old_solution () const;
453 
463  get_old_old_solution () const;
464 
473  get_reaction_vector () const;
474 
483  get_mesh_velocity () const;
484 
489  const DoFHandler<dim> &
490  get_dof_handler () const;
491 
500  const FiniteElement<dim> &
501  get_fe () const;
502 
507  get_system_matrix () const;
508 
513  get_system_preconditioner_matrix () const;
514 
526  get_material_model () const;
527 
532  void
533  compute_material_model_input_values (const LinearAlgebra::BlockVector &input_solution,
534  const FEValuesBase<dim,dim> &input_finite_element_values,
535  const typename DoFHandler<dim>::active_cell_iterator &cell,
536  const bool compute_strainrate,
537  MaterialModel::MaterialModelInputs<dim> &material_model_inputs) const;
538 
543  get_gravity_model () const;
544 
549  get_initial_topography_model () const;
550 
555  get_geometry_model () const;
556 
557 
563  get_adiabatic_conditions () const;
564 
573  bool has_boundary_temperature () const;
574 
583  get_boundary_temperature () const;
584 
592  get_boundary_temperature_manager () const;
593 
599  get_boundary_heat_flux () const;
600 
608  bool has_boundary_composition () const;
609 
618  get_boundary_composition () const;
619 
627  get_boundary_composition_manager () const;
628 
633  const std::map<types::boundary_id,std::unique_ptr<BoundaryTraction::Interface<dim>>> &
634  get_boundary_traction () const;
635 
644  get_initial_temperature () const;
645 
653  get_initial_temperature_manager () const;
654 
661  get_initial_composition () const;
662 
669  get_initial_composition_manager () const;
670 
675  const std::set<types::boundary_id> &
676  get_fixed_temperature_boundary_indicators () const;
677 
682  const std::set<types::boundary_id> &
683  get_fixed_heat_flux_boundary_indicators () const;
684 
689  const std::set<types::boundary_id> &
690  get_fixed_composition_boundary_indicators () const;
691 
698  const std::set<types::boundary_id> &
699  get_mesh_deformation_boundary_indicators () const;
700 
708  get_boundary_velocity_manager () const;
709 
716  get_heating_model_manager () const;
717 
725  get_mesh_refinement_manager () const;
726 
730  const MeltHandler<dim> &
731  get_melt_handler () const;
732 
737  get_volume_of_fluid_handler () const;
738 
743  const NewtonHandler<dim> &
744  get_newton_handler () const;
745 #ifdef ASPECT_WITH_WORLD_BUILDER
746 
753  const WorldBuilder::World &
754  get_world_builder () const;
755 #endif
756 
761  get_mesh_deformation_handler () const;
762 
768  const LateralAveraging<dim> &
769  get_lateral_averaging () const;
770 
775  const AffineConstraints<double> &
776  get_current_constraints () const;
777 
810  bool simulator_is_past_initialization () const;
811 
816  double
817  get_pressure_scaling () const;
818 
825  bool
826  pressure_rhs_needs_compatibility_modification() const;
827 
831  bool
832  model_has_prescribed_stokes_solution () const;
833 
839  static
840  void
841  get_composition_values_at_q_point (const std::vector<std::vector<double>> &composition_values,
842  const unsigned int q,
843  std::vector<double> &composition_values_at_q_point);
844 
856  TableHandler &get_statistics_object() const;
857 
858 
872  template <typename PostprocessorType>
874  PostprocessorType *
875  find_postprocessor () const;
876 
881  get_postprocess_manager () const;
882 
887  const Particle::World<dim> &
888  get_particle_world() const;
889 
897  get_particle_world();
898 
902  bool is_stokes_matrix_free();
903 
909  get_stokes_matrix_free () const;
910 
923  compute_net_angular_momentum(const bool use_constant_density,
924  const LinearAlgebra::BlockVector &solution,
925  const bool limit_to_top_faces = false) const;
926 
929  private:
934  };
935 
936  template <int dim>
937  template <typename PostprocessorType>
938  inline
939  PostprocessorType *
941  {
942  if (get_postprocess_manager().template has_matching_postprocessor<PostprocessorType>())
943  return &get_postprocess_manager().template get_matching_postprocessor<PostprocessorType>();
944 
945  return nullptr;
946  }
947 }
948 
949 
950 #endif
const Simulator< dim > * simulator
::TrilinosWrappers::MPI::BlockVector BlockVector
Definition: global.h:229
StructuredDataLookup< dim > DEAL_II_DEPRECATED
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
Definition: global.h:240