ASPECT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 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 
22 #ifndef _aspect_mesh_deformation_interface_h
23 #define _aspect_mesh_deformation_interface_h
24 
25 #include <aspect/plugins.h>
27 #include <aspect/global.h>
28 
29 #include <deal.II/fe/fe_system.h>
30 #include <deal.II/dofs/dof_handler.h>
31 #include <deal.II/lac/affine_constraints.h>
32 #include <deal.II/base/index_set.h>
33 #include <deal.II/base/mg_level_object.h>
34 #include <deal.II/lac/la_parallel_vector.h>
35 #include <deal.II/multigrid/mg_constrained_dofs.h>
36 #include <deal.II/multigrid/mg_transfer_matrix_free.h>
37 #include <deal.II/multigrid/mg_transfer_global_coarsening.templates.h>
39 
40 namespace aspect
41 {
42  namespace Assemblers
43  {
50  template <int dim>
52  public SimulatorAccess<dim>
53  {
54  public:
55  ApplyStabilization(const double stabilization_theta);
56 
57  void
60 
61  private:
67  const double free_surface_theta;
68  };
69  }
70 
71  template <int dim> class Simulator;
72 
77  namespace MeshDeformation
78  {
88  template <int dim>
90  {
91  public:
95  virtual bool needs_surface_stabilization() const;
96 
97 
105  virtual
106  Tensor<1,dim>
107  compute_initial_deformation_on_boundary(const types::boundary_id boundary_indicator,
108  const Point<dim> &position) const;
109 
117  virtual
118  void
119  compute_velocity_constraints_on_boundary(const DoFHandler<dim> &mesh_deformation_dof_handler,
120  AffineConstraints<double> &mesh_velocity_constraints,
121  const std::set<types::boundary_id> &boundary_ids) const;
122 
123 
145  virtual
146  void save (std::map<std::string, std::string> &status_strings) const;
147 
159  virtual
160  void load (const std::map<std::string, std::string> &status_strings);
161  };
162 
163 
164 
170  template <int dim>
171  class MeshDeformationHandler: public SimulatorAccess<dim>
172  {
173  public:
181 
185  ~MeshDeformationHandler() override;
186 
192  void initialize();
193 
198  void set_assemblers(const SimulatorAccess<dim> &simulator_access,
199  aspect::Assemblers::Manager<dim> &assemblers) const;
200 
205  void update();
206 
214  void execute();
215 
220  void setup_dofs();
221 
225  static
226  void declare_parameters (ParameterHandler &prm);
227 
231  void parse_parameters (ParameterHandler &prm);
232 
237  template <class Archive>
238  void save (Archive &ar,
239  const unsigned int version) const;
240 
245  template <class Archive>
246  void load (Archive &ar,
247  const unsigned int version);
248 
249  BOOST_SERIALIZATION_SPLIT_MEMBER()
250 
251 
268  static
269  void
270  register_mesh_deformation
271  (const std::string &name,
272  const std::string &description,
273  void (*declare_parameters_function) (ParameterHandler &),
274  std::unique_ptr<Interface<dim>> (*factory_function) ());
275 
280  const std::map<types::boundary_id, std::vector<std::string>> &
281  get_active_mesh_deformation_names () const;
282 
287  const std::map<types::boundary_id,std::vector<std::unique_ptr<Interface<dim>>>> &
288  get_active_mesh_deformation_models () const;
289 
294  const std::set<types::boundary_id> &
295  get_active_mesh_deformation_boundary_indicators () const;
296 
301  const std::set<types::boundary_id> &
302  get_boundary_indicators_requiring_stabilization () const;
303 
309  const std::set<types::boundary_id> &
310  get_free_surface_boundary_indicators () const;
311 
315  double get_free_surface_theta () const;
316 
331  const LinearAlgebra::Vector &
332  get_initial_topography () const;
333 
338  const LinearAlgebra::Vector &
339  get_mesh_displacements () const;
340 
344  const DoFHandler<dim> &
345  get_mesh_deformation_dof_handler () const;
346 
356  template <typename MeshDeformationType,
357  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshDeformationType>::value>>
358  bool
359  has_matching_mesh_deformation_object () const;
360 
372  template <typename MeshDeformationType,
373  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshDeformationType>::value>>
374  const MeshDeformationType &
375  get_matching_mesh_deformation_object () const;
376 
377 
383  const Mapping<dim> &
384  get_level_mapping(const unsigned int level) const;
385 
395  static
396  void
397  write_plugin_graph (std::ostream &output_stream);
398 
402  DeclException1 (ExcMeshDeformationNameNotFound,
403  std::string,
404  << "Could not find entry <"
405  << arg1
406  << "> among the names of registered mesh deformation objects.");
407 
408  private:
416  void make_initial_constraints ();
417 
430  void make_constraints ();
431 
436  void compute_mesh_displacements ();
437 
442  void compute_mesh_displacements_gmg ();
443 
459  void set_initial_topography ();
460 
464  void interpolate_mesh_velocity ();
465 
469  void update_multilevel_deformation ();
470 
476 
481  const FESystem<dim> mesh_deformation_fe;
482 
487 
494 
502 
507 
517 
527 
532 
537 
542  AffineConstraints<double> mesh_velocity_constraints;
543 
548  AffineConstraints<double> mesh_vertex_constraints;
549 
554  std::map<types::boundary_id,std::vector<std::unique_ptr<Interface<dim>>>> mesh_deformation_objects;
555 
560  std::map<types::boundary_id, std::vector<std::string>> mesh_deformation_object_names;
561 
568 
574 
581  std::set<types::boundary_id> zero_mesh_deformation_boundary_indicators;
582 
588  std::set<types::boundary_id> free_surface_boundary_indicators;
589 
594  std::set<types::boundary_id> boundary_indicators_requiring_stabilization;
595 
597 
604 
608  MGLevelObject<std::unique_ptr<Mapping<dim>>> level_mappings;
609 
613  MGLevelObject<::LinearAlgebra::distributed::Vector<double>> level_displacements;
614 
619 
623  MGConstrainedDoFs mg_constrained_dofs;
624 
625  friend class Simulator<dim>;
626  friend class SimulatorAccess<dim>;
627  };
628 
629 
630  template <int dim>
631  template <class Archive>
633  const unsigned int) const
634  {
635  // let all the mesh deformation plugins save their data in a map and then
636  // serialize that
637  //TODO: for now we assume the same plugins are active before and after restart.
638  std::map<std::string,std::string> saved_text;
639  for (const auto &boundary_id_and_mesh_deformation_objects : mesh_deformation_objects)
640  for (const auto &p : boundary_id_and_mesh_deformation_objects.second)
641  p->save (saved_text);
642 
643  ar &saved_text;
644  }
645 
646 
647  template <int dim>
648  template <class Archive>
650  const unsigned int)
651  {
652  // get the map back out of the stream; then let the mesh deformation plugins
653  // that we currently have get their data from there. note that this
654  // may not be the same set ofmesh deformation plugins we had when we saved
655  // their data
656  std::map<std::string,std::string> saved_text;
657  ar &saved_text;
658 
659  for (const auto &boundary_id_and_mesh_deformation_objects : mesh_deformation_objects)
660  for (const auto &p : boundary_id_and_mesh_deformation_objects.second)
661  p->load (saved_text);
662  }
663 
664 
665 
666 
667  template <int dim>
668  template <typename MeshDeformationType, typename>
669  inline
670  bool
672  {
673  for (const auto &object_iterator : mesh_deformation_objects)
674  for (const auto &p : object_iterator.second)
675  if (Plugins::plugin_type_matches<MeshDeformationType>(*p))
676  return true;
677 
678  return false;
679  }
680 
681 
682 
683  template <int dim>
684  template <typename MeshDeformationType, typename>
685  inline
686  const MeshDeformationType &
688  {
689  AssertThrow(has_matching_mesh_deformation_object<MeshDeformationType> (),
690  ExcMessage("You asked MeshDeformation::MeshDeformationHandler::get_matching_mesh_deformation_object() for a "
691  "mesh deformation object of type <" + boost::core::demangle(typeid(MeshDeformationType).name()) + "> "
692  "that could not be found in the current model. Activate this "
693  "mesh deformation in the input file."));
694 
695  for (const auto &object_iterator : mesh_deformation_objects)
696  for (const auto &p : object_iterator.second)
697  if (Plugins::plugin_type_matches<MeshDeformationType>(*p))
698  return Plugins::get_plugin_as_type<MeshDeformationType>(*p);
699 
700  // We will never get here, because we had the Assert above. Just to avoid warnings.
701  typename std::vector<std::unique_ptr<Interface<dim>>>::const_iterator mesh_def;
702  return Plugins::get_plugin_as_type<MeshDeformationType>(*(*mesh_def));
703 
704  }
705 
706 
713  template <int dim>
714  std::string
716 
717 
718 
726 #define ASPECT_REGISTER_MESH_DEFORMATION_MODEL(classname,name,description) \
727  template class classname<2>; \
728  template class classname<3>; \
729  namespace ASPECT_REGISTER_MESH_DEFORMATION_MODEL_ ## classname \
730  { \
731  aspect::internal::Plugins::RegisterHelper<aspect::MeshDeformation::Interface<2>,classname<2>> \
732  dummy_ ## classname ## _2d (&aspect::MeshDeformation::MeshDeformationHandler<2>::register_mesh_deformation, \
733  name, description); \
734  aspect::internal::Plugins::RegisterHelper<aspect::MeshDeformation::Interface<3>,classname<3>> \
735  dummy_ ## classname ## _3d (&aspect::MeshDeformation::MeshDeformationHandler<3>::register_mesh_deformation, \
736  name, description); \
737  }
738  }
739 }
740 
741 #endif
virtual void parse_parameters(ParameterHandler &prm)
const Simulator< dim > * simulator
::TrilinosWrappers::MPI::BlockVector BlockVector
Definition: global.h:269
ApplyStabilization(const double stabilization_theta)
std::set< types::boundary_id > prescribed_mesh_deformation_boundary_indicators
Definition: interface.h:567
MGLevelObject<::LinearAlgebra::distributed::Vector< double > > level_displacements
Definition: interface.h:613
void write_plugin_graph(std::ostream &output_stream)
std::set< types::boundary_id > boundary_indicators_requiring_stabilization
Definition: interface.h:594
::TrilinosWrappers::MPI::Vector Vector
Definition: global.h:263
AffineConstraints< double > mesh_velocity_constraints
Definition: interface.h:542
std::set< types::boundary_id > tangential_mesh_deformation_boundary_indicators
Definition: interface.h:573
AffineConstraints< double > mesh_vertex_constraints
Definition: interface.h:548
MGLevelObject< std::unique_ptr< Mapping< dim > > > level_mappings
Definition: interface.h:608
std::map< types::boundary_id, std::vector< std::string > > mesh_deformation_object_names
Definition: interface.h:560
LinearAlgebra::BlockVector mesh_velocity
Definition: interface.h:493
std::map< types::boundary_id, std::vector< std::unique_ptr< Interface< dim > > > > mesh_deformation_objects
Definition: interface.h:554
MGTransferMF< dim, double > mg_transfer
Definition: interface.h:618
std::string get_valid_model_names_pattern()
::MGTransferMatrixFree< dim, NumberType > MGTransferMF
Definition: compat.h:64
void load(Archive &ar, const unsigned int version)
Definition: interface.h:649
void execute(internal::Assembly::Scratch::ScratchBase< dim > &scratch, internal::Assembly::CopyData::CopyDataBase< dim > &data) const override
std::set< types::boundary_id > free_surface_boundary_indicators
Definition: interface.h:588
void save(Archive &ar, const unsigned int version) const
Definition: interface.h:632
std::set< types::boundary_id > zero_mesh_deformation_boundary_indicators
Definition: interface.h:581
static void declare_parameters(ParameterHandler &prm)
DeclException1(ProbabilityFunctionNegative, Point< dim >,<< "Your probability density function in the particle generator " "returned a negative probability density for the following position: "<< arg1<< ". Please check your function expression.")