ASPECT
simulator_access.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2020 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 
30 #include <deal.II/base/timer.h>
34 #include <deal.II/fe/fe.h>
35 #include <deal.II/fe/mapping_q.h>
36 
37 #if !DEAL_II_VERSION_GTE(9,1,0)
39 #else
41 #endif
42 
43 namespace WorldBuilder
44 {
45  class World;
46 }
47 
48 namespace aspect
49 {
50  using namespace dealii;
51 
52  // forward declarations:
53  template <int dim> class Simulator;
54  template <int dim> struct SimulatorSignals;
55  template <int dim> class LateralAveraging;
56  template <int dim> struct RotationProperties;
57 
58  namespace GravityModel
59  {
60  template <int dim> class Interface;
61  }
62 
63  namespace HeatingModel
64  {
65  template <int dim> class Manager;
66  }
67 
68  namespace MaterialModel
69  {
70  template <int dim> class Interface;
71  }
72 
73  namespace InitialTemperature
74  {
75  template <int dim> class Manager;
76  template <int dim> class Interface;
77  }
78 
79  namespace BoundaryTemperature
80  {
81  template <int dim> class Manager;
82  template <int dim> class Interface;
83  }
84 
85  namespace BoundaryHeatFlux
86  {
87  template <int dim> class Interface;
88  }
89 
90  namespace BoundaryComposition
91  {
92  template <int dim> class Manager;
93  template <int dim> class Interface;
94  }
95 
96  namespace BoundaryTraction
97  {
98  template <int dim> class Interface;
99  }
100 
101  namespace BoundaryVelocity
102  {
103  template <int dim> class Manager;
104  template <int dim> class Interface;
105  }
106 
107  namespace InitialComposition
108  {
109  template <int dim> class Manager;
110  template <int dim> class Interface;
111  }
112 
113  namespace InitialTopographyModel
114  {
115  template <int dim> class Interface;
116  }
117 
118  namespace MeshRefinement
119  {
120  template <int dim> class Manager;
121  }
122 
123  namespace AdiabaticConditions
124  {
125  template <int dim> class Interface;
126  }
127 
128  namespace Postprocess
129  {
130  template <int dim> class Manager;
131  }
132 
133  template <int dim> class MeltHandler;
134  template <int dim> class VolumeOfFluidHandler;
135 
136  namespace MeshDeformation
137  {
138  template <int dim> class MeshDeformationHandler;
139  }
140 
141  template <int dim> class NewtonHandler;
142 
143  template <int dim> class StokesMatrixFreeHandler;
144 
145  namespace Particle
146  {
147  template <int dim> class World;
148  }
149 
171  template <int dim>
172  class SimulatorAccess
173  {
174  public:
181  SimulatorAccess ();
182 
187  SimulatorAccess (const Simulator<dim> &simulator_object);
188 
193  virtual
194  ~SimulatorAccess ();
195 
205  virtual void initialize_simulator (const Simulator<dim> &simulator_object);
206 
217  const Introspection<dim> &
218  introspection () const;
219 
227  const Simulator<dim> &
228  get_simulator () const;
229 
234  const Parameters<dim> &
235  get_parameters () const;
236 
241  get_signals() const;
242 
246  MPI_Comm
247  get_mpi_communicator () const;
248 
255  TimerOutput &
256  get_computing_timer () const;
257 
263  const ConditionalOStream &
264  get_pcout () const;
265 
269  double get_time () const;
270 
274  double
275  get_timestep () const;
276 
280  double
281  get_old_timestep () const;
282 
286  unsigned int
287  get_timestep_number () const;
288 
292  unsigned int
293  get_nonlinear_iteration () const;
294 
300  get_triangulation () const;
301 
305  double
306  get_volume () const;
307 
312  const Mapping<dim> &
313  get_mapping () const;
314 
320  std::string
321  get_output_directory () const;
322 
326  bool
327  include_adiabatic_heating () const;
328 
332  bool
333  include_latent_heat () const;
334 
338  bool
339  include_melt_transport () const;
340 
344  int
345  get_stokes_velocity_degree () const;
346 
350  double
351  get_adiabatic_surface_temperature () const;
352 
356  double
357  get_surface_pressure () const;
358 
364  bool
365  convert_output_to_years () const;
366 
376  unsigned int
377  get_pre_refinement_step () const;
378 
383  unsigned int
384  n_compositional_fields () const;
385 
392  void
393  get_refinement_criteria(Vector<float> &estimated_error_per_cell) const;
394 
405  void
406  get_artificial_viscosity(Vector<float> &viscosity_per_cell,
407  const bool skip_interior_cells = false) const;
408 
413  void
414  get_artificial_viscosity_composition(Vector<float> &viscosity_per_cell,
415  const unsigned int compositional_variable) const;
434  get_current_linearization_point () const;
435 
447  get_solution () const;
448 
458  get_old_solution () const;
459 
469  get_old_old_solution () const;
470 
479  get_reaction_vector () const;
480 
489  get_mesh_velocity () const;
490 
495  const DoFHandler<dim> &
496  get_dof_handler () const;
497 
506  const FiniteElement<dim> &
507  get_fe () const;
508 
513  get_system_matrix () const;
514 
519  get_system_preconditioner_matrix () const;
520 
532  get_material_model () const;
533 
538  void
539  compute_material_model_input_values (const LinearAlgebra::BlockVector &input_solution,
540  const FEValuesBase<dim,dim> &input_finite_element_values,
541  const typename DoFHandler<dim>::active_cell_iterator &cell,
542  const bool compute_strainrate,
543  MaterialModel::MaterialModelInputs<dim> &material_model_inputs) const;
544 
549  get_gravity_model () const;
550 
555  get_initial_topography_model () const;
556 
561  get_geometry_model () const;
562 
563 
569  get_adiabatic_conditions () const;
570 
579  bool has_boundary_temperature () const;
580 
587  DEAL_II_DEPRECATED
589  get_boundary_temperature () const;
590 
598  get_boundary_temperature_manager () const;
599 
605  get_boundary_heat_flux () const;
606 
614  bool has_boundary_composition () const;
615 
622  DEAL_II_DEPRECATED
624  get_boundary_composition () const;
625 
633  get_boundary_composition_manager () const;
634 
639  const std::map<types::boundary_id,std::unique_ptr<BoundaryTraction::Interface<dim> > > &
640  get_boundary_traction () const;
641 
648  DEAL_II_DEPRECATED
650  get_initial_temperature () const;
651 
659  get_initial_temperature_manager () const;
660 
665  DEAL_II_DEPRECATED
667  get_initial_composition () const;
668 
675  get_initial_composition_manager () const;
676 
681  const std::set<types::boundary_id> &
682  get_fixed_temperature_boundary_indicators () const;
683 
688  const std::set<types::boundary_id> &
689  get_fixed_heat_flux_boundary_indicators () const;
690 
695  const std::set<types::boundary_id> &
696  get_fixed_composition_boundary_indicators () const;
697 
704  const std::set<types::boundary_id> &
705  get_mesh_deformation_boundary_indicators () const;
706 
714  get_boundary_velocity_manager () const;
715 
722  get_heating_model_manager () const;
723 
731  get_mesh_refinement_manager () const;
732 
736  const MeltHandler<dim> &
737  get_melt_handler () const;
738 
743  get_volume_of_fluid_handler () const;
744 
749  const NewtonHandler<dim> &
750  get_newton_handler () const;
751 #ifdef ASPECT_WITH_WORLD_BUILDER
752 
759  const WorldBuilder::World &
760  get_world_builder () const;
761 #endif
762 
767  get_mesh_deformation_handler () const;
768 
774  const LateralAveraging<dim> &
775  get_lateral_averaging () const;
776 
782  get_current_constraints () const;
783 
816  bool simulator_is_past_initialization () const;
817 
822  double
823  get_pressure_scaling () const;
824 
831  bool
832  pressure_rhs_needs_compatibility_modification() const;
833 
837  bool
838  model_has_prescribed_stokes_solution () const;
839 
845  static
846  void
847  get_composition_values_at_q_point (const std::vector<std::vector<double> > &composition_values,
848  const unsigned int q,
849  std::vector<double> &composition_values_at_q_point);
850 
862  TableHandler &get_statistics_object() const;
863 
864 
878  template <typename PostprocessorType>
879  DEAL_II_DEPRECATED
880  PostprocessorType *
881  find_postprocessor () const;
882 
887  get_postprocess_manager () const;
888 
893  const Particle::World<dim> &
894  get_particle_world() const;
895 
903  get_particle_world();
904 
908  bool is_stokes_matrix_free();
909 
915  get_stokes_matrix_free () const;
916 
929  compute_net_angular_momentum(const bool use_constant_density,
930  const LinearAlgebra::BlockVector &solution,
931  const bool limit_to_top_faces = false) const;
932 
935  private:
940  };
941 
942  template <int dim>
943  template <typename PostprocessorType>
944  inline
945  PostprocessorType *
947  {
948  if (get_postprocess_manager().template has_matching_postprocessor<PostprocessorType>())
949  return &get_postprocess_manager().template get_matching_postprocessor<PostprocessorType>();
950 
951  return nullptr;
952  }
953 }
954 
955 
956 #endif
const Simulator< dim > * simulator
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const Triangulation< dim, spacedim > & get_triangulation() const
std::string get_time()