ASPECT
melt.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2016 - 2021 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 
152  virtual ~MeltFractionModel () = default;
153  };
154 
158  template <int dim>
160  {
161  public:
166  virtual double reference_darcy_coefficient () const = 0;
167 
185  double p_c_scale (const MaterialModel::MaterialModelInputs<dim> &inputs,
187  const MeltHandler<dim> &melt_handler,
188  const bool consider_is_melt_cell) const;
189  };
190 
191 
192  }
193 
194 
195 
196  namespace Assemblers
197  {
203  template <int dim>
205  public SimulatorAccess<dim>
206  {
207  public:
213  void
214  create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
215  };
216 
221  template <int dim>
223  {
224  public:
225  void
226  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
227  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
228  };
229 
234  template <int dim>
235  class MeltStokesSystem : public MeltInterface<dim>
236  {
237  public:
238  void
239  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
240  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
241  };
242 
243 
249  template <int dim>
251  {
252  public:
253  void
254  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
255  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
256  };
257 
262  template <int dim>
264  {
265  public:
266  void
267  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
268  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
269 
274  std::vector<double>
275  compute_residual(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base) const override;
276  };
277 
283  template <int dim>
285  {
286  public:
287  void
288  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
289  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
290  };
291 
295  template <int dim>
297  {
298  public:
299  void
300  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
301  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
302  };
303  }
304 
305 
306  namespace Melt
307  {
308  template <int dim>
309  struct Parameters
310  {
315  static void declare_parameters (ParameterHandler &prm);
316 
326  void parse_parameters (ParameterHandler &prm);
327 
339 
349 
354 
364  };
365  }
366 
367 
374  template <int dim>
375  class MeltHandler: public SimulatorAccess<dim>
376  {
377  public:
378  MeltHandler(ParameterHandler &prm);
379 
386  static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
387 
393  void edit_finite_element_variables(const Parameters<dim> &parameters,
394  std::vector<VariableDeclaration<dim>> &variables);
395 
401  void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
402 
403 
410  void initialize() const;
411 
415  void initialize_simulator (const Simulator<dim> &simulator_object) override;
416 
431  void compute_melt_variables(LinearAlgebra::BlockSparseMatrix &system_matrix,
432  LinearAlgebra::BlockVector &solution,
433  LinearAlgebra::BlockVector &system_rhs);
434 
438  bool is_porosity (const typename Simulator<dim>::AdvectionField &advection_field) const;
439 
444  void apply_free_surface_stabilization_with_melt (const double free_surface_theta,
445  const typename DoFHandler<dim>::active_cell_iterator &cell,
448 
455  void add_current_constraints(AffineConstraints<double> &constraints);
456 
462  bool is_melt_cell(const typename DoFHandler<dim>::active_cell_iterator &cell) const;
463 
471  double limited_darcy_coefficient(const double K_D,
472  const bool is_melt_cell) const;
473 
478  get_boundary_fluid_pressure () const;
479 
485 
486  private:
492  const std::unique_ptr<aspect::BoundaryFluidPressure::Interface<dim>> boundary_fluid_pressure;
493 
501  std::vector<bool> is_melt_cell_vector;
502 
509  AffineConstraints<double> current_constraints;
510 
511  };
512 
513 }
514 
515 #endif
bool average_melt_velocity
Definition: melt.h:363
::TrilinosWrappers::MPI::BlockVector BlockVector
Definition: global.h:239
std::vector< bool > is_melt_cell_vector
Definition: melt.h:501
Melt::Parameters< dim > melt_parameters
Definition: melt.h:484
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
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
Definition: global.h:250
void declare_parameters(ParameterHandler &prm)
bool heat_advection_by_melt
Definition: melt.h:348
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:338
const std::unique_ptr< aspect::BoundaryFluidPressure::Interface< dim > > boundary_fluid_pressure
Definition: melt.h:492
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:353
MeltOutputs(const unsigned int n_points, const unsigned int)
Definition: melt.h:79
Definition: compat.h:42
AffineConstraints< double > current_constraints
Definition: melt.h:509