ASPECT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
newton.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_newton_h
23 #define _aspect_newton_h
24 
27 #include <aspect/global.h>
29 
30 #include <deal.II/base/signaling_nan.h>
31 
32 namespace aspect
33 {
34  namespace MaterialModel
35  {
39  template <int dim>
41  {
42  public:
47  MaterialModelDerivatives (const unsigned int n_points);
48 
52  std::vector<double> viscosity_derivative_wrt_pressure;
53 
57  std::vector<SymmetricTensor<2,dim>> viscosity_derivative_wrt_strain_rate;
58 
64 
71  std::vector<double> dilation_derivative_wrt_pressure;
72 
76  std::vector<SymmetricTensor<2,dim>> dilation_derivative_wrt_strain_rate;
77  };
78  }
79 
80 
81  namespace Newton
82  {
83  struct Parameters
84  {
85 
95  {
96  none = 0,
97  symmetric = 1,
98  PD = 2,
99  SPD = symmetric | PD
100  };
101 
102 
108  friend
111  const Stabilization b)
112  {
113  return static_cast<Stabilization>(
114  static_cast<int>(a) | static_cast<int>(b));
115  }
116 
117 
123  friend
125  operator& (const Stabilization a,
126  const Stabilization b)
127  {
128  return static_cast<Stabilization>(
129  static_cast<int>(a) & static_cast<int>(b));
130  }
131 
132 
137  static void declare_parameters (ParameterHandler &prm);
138 
143  void parse_parameters (ParameterHandler &prm);
144 
161 
164 
171 
178 
185  };
186 
187 
192  std::string
193  to_string(const Newton::Parameters::Stabilization preconditioner_stabilization);
194  }
195 
196 
197 
201  template <int dim>
202  class NewtonHandler : public SimulatorAccess<dim>
203  {
204  public:
210  void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
211 
217  static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
218 
224  };
225 
226 
227  namespace Assemblers
228  {
233  template <int dim>
235  public SimulatorAccess<dim>
236  {
237  public:
243  void
244  create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
245  };
246 
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 
262  void create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
263  };
264 
269  template <int dim>
271  {
272  public:
273  void
274  execute (internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
275  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
276 
282  void create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
283  };
284 
289  template <int dim>
291  public SimulatorAccess<dim>
292  {
293  public:
294  void
295  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
296  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
297  };
298 
305  template <int dim>
307  public SimulatorAccess<dim>
308  {
309  public:
310  void
311  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
312  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
313  };
314 
323  template <int dim>
325  public SimulatorAccess<dim>
326  {
327  public:
328  void
329  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
330  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
331  };
332 
339  template <int dim>
341  public SimulatorAccess<dim>
342  {
343  public:
344  void
345  execute(internal::Assembly::Scratch::ScratchBase<dim> &scratch_base,
346  internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
347  };
348 
357  template <int dim>
359  public SimulatorAccess<dim>
360  {
361  public:
362  void
365  };
366  }
367 }
368 
369 #endif
unsigned int max_pre_newton_nonlinear_iterations
Definition: newton.h:180
bool use_newton_residual_scaling_method
Definition: newton.h:182
double maximum_linear_stokes_solver_tolerance
Definition: newton.h:183
double nonlinear_switch_tolerance
Definition: newton.h:177
std::string to_string(const Newton::Parameters::Stabilization preconditioner_stabilization)
bool use_Eisenstat_Walker_method_for_Picard_iterations
Definition: newton.h:179
void declare_parameters(ParameterHandler &prm)
Stabilization preconditioner_stabilization
Definition: newton.h:162
std::vector< SymmetricTensor< 2, dim > > viscosity_derivative_wrt_strain_rate
Definition: newton.h:57
std::vector< double > dilation_derivative_wrt_pressure
Definition: newton.h:71
Stabilization velocity_block_stabilization
Definition: newton.h:163
std::vector< SymmetricTensor< 2, dim > > dilation_derivative_wrt_strain_rate
Definition: newton.h:76
double newton_derivative_scaling_factor
Definition: newton.h:160
Newton::Parameters parameters
Definition: newton.h:223
MaterialModelDerivatives(const unsigned int n_points)
std::vector< double > viscosity_derivative_averaging_weights
Definition: newton.h:63
Dependence operator|(const Dependence d1, const Dependence d2)
Definition: interface.h:99
unsigned int max_newton_line_search_iterations
Definition: newton.h:181
std::vector< double > viscosity_derivative_wrt_pressure
Definition: newton.h:52