ASPECT
melt.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2016 - 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_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) override;
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) override;
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  void
215  create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
216  };
217 
222  template <int dim>
224  {
225  public:
226  void
227  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
228  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
229  };
230 
235  template <int dim>
236  class MeltStokesSystem : public MeltInterface<dim>
237  {
238  public:
239  void
240  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
241  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
242  };
243 
244 
250  template <int dim>
252  {
253  public:
254  void
255  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
256  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
257  };
258 
263  template <int dim>
264  class MeltAdvectionSystem : public MeltInterface<dim>
265  {
266  public:
267  void
268  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
269  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
270 
275  std::vector<double>
276  compute_residual(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base) const override;
277  };
278 
284  template <int dim>
286  {
287  public:
288  void
289  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
290  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
291  };
292 
296  template <int dim>
298  {
299  public:
300  void
301  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
302  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
303  };
304  }
305 
306 
307  namespace Melt
308  {
309  template <int dim>
310  struct Parameters
311  {
316  static void declare_parameters (ParameterHandler &prm);
317 
327  void parse_parameters (ParameterHandler &prm);
328 
340 
350 
355 
365  };
366  }
367 
368 
375  template <int dim>
376  class MeltHandler: public SimulatorAccess<dim>
377  {
378  public:
380 
387  static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
388 
394  void edit_finite_element_variables(const Parameters<dim> &parameters,
395  std::vector<VariableDeclaration<dim> > &variables);
396 
402  void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
403 
404 
411  void initialize() const;
412 
416  void initialize_simulator (const Simulator<dim> &simulator_object) override;
417 
432  void compute_melt_variables(LinearAlgebra::BlockSparseMatrix &system_matrix,
433  LinearAlgebra::BlockVector &solution,
434  LinearAlgebra::BlockVector &system_rhs);
435 
439  bool is_porosity (const typename Simulator<dim>::AdvectionField &advection_field) const;
440 
445  void apply_free_surface_stabilization_with_melt (const double free_surface_theta,
446  const typename DoFHandler<dim>::active_cell_iterator &cell,
449 
456  void add_current_constraints(AffineConstraints<double> &constraints);
457 
463  bool is_melt_cell(const typename DoFHandler<dim>::active_cell_iterator &cell) const;
464 
472  double limited_darcy_coefficient(const double K_D,
473  const bool is_melt_cell) const;
474 
479  get_boundary_fluid_pressure () const;
480 
486 
487  private:
493  const std::unique_ptr<aspect::BoundaryFluidPressure::Interface<dim> > boundary_fluid_pressure;
494 
502  std::vector<bool> is_melt_cell_vector;
503 
511 
512  };
513 
514 }
515 
516 #endif
bool average_melt_velocity
Definition: melt.h:364
std::vector< bool > is_melt_cell_vector
Definition: melt.h:502
Melt::Parameters< dim > melt_parameters
Definition: melt.h:485
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:349
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
double melt_scaling_factor_threshold
Definition: melt.h:339
const std::unique_ptr< aspect::BoundaryFluidPressure::Interface< dim > > boundary_fluid_pressure
Definition: melt.h:493
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:354
MeltOutputs(const unsigned int n_points, const unsigned int)
Definition: melt.h:79
AffineConstraints< double > current_constraints
Definition: melt.h:510