ASPECT
melt.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2016 - 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_melt_h
23 #define _aspect_melt_h
24 
26 #include <aspect/global.h>
30 
31 namespace aspect
32 {
33  using namespace dealii;
34 
35  namespace MaterialModel
36  {
42  template <int dim>
44  {
45  public:
54  MeltInputs (const unsigned int n_points);
55 
59  std::vector<double> compaction_pressures;
60 
65  std::vector<Tensor<1,dim> > fluid_velocities;
66 
70  void fill (const LinearAlgebra::BlockVector &solution,
71  const FEValuesBase<dim> &fe_values,
72  const Introspection<dim> &introspection);
73  };
74 
75  template <int dim>
77  {
78  public:
79  MeltOutputs (const unsigned int n_points,
80  const unsigned int /*n_comp*/)
81  {
82  compaction_viscosities.resize(n_points);
83  fluid_viscosities.resize(n_points);
84  permeabilities.resize(n_points);
85  fluid_densities.resize(n_points);
86  fluid_density_gradients.resize(n_points, Tensor<1,dim>());
87  }
88 
94  std::vector<double> compaction_viscosities;
95 
99  std::vector<double> fluid_viscosities;
100 
104  std::vector<double> permeabilities;
105 
109  std::vector<double> fluid_densities;
110 
117  std::vector<Tensor<1,dim> > fluid_density_gradients;
118 
124  void average (const MaterialAveraging::AveragingOperation operation,
125  const FullMatrix<double> &projection_matrix,
126  const FullMatrix<double> &expansion_matrix);
127  };
128 
133  template <int dim>
135  {
136  public:
145  virtual void melt_fractions (const MaterialModel::MaterialModelInputs<dim> &in,
146  std::vector<double> &melt_fractions) const = 0;
147 
153  {};
154  };
155 
159  template <int dim>
161  {
162  public:
167  virtual double reference_darcy_coefficient () const = 0;
168 
186  double p_c_scale (const MaterialModel::MaterialModelInputs<dim> &inputs,
188  const MeltHandler<dim> &melt_handler,
189  const bool consider_is_melt_cell) const;
190  };
191 
192 
193  }
194 
195 
196 
197  namespace Assemblers
198  {
204  template <int dim>
206  public SimulatorAccess<dim>
207  {
208  public:
214  virtual
215  void
216  create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const;
217  };
218 
223  template <int dim>
225  {
226  public:
227  virtual
228  void
229  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
231  };
232 
237  template <int dim>
238  class MeltStokesSystem : public MeltInterface<dim>
239  {
240  public:
241  virtual
242  void
243  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
245  };
246 
247 
253  template <int dim>
255  {
256  public:
257  virtual
258  void
259  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
261  };
262 
267  template <int dim>
268  class MeltAdvectionSystem : public MeltInterface<dim>
269  {
270  public:
271  virtual
272  void
273  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
275 
280  virtual
281  std::vector<double>
282  compute_residual(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base) const;
283  };
284 
290  template <int dim>
292  {
293  public:
294  virtual
295  void
296  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
298  };
299 
303  template <int dim>
305  {
306  public:
307  virtual
308  void
309  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
311  };
312  }
313 
314 
315  namespace Melt
316  {
317  template <int dim>
318  struct Parameters
319  {
324  static void declare_parameters (ParameterHandler &prm);
325 
335  void parse_parameters (ParameterHandler &prm);
336 
348 
358 
363 
373  };
374  }
375 
376 
383  template <int dim>
384  class MeltHandler: public SimulatorAccess<dim>
385  {
386  public:
388 
395  static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
396 
402  void edit_finite_element_variables(const Parameters<dim> &parameters,
403  std::vector<VariableDeclaration<dim> > &variables);
404 
410  void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
411 
412 
419  void initialize() const;
420 
424  void initialize_simulator (const Simulator<dim> &simulator_object);
425 
436  void compute_melt_variables(LinearAlgebra::BlockVector &solution);
437 
441  bool is_porosity (const typename Simulator<dim>::AdvectionField &advection_field) const;
442 
447  void apply_free_surface_stabilization_with_melt (const double free_surface_theta,
448  const typename DoFHandler<dim>::active_cell_iterator &cell,
451 
458  void add_current_constraints(ConstraintMatrix &constraints);
459 
465  bool is_melt_cell(const typename DoFHandler<dim>::active_cell_iterator &cell) const;
466 
474  double limited_darcy_coefficient(const double K_D,
475  const bool is_melt_cell) const;
476 
481  get_boundary_fluid_pressure () const;
482 
488 
489  private:
495  const std::unique_ptr<aspect::BoundaryFluidPressure::Interface<dim> > boundary_fluid_pressure;
496 
504  std::vector<bool> is_melt_cell_vector;
505 
512  ConstraintMatrix current_constraints;
513 
514  };
515 
516 }
517 
518 #endif
bool average_melt_velocity
Definition: melt.h:372
std::vector< bool > is_melt_cell_vector
Definition: melt.h:504
Melt::Parameters< dim > melt_parameters
Definition: melt.h:487
void average(const AveragingOperation operation, const typename DoFHandler< dim >::active_cell_iterator &cell, const Quadrature< dim > &quadrature_formula, const Mapping< dim > &mapping, MaterialModelOutputs< dim > &values_out)
std::vector< Tensor< 1, dim > > fluid_density_gradients
Definition: melt.h:117
void declare_parameters(ParameterHandler &prm)
bool heat_advection_by_melt
Definition: melt.h:357
std::vector< Tensor< 1, dim > > fluid_velocities
Definition: melt.h:65
std::vector< double > fluid_viscosities
Definition: melt.h:99
std::vector< double > compaction_viscosities
Definition: melt.h:94
std::vector< double > permeabilities
Definition: melt.h:104
Definition: compat.h:37
double melt_scaling_factor_threshold
Definition: melt.h:347
const std::unique_ptr< aspect::BoundaryFluidPressure::Interface< dim > > boundary_fluid_pressure
Definition: melt.h:495
std::vector< double > compaction_pressures
Definition: melt.h:59
std::vector< double > fluid_densities
Definition: melt.h:109
bool use_discontinuous_p_c
Definition: melt.h:362
MeltOutputs(const unsigned int n_points, const unsigned int)
Definition: melt.h:79
ConstraintMatrix current_constraints
Definition: melt.h:512