ASPECT
free_surface.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_mesh_deformation_free_surface_h
23 #define _aspect_mesh_deformation_free_surface_h
24 
26 
29 
30 
31 namespace aspect
32 {
33  using namespace dealii;
34 
35  namespace Assemblers
36  {
43  template <int dim>
45  public SimulatorAccess<dim>
46  {
47  public:
48  ApplyStabilization(const double stabilization_theta);
49 
50  virtual ~ApplyStabilization () {};
51 
52  virtual
53  void
56 
57  private:
63  const double free_surface_theta;
64  };
65  }
66 
67  namespace MeshDeformation
68  {
75  template<int dim>
76  class FreeSurface : public Interface<dim>, public SimulatorAccess<dim>
77  {
78  public:
83  virtual void initialize();
84 
89  void set_assemblers(const SimulatorAccess<dim> &simulator_access,
90  aspect::Assemblers::Manager<dim> &assemblers) const;
91 
98  virtual
99  void
100  compute_velocity_constraints_on_boundary(const DoFHandler<dim> &mesh_deformation_dof_handler,
101  ConstraintMatrix &mesh_velocity_constraints,
102  const std::set<types::boundary_id> &boundary_id) const;
103 
107  static
109 
113  void parse_parameters (ParameterHandler &prm);
114 
115  private:
120  void project_velocity_onto_boundary (const DoFHandler<dim> &free_surface_dof_handler,
121  const IndexSet &mesh_locally_owned,
122  const IndexSet &mesh_locally_relevant,
123  LinearAlgebra::Vector &output) const;
124 
131 
136  {
137  enum Direction { normal, vertical };
138  };
139 
145  };
146  }
147 }
148 
149 
150 #endif
void declare_parameters(ParameterHandler &prm)
Definition: compat.h:37
SurfaceAdvection::Direction advection_direction
Definition: free_surface.h:144