ASPECT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
melt.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2016 - 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_melt_h
23 #define _aspect_melt_h
24 
26 #include <aspect/global.h>
30 
31 namespace aspect
32 {
33  namespace MaterialModel
34  {
40  template <int dim>
42  {
43  public:
52  MeltInputs (const unsigned int n_points);
53 
57  std::vector<double> compaction_pressures;
58 
63  std::vector<Tensor<1,dim>> fluid_velocities;
64 
68  void fill (const LinearAlgebra::BlockVector &solution,
69  const FEValuesBase<dim> &fe_values,
70  const Introspection<dim> &introspection) override;
71  };
72 
73  template <int dim>
75  {
76  public:
85  MeltOutputs (const unsigned int n_points,
86  const unsigned int n_comp);
87 
93  std::vector<double> compaction_viscosities;
94 
98  std::vector<double> fluid_viscosities;
99 
103  std::vector<double> permeabilities;
104 
108  std::vector<double> fluid_densities;
109 
116  std::vector<Tensor<1,dim>> fluid_density_gradients;
117 
123  void average (const MaterialAveraging::AveragingOperation operation,
124  const FullMatrix<double> &projection_matrix,
125  const FullMatrix<double> &expansion_matrix) override;
126  };
127 
139  template <int dim>
141  {
142  public:
147  virtual ~MeltFractionModel () = default;
148 
161  virtual void melt_fractions (const MaterialModel::MaterialModelInputs<dim> &in,
162  std::vector<double> &melt_fractions,
163  const MaterialModel::MaterialModelOutputs<dim> *out = nullptr) const = 0;
164 
172  template <typename ModelType>
173  static
174  bool is_melt_fraction_model (const ModelType &model_object);
175 
181  template <typename ModelType>
182  static
183  const MeltFractionModel<dim> &
184  as_melt_fraction_model (const ModelType &model_object);
185  };
186 
187 
188 
189  template <int dim>
190  template <typename ModelType>
191  inline
192  bool
193  MeltFractionModel<dim>::is_melt_fraction_model (const ModelType &model_object)
194  {
195  return (dynamic_cast<const MeltFractionModel<dim>*>(&model_object)
196  != nullptr);
197  }
198 
199 
200  template <int dim>
201  template <typename ModelType>
202  inline
203  const MeltFractionModel<dim> &
204  MeltFractionModel<dim>::as_melt_fraction_model (const ModelType &model_object)
205  {
206  Assert (is_melt_fraction_model(model_object) == true,
207  ExcMessage ("This function can only be called for model objects "
208  "whose types are derived from MeltFractionModel."));
209 
210  return dynamic_cast<const MeltFractionModel<dim>&>(model_object);
211  }
212 
213 
214 
218  template <int dim>
220  {
221  public:
226  virtual double reference_darcy_coefficient () const = 0;
227 
245  double p_c_scale (const MaterialModel::MaterialModelInputs<dim> &inputs,
247  const MeltHandler<dim> &melt_handler,
248  const bool consider_is_melt_cell) const;
249  };
250 
251 
252  }
253 
254 
255 
256  namespace Assemblers
257  {
263  template <int dim>
265  public SimulatorAccess<dim>
266  {
267  public:
273  void
274  create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
275  };
276 
281  template <int dim>
283  {
284  public:
285  void
286  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
287  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
288  };
289 
294  template <int dim>
295  class MeltStokesSystem : public MeltInterface<dim>
296  {
297  public:
298  void
299  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
300  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
301  };
302 
303 
309  template <int dim>
311  {
312  public:
313  void
314  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
315  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
316  };
317 
322  template <int dim>
324  {
325  public:
326  void
327  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
328  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
329 
334  std::vector<double>
335  compute_residual(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base) const override;
336  };
337 
343  template <int dim>
345  {
346  public:
347  void
348  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
349  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
350  };
351 
355  template <int dim>
357  {
358  public:
359  void
360  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
361  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
362  };
363  }
364 
365 
366  namespace Melt
367  {
368  template <int dim>
369  struct Parameters
370  {
375  static void declare_parameters (ParameterHandler &prm);
376 
386  void parse_parameters (ParameterHandler &prm);
387 
399 
409 
414 
424  };
425  }
426 
427 
434  template <int dim>
435  class MeltHandler: public SimulatorAccess<dim>
436  {
437  public:
438  MeltHandler(ParameterHandler &prm);
439 
446  static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
447 
453  void edit_finite_element_variables(const Parameters<dim> &parameters,
454  std::vector<VariableDeclaration<dim>> &variables);
455 
461  void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
462 
463 
470  void initialize() const;
471 
475  void initialize_simulator (const Simulator<dim> &simulator_object) override;
476 
491  void compute_melt_variables(LinearAlgebra::BlockSparseMatrix &system_matrix,
492  LinearAlgebra::BlockVector &solution,
493  LinearAlgebra::BlockVector &system_rhs) const;
494 
498  bool is_porosity (const AdvectionField &advection_field) const;
499 
504  void apply_free_surface_stabilization_with_melt (const double free_surface_theta,
505  const typename DoFHandler<dim>::active_cell_iterator &cell,
508 
515  void add_current_constraints(AffineConstraints<double> &constraints);
516 
522  bool is_melt_cell(const typename DoFHandler<dim>::active_cell_iterator &cell) const;
523 
531  double limited_darcy_coefficient(const double K_D,
532  const bool is_melt_cell) const;
533 
538  get_boundary_fluid_pressure () const;
539 
545 
546  private:
552  const std::unique_ptr<aspect::BoundaryFluidPressure::Interface<dim>> boundary_fluid_pressure;
553 
561  std::vector<bool> is_melt_cell_vector;
562 
569  AffineConstraints<double> current_constraints;
570 
571  };
572 
573 }
574 
575 #endif
bool average_melt_velocity
Definition: melt.h:423
::TrilinosWrappers::MPI::BlockVector BlockVector
Definition: global.h:269
std::vector< bool > is_melt_cell_vector
Definition: melt.h:561
Melt::Parameters< dim > melt_parameters
Definition: melt.h:544
void declare_parameters(ParameterHandler &prm)
static const MeltFractionModel< dim > & as_melt_fraction_model(const ModelType &model_object)
Definition: melt.h:204
std::vector< Tensor< 1, dim > > fluid_density_gradients
Definition: melt.h:116
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
Definition: global.h:280
bool heat_advection_by_melt
Definition: melt.h:408
std::vector< Tensor< 1, dim > > fluid_velocities
Definition: melt.h:63
std::vector< double > fluid_viscosities
Definition: melt.h:98
std::vector< double > compaction_viscosities
Definition: melt.h:93
std::vector< double > permeabilities
Definition: melt.h:103
double melt_scaling_factor_threshold
Definition: melt.h:398
const std::unique_ptr< aspect::BoundaryFluidPressure::Interface< dim > > boundary_fluid_pressure
Definition: melt.h:552
MeltInputs(const unsigned int n_points)
void fill(const LinearAlgebra::BlockVector &solution, const FEValuesBase< dim > &fe_values, const Introspection< dim > &introspection) override
std::vector< double > compaction_pressures
Definition: melt.h:57
static bool is_melt_fraction_model(const ModelType &model_object)
Definition: melt.h:193
void average(const AveragingOperation operation, const typename DoFHandler< dim >::active_cell_iterator &cell, const Quadrature< dim > &quadrature_formula, const Mapping< dim > &mapping, const MaterialProperties::Property &requested_properties, MaterialModelOutputs< dim > &values_out)
std::vector< double > fluid_densities
Definition: melt.h:108
bool use_discontinuous_p_c
Definition: melt.h:413
AffineConstraints< double > current_constraints
Definition: melt.h:569