newton.h

Go to the documentation of this file.

/*
  Copyright (C) 2016 - 2024 by the authors of the ASPECT code.
 
  This file is part of ASPECT.
 
  ASPECT is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; either version 2, or (at your option)
  any later version.
 
  ASPECT is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.
 
  You should have received a copy of the GNU General Public License
  along with ASPECT; see the file LICENSE.  If not see
  <http://www.gnu.org/licenses/>.
 */
 
 
#ifndef _aspect_newton_h
#define _aspect_newton_h
 
#include <aspect/simulator/assemblers/interface.h>
#include <aspect/simulator_access.h>
#include <aspect/global.h>
#include <aspect/material_model/interface.h>
 
#include <deal.II/base/signaling_nan.h>
 
namespace aspect
{
  namespace MaterialModel
  {
    template <int dim>
    class MaterialModelDerivatives : public AdditionalMaterialOutputs<dim>
    {
      public:
        MaterialModelDerivatives (const unsigned int n_points);
 
        std::vector<double> viscosity_derivative_wrt_pressure;
 
        std::vector<SymmetricTensor<2,dim>> viscosity_derivative_wrt_strain_rate;
 
        std::vector<double> viscosity_derivative_averaging_weights;
 
        std::vector<double> dilation_derivative_wrt_pressure;
 
        std::vector<SymmetricTensor<2,dim>> dilation_derivative_wrt_strain_rate;
    };
  }
 
 
  namespace Newton
  {
    struct Parameters
    {
 
      enum Stabilization
      {
        none = 0,
        symmetric = 1,
        PD = 2,
        SPD = symmetric | PD
      };
 
 
      friend
      Stabilization
      operator| (const Stabilization a,
                 const Stabilization b)
      {
        return static_cast<Stabilization>(
                 static_cast<int>(a) | static_cast<int>(b));
      }
 
 
      friend
      Stabilization
      operator& (const Stabilization a,
                 const Stabilization b)
      {
        return static_cast<Stabilization>(
                 static_cast<int>(a) & static_cast<int>(b));
      }
 
 
      static void declare_parameters (ParameterHandler &prm);
 
      void parse_parameters (ParameterHandler &prm);
 
      double              newton_derivative_scaling_factor;
 
      Stabilization       preconditioner_stabilization;
      Stabilization       velocity_block_stabilization;
 
      bool                use_Newton_failsafe;
 
      double              nonlinear_switch_tolerance;
 
      bool                use_Eisenstat_Walker_method_for_Picard_iterations;
      unsigned int        max_pre_newton_nonlinear_iterations;
      unsigned int        max_newton_line_search_iterations;
      bool                use_newton_residual_scaling_method;
      double              maximum_linear_stokes_solver_tolerance;
      double              SPD_safety_factor;
    };
 
 
    std::string
    to_string(const Newton::Parameters::Stabilization preconditioner_stabilization);
  }
 
 
 
  template <int dim>
  class NewtonHandler : public SimulatorAccess<dim>
  {
    public:
      void set_assemblers (Assemblers::Manager<dim> &assemblers) const;
 
      static void create_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &output);
 
      Newton::Parameters parameters;
  };
 
 
  namespace Assemblers
  {
    template <int dim>
    class NewtonInterface : public aspect::Assemblers::Interface<dim>,
      public SimulatorAccess<dim>
    {
      public:
        void
        create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
    };
 
    template <int dim>
    class NewtonStokesPreconditioner : public NewtonInterface<dim>
    {
      public:
        void
        execute (internal::Assembly::Scratch::ScratchBase<dim>  &scratch_base,
                 internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
 
        void create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
    };
 
    template <int dim>
    class NewtonStokesCompressiblePreconditioner : public NewtonInterface<dim>
    {
      public:
        void
        execute (internal::Assembly::Scratch::ScratchBase<dim>  &scratch_base,
                 internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
    };
 
    template <int dim>
    class NewtonStokesIncompressibleTerms : public NewtonInterface<dim>
    {
      public:
        void
        execute (internal::Assembly::Scratch::ScratchBase<dim>  &scratch_base,
                 internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
 
        void create_additional_material_model_outputs(MaterialModel::MaterialModelOutputs<dim> &outputs) const override;
    };
 
    template <int dim>
    class NewtonStokesCompressibleStrainRateViscosityTerm : public Assemblers::Interface<dim>,
      public SimulatorAccess<dim>
    {
      public:
        void
        execute(internal::Assembly::Scratch::ScratchBase<dim>  &scratch_base,
                internal::Assembly::CopyData::CopyDataBase<dim> &data_base) const override;
    };
  }
}
 
#endif