ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2017 - 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 #ifndef _aspect_simulator_assemblers_interface_h
22 #define _aspect_simulator_assemblers_interface_h
23 
24 #include <aspect/global.h>
27 
28 #include <deal.II/fe/fe_values.h>
29 
30 namespace aspect
31 {
32  using namespace dealii;
33 
34  template <int dim>
35  class Simulator;
36 
37  struct AdvectionField;
38 
46  namespace internal
47  {
48  namespace Assembly
49  {
50  namespace Scratch
51  {
58  template <int dim>
59  struct ScratchBase
60  {
62  :
63  cell(),
64  face_number(numbers::invalid_unsigned_int)
65  {}
66 
67  ScratchBase(const ScratchBase &scratch)
68  :
69  cell(scratch.cell),
70  face_number(scratch.face_number)
71  {}
72 
73  virtual ~ScratchBase () = default;
74 
79 
86  unsigned face_number;
87  };
88 
93  template <int dim>
94  struct StokesPreconditioner: public ScratchBase<dim>
95  {
96  StokesPreconditioner (const FiniteElement<dim> &finite_element,
97  const Quadrature<dim> &quadrature,
98  const Mapping<dim> &mapping,
99  const UpdateFlags update_flags,
100  const unsigned int n_compositional_fields,
101  const unsigned int stokes_dofs_per_cell,
102  const bool add_compaction_pressure,
103  const bool rebuild_matrix);
104  StokesPreconditioner (const StokesPreconditioner &scratch);
105 
106  ~StokesPreconditioner () override;
107 
109 
110  void reinit (const typename DoFHandler<dim>::active_cell_iterator &cell_ref);
111 
112  std::vector<types::global_dof_index> local_dof_indices;
113  std::vector<unsigned int> dof_component_indices;
114  std::vector<SymmetricTensor<2,dim> > grads_phi_u;
115  std::vector<double> div_phi_u;
116  std::vector<double> phi_p;
117  std::vector<double> phi_p_c;
118  std::vector<Tensor<1,dim> > grad_phi_p;
119 
126 
133  };
134 
145  template <int dim>
146  struct StokesSystem : public StokesPreconditioner<dim>
147  {
148  StokesSystem (const FiniteElement<dim> &finite_element,
149  const Mapping<dim> &mapping,
150  const Quadrature<dim> &quadrature,
151  const Quadrature<dim-1> &face_quadrature,
152  const UpdateFlags update_flags,
153  const UpdateFlags face_update_flags,
154  const unsigned int n_compositional_fields,
155  const unsigned int stokes_dofs_per_cell,
156  const bool add_compaction_pressure,
157  const bool use_reference_density_profile,
158  const bool rebuild_stokes_matrix,
159  const bool rebuild_newton_stokes_matrix);
160 
161  StokesSystem (const StokesSystem<dim> &scratch);
162 
164 
166 
167  void reinit (const typename DoFHandler<dim>::active_cell_iterator &cell_ref,
168  const unsigned face_number_ref);
169 
170  std::vector<Tensor<1,dim> > phi_u;
171  std::vector<Tensor<1,dim> > velocity_values;
172  std::vector<double> velocity_divergence;
173  std::vector<Tensor<1,dim> > temperature_gradients;
174 
185 
195  std::vector<double> reference_densities;
197 
204  };
205 
210  template <int dim>
211  struct AdvectionSystem: public ScratchBase<dim>
212  {
213  AdvectionSystem (const FiniteElement<dim> &finite_element,
214  const FiniteElement<dim> &advection_element,
215  const Mapping<dim> &mapping,
216  const Quadrature<dim> &quadrature,
217  const Quadrature<dim-1> &face_quadrature,
218  const UpdateFlags update_flags,
219  const UpdateFlags face_update_flags,
220  const unsigned int n_compositional_fields,
221  const typename Simulator<dim>::AdvectionField &field);
222  AdvectionSystem (const AdvectionSystem &scratch);
223 
225 
226  void reinit (const typename DoFHandler<dim>::active_cell_iterator &cell_ref);
227 
228  std::unique_ptr<FEFaceValues<dim> > face_finite_element_values;
229  std::unique_ptr<FEFaceValues<dim> > neighbor_face_finite_element_values;
230  std::unique_ptr<FESubfaceValues<dim> > subface_finite_element_values;
231 
232  std::vector<types::global_dof_index> local_dof_indices;
233 
242  std::vector<double> phi_field;
243  std::vector<Tensor<1,dim> > grad_phi_field;
244  std::vector<double> laplacian_phi_field;
245  std::vector<double> face_phi_field;
246  std::vector<Tensor<1,dim> > face_grad_phi_field;
247  std::vector<double> neighbor_face_phi_field;
248  std::vector<Tensor<1,dim> > neighbor_face_grad_phi_field;
249 
250  std::vector<Tensor<1,dim> > old_velocity_values;
251  std::vector<Tensor<1,dim> > old_old_velocity_values;
252 
253  std::vector<double> old_pressure;
254  std::vector<double> old_old_pressure;
255  std::vector<Tensor<1,dim> > old_pressure_gradients;
256  std::vector<Tensor<1,dim> > old_old_pressure_gradients;
257 
258  std::vector<SymmetricTensor<2,dim> > old_strain_rates;
259  std::vector<SymmetricTensor<2,dim> > old_old_strain_rates;
260 
261  std::vector<double> old_temperature_values;
262  std::vector<double> old_old_temperature_values;
263 
264  std::vector<double> old_field_values;
265  std::vector<double> old_old_field_values;
266  std::vector<Tensor<1,dim> > old_field_grads;
267  std::vector<Tensor<1,dim> > old_old_field_grads;
268  std::vector<double> old_field_laplacians;
269  std::vector<double> old_old_field_laplacians;
270 
271  std::vector<std::vector<double> > old_composition_values;
272  std::vector<std::vector<double> > old_old_composition_values;
273 
274  std::vector<double> current_temperature_values;
275  std::vector<Tensor<1,dim> > current_velocity_values;
276  std::vector<Tensor<1,dim> > face_current_velocity_values;
277  std::vector<Tensor<1,dim> > mesh_velocity_values;
278  std::vector<Tensor<1,dim> > face_mesh_velocity_values;
279 
280  std::vector<SymmetricTensor<2,dim> > current_strain_rates;
281  std::vector<std::vector<double> > current_composition_values;
282  std::vector<double> current_velocity_divergences;
283 
290 
293 
296 
307 
316 
322  };
323  }
324 
337  namespace CopyData
338  {
344  template <int dim>
346  {
347  virtual ~CopyDataBase () = default;
348  };
349 
355  template <int dim>
356  struct StokesPreconditioner: public CopyDataBase<dim>
357  {
358  StokesPreconditioner (const unsigned int stokes_dofs_per_cell);
359 
360  StokesPreconditioner (const StokesPreconditioner &data);
361 
362  ~StokesPreconditioner () override = default;
363  StokesPreconditioner<dim> &operator= (const StokesPreconditioner<dim> &data) = default;
364 
366  std::vector<types::global_dof_index> local_dof_indices;
367 
373  void extract_stokes_dof_indices(const std::vector<types::global_dof_index> &all_dof_indices,
374  const Introspection<dim> &introspection,
375  const FiniteElement<dim> &finite_element);
376  };
377 
383  template <int dim>
384  struct StokesSystem : public StokesPreconditioner<dim>
385  {
386  StokesSystem (const unsigned int stokes_dofs_per_cell,
387  const bool do_pressure_rhs_compatibility_modification);
388  StokesSystem (const StokesSystem<dim> &data);
389 
390  ~StokesSystem () override = default;
391  StokesSystem<dim> &operator= (const StokesSystem<dim> &data) = default;
392 
395  };
396 
402  template <int dim>
403  struct AdvectionSystem: public CopyDataBase<dim>
404  {
413  AdvectionSystem (const FiniteElement<dim> &finite_element,
414  const bool field_is_discontinuous);
415 
421 
432  std::vector<FullMatrix<double> > local_matrices_int_ext;
433  std::vector<FullMatrix<double> > local_matrices_ext_int;
434  std::vector<FullMatrix<double> > local_matrices_ext_ext;
435 
441 
447  std::vector<bool> assembled_matrices;
448 
458  std::vector<types::global_dof_index> local_dof_indices;
459 
467  std::vector<std::vector<types::global_dof_index> > neighbor_dof_indices;
468  };
469  }
470  }
471  }
472 
477  namespace Assemblers
478  {
493  template <int dim>
494  class Interface
495  {
496  public:
497  virtual ~Interface ();
498 
512  virtual
513  void
516 
550  virtual void create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &) const;
551 
565  virtual
566  std::vector<double>
567  compute_residual(internal::Assembly::Scratch::ScratchBase<dim> &) const;
568  };
569 
596  template <int dim>
597  class Manager
598  {
599  public:
600 
604  void reset ();
605 
610  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > stokes_preconditioner;
611 
617  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > stokes_system;
618 
626  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > stokes_system_on_boundary_face;
627 
632  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > advection_system;
633 
641  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > advection_system_on_boundary_face;
642 
650  std::vector<std::unique_ptr<Assemblers::Interface<dim> > > advection_system_on_interior_face;
651 
671  struct Properties
672  {
677  Properties ();
678 
687 
696 
702 
710  };
711 
721  std::vector<Properties> advection_system_assembler_properties;
723  };
724  }
725 }
726 
727 
728 #endif
std::vector< SymmetricTensor< 2, dim > > old_strain_rates
Definition: interface.h:258
std::vector< Tensor< 1, dim > > old_field_grads
Definition: interface.h:266
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
std::vector< types::global_dof_index > local_dof_indices
Definition: interface.h:458
std::vector< Tensor< 1, dim > > old_velocity_values
Definition: interface.h:250
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > stokes_system
Definition: interface.h:617
Properties stokes_preconditioner_assembler_properties
Definition: interface.h:718
MaterialModel::MaterialModelInputs< dim > face_material_model_inputs
Definition: interface.h:183
std::vector< Tensor< 1, dim > > face_mesh_velocity_values
Definition: interface.h:278
std::vector< Tensor< 1, dim > > grad_phi_field
Definition: interface.h:243
std::vector< Tensor< 1, dim > > face_grad_phi_field
Definition: interface.h:246
MaterialModel::MaterialModelOutputs< dim > material_model_outputs
Definition: interface.h:289
MaterialModel::MaterialModelInputs< dim > material_model_inputs
Definition: interface.h:288
HeatingModel::HeatingModelOutputs face_heating_model_outputs
Definition: interface.h:305
MaterialModel::MaterialModelOutputs< dim > neighbor_face_material_model_outputs
Definition: interface.h:295
std::vector< Tensor< 1, dim > > mesh_velocity_values
Definition: interface.h:277
const Simulator< dim >::AdvectionField * advection_field
Definition: interface.h:315
std::vector< Tensor< 1, dim > > old_pressure_gradients
Definition: interface.h:255
std::vector< SymmetricTensor< 2, dim > > old_old_strain_rates
Definition: interface.h:259
std::vector< Tensor< 1, dim > > old_old_field_grads
Definition: interface.h:267
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > advection_system_on_interior_face
Definition: interface.h:650
std::vector< Tensor< 1, dim > > temperature_gradients
Definition: interface.h:173
std::vector< SymmetricTensor< 2, dim > > current_strain_rates
Definition: interface.h:280
ScratchBase(const ScratchBase &scratch)
Definition: interface.h:67
UpdateFlags
std::vector< Tensor< 1, dim > > current_velocity_values
Definition: interface.h:275
std::vector< Properties > advection_system_assembler_properties
Definition: interface.h:721
MaterialModel::MaterialModelOutputs< dim > face_material_model_outputs
Definition: interface.h:292
HeatingModel::HeatingModelOutputs neighbor_face_heating_model_outputs
Definition: interface.h:306
std::vector< Tensor< 1, dim > > neighbor_face_grad_phi_field
Definition: interface.h:248
std::unique_ptr< FEFaceValues< dim > > neighbor_face_finite_element_values
Definition: interface.h:229
std::vector< std::vector< types::global_dof_index > > neighbor_dof_indices
Definition: interface.h:467
std::vector< FullMatrix< double > > local_matrices_int_ext
Definition: interface.h:432
std::vector< Tensor< 1, dim > > old_old_velocity_values
Definition: interface.h:251
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > advection_system_on_boundary_face
Definition: interface.h:641
std::vector< FullMatrix< double > > local_matrices_ext_int
Definition: interface.h:433
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > advection_system
Definition: interface.h:632
std::vector< std::vector< double > > old_composition_values
Definition: interface.h:271
std::vector< SymmetricTensor< 2, dim > > grads_phi_u
Definition: interface.h:114
Properties stokes_system_assembler_on_boundary_face_properties
Definition: interface.h:720
DoFHandler< dim >::active_cell_iterator cell
Definition: interface.h:78
std::vector< std::vector< double > > old_old_composition_values
Definition: interface.h:272
MaterialModel::MaterialModelOutputs< dim > face_material_model_outputs
Definition: interface.h:184
std::vector< Tensor< 1, dim > > old_old_pressure_gradients
Definition: interface.h:256
std::vector< types::global_dof_index > local_dof_indices
Definition: interface.h:366
std::vector< Tensor< 1, dim > > velocity_values
Definition: interface.h:171
HeatingModel::HeatingModelOutputs heating_model_outputs
Definition: interface.h:304
std::vector< Properties > advection_system_assembler_on_face_properties
Definition: interface.h:722
std::vector< double > reference_densities_depth_derivative
Definition: interface.h:196
MaterialModel::MaterialModelInputs< dim > neighbor_face_material_model_inputs
Definition: interface.h:294
std::unique_ptr< FESubfaceValues< dim > > subface_finite_element_values
Definition: interface.h:230
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > stokes_system_on_boundary_face
Definition: interface.h:626
std::vector< Tensor< 1, dim > > face_current_velocity_values
Definition: interface.h:276
std::vector< FullMatrix< double > > local_matrices_ext_ext
Definition: interface.h:434
std::vector< Tensor< 1, dim > > phi_u
Definition: interface.h:170
MaterialModel::MaterialModelInputs< dim > face_material_model_inputs
Definition: interface.h:291
MaterialModel::MaterialModelInputs< dim > material_model_inputs
Definition: interface.h:124
std::vector< std::vector< double > > current_composition_values
Definition: interface.h:281
std::unique_ptr< FEFaceValues< dim > > face_finite_element_values
Definition: interface.h:228
std::vector< std::unique_ptr< Assemblers::Interface< dim > > > stokes_preconditioner
Definition: interface.h:610
MaterialModel::MaterialModelOutputs< dim > material_model_outputs
Definition: interface.h:125
Properties stokes_system_assembler_properties
Definition: interface.h:719
std::vector< types::global_dof_index > local_dof_indices
Definition: interface.h:232
std::vector< types::global_dof_index > local_dof_indices
Definition: interface.h:112