parameters.h

Go to the documentation of this file.

// ------------------------------------------------------
//
// SPDX-License-Identifier: GPL-2.0-or-later
//
// SPDX-FileCopyrightText: Copyright (C) 2011-2026 by the ASPECT authors.
//
// This file is part of ASPECT.
//
// Detailed license information governing the source code and
// contributions can be found in the folder LICENSES and in
// LICENSE and CONTRIBUTING.md at the top level directory.
//
// ------------------------------------------------------


#ifndef _aspect_parameters_h
#define _aspect_parameters_h

#include <deal.II/base/parameter_handler.h>

#include <aspect/global.h>
#include <aspect/material_model/interface.h>


namespace aspect
{
  struct CompositionalFieldDescription;

  // forward declaration:
  namespace GeometryModel
  {
    template <int dim>
    class Interface;
  }

  template <int dim>
  struct Parameters
  {

    struct NonlinearSolver
    {
      enum Kind
      {
        no_Advection_no_Stokes,
        no_Advection_single_Stokes,
        no_Advection_single_Stokes_first_timestep_only,
        no_Advection_iterated_Stokes,
        no_Advection_iterated_defect_correction_Stokes,
        single_Advection_no_Stokes,
        single_Advection_single_Stokes,
        single_Advection_iterated_Stokes,
        single_Advection_iterated_defect_correction_Stokes,
        single_Advection_iterated_Newton_Stokes,
        iterated_Advection_no_Stokes,
        iterated_Advection_and_Stokes,
        iterated_Advection_and_defect_correction_Stokes,
        iterated_Advection_and_Newton_Stokes
      };
    };

    static
    bool
    is_defect_correction(const typename NonlinearSolver::Kind &input)
    {
      return input == NonlinearSolver::iterated_Advection_and_Newton_Stokes ||
             input == NonlinearSolver::single_Advection_iterated_Newton_Stokes ||
             input == NonlinearSolver::no_Advection_iterated_defect_correction_Stokes ||
             input == NonlinearSolver::single_Advection_iterated_defect_correction_Stokes ||
             input == NonlinearSolver::iterated_Advection_and_defect_correction_Stokes
             ?
             true
             :
             false;
    }

    static
    bool
    solve_Stokes_iteratively(const typename NonlinearSolver::Kind &input)
    {
      return input == NonlinearSolver::no_Advection_iterated_Stokes ||
             input == NonlinearSolver::no_Advection_iterated_defect_correction_Stokes ||
             input == NonlinearSolver::single_Advection_iterated_Stokes ||
             input == NonlinearSolver::single_Advection_iterated_defect_correction_Stokes ||
             input == NonlinearSolver::single_Advection_iterated_Newton_Stokes ||
             input == NonlinearSolver::iterated_Advection_and_Stokes ||
             input == NonlinearSolver::iterated_Advection_and_defect_correction_Stokes ||
             input == NonlinearSolver::iterated_Advection_and_Newton_Stokes
             ?
             true
             :
             false;
    }

    struct NonlinearSolverFailureStrategy
    {
      enum Kind
      {
        continue_with_next_timestep,
        cut_timestep_size,
        abort_program
      };

      static
      Kind
      parse(const std::string &input)
      {
        if (input == "continue with next timestep")
          return continue_with_next_timestep;
        else if (input == "cut timestep size")
          return cut_timestep_size;
        else if (input == "abort program")
          return abort_program;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }
    };

    struct LinearSolverFailureStrategy
    {
      enum Kind
      {
        continue_with_nonlinear_solver,
        abort
      };

      static
      Kind
      parse(const std::string &input)
      {
        if (input == "continue with nonlinear solver")
          return continue_with_nonlinear_solver;
        else if (input == "abort")
          return abort;
        else
          AssertThrow(false, ExcNotImplemented());
      }
    };

    struct NullspaceRemoval
    {
      enum Kind
      {
        none = 0,
        net_translation_x = 0x1,
        net_translation_y = 0x2,
        net_translation_z = 0x4,
        net_translation   = 0x1+0x2+0x4,
        linear_momentum_x = 0x8,
        linear_momentum_y = 0x10,
        linear_momentum_z = 0x20,
        linear_momentum   = linear_momentum_x+linear_momentum_y+linear_momentum_z,
        net_rotation      = 0x40,
        angular_momentum  = 0x80,
        net_surface_rotation = 0x100,
        any_translation = net_translation+linear_momentum,
        any_rotation = net_rotation+angular_momentum+net_surface_rotation,
      };
    };

    struct AdvectionFieldMethod
    {
      enum Kind
      {
        fem_field,
        particles,
        volume_of_fluid,
        static_field,
        fem_melt_field,
        fem_darcy_field,
        prescribed_field,
        prescribed_field_with_diffusion
      };
    };

    struct Formulation
    {
      enum Kind
      {
        boussinesq_approximation,
        anelastic_liquid_approximation,
        isentropic_compression,
        custom
      };

      static
      Kind
      parse(const std::string &input)
      {
        if (input == "isentropic compression")
          return Formulation::isentropic_compression;
        else if (input == "anelastic liquid approximation")
          return Formulation::anelastic_liquid_approximation;
        else if (input == "Boussinesq approximation")
          return Formulation::boussinesq_approximation;
        else if (input == "custom")
          return Formulation::custom;
        else
          AssertThrow(false, ExcNotImplemented());

        return Formulation::Kind();
      }

      struct MassConservation
      {
        enum Kind
        {
          isentropic_compression,
          hydrostatic_compression,
          reference_density_profile,
          implicit_reference_density_profile,
          incompressible,
          projected_density_field,
          ask_material_model
        };

        static Kind
        parse(const std::string &input)
        {
          if (input == "isentropic compression")
            return Formulation::MassConservation::isentropic_compression;
          else if (input == "hydrostatic compression")
            return Formulation::MassConservation::hydrostatic_compression;
          else if (input == "reference density profile")
            return Formulation::MassConservation::reference_density_profile;
          else if (input == "implicit reference density profile")
            return Formulation::MassConservation::implicit_reference_density_profile;
          else if (input == "incompressible")
            return Formulation::MassConservation::incompressible;
          else if (input == "projected density field")
            return Formulation::MassConservation::projected_density_field;
          else if (input == "ask material model")
            return Formulation::MassConservation::ask_material_model;
          else
            AssertThrow(false, ExcNotImplemented());

          return Formulation::MassConservation::Kind();
        }
      };

      struct TemperatureEquation
      {
        enum Kind
        {
          real_density,
          reference_density_profile
        };

        static
        Kind
        parse(const std::string &input)
        {
          if (input == "real density")
            return Formulation::TemperatureEquation::real_density;
          else if (input == "reference density profile")
            return Formulation::TemperatureEquation::reference_density_profile;
          else
            AssertThrow(false, ExcNotImplemented());

          return Formulation::TemperatureEquation::Kind();
        }
      };
    };

    struct AdvectionStabilizationMethod
    {
      enum Kind
      {
        entropy_viscosity,
        supg
      };

      static
      Kind
      parse(const std::string &input)
      {
        if (input == "entropy viscosity")
          return entropy_viscosity;
        else if (input == "SUPG")
          return supg;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }

      static std::string get_options_string()
      {
        return "entropy viscosity|SUPG";
      }
    };

    struct StokesSolverType
    {
      enum Kind
      {
        block_amg,
        direct_solver,
        block_gmg,
        default_solver
      };

      static const std::string pattern()
      {
        return "default solver|block AMG|direct solver|block GMG";
      }

      static Kind
      parse(const std::string &input)
      {
        if (input == "block AMG")
          return block_amg;
        else if (input == "direct solver")
          return direct_solver;
        else if (input == "block GMG")
          return block_gmg;
        else if (input == "default solver")
          return default_solver;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }
    };

    struct StokesGMGType
    {
      enum Kind
      {
        local_smoothing,
        global_coarsening
      };

      static const std::string pattern()
      {
        return "local smoothing|global coarsening";
      }

      static Kind
      parse(const std::string &input)
      {
        if (input == "local smoothing")
          return local_smoothing;
        else if (input == "global coarsening")
          return global_coarsening;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }
    };

    struct StokesKrylovType
    {
      enum Kind
      {
        gmres,
        idr_s
      };

      static const std::string pattern()
      {
        return "GMRES|IDR(s)";
      }

      static Kind
      parse(const std::string &input)
      {
        if (input == "GMRES")
          return gmres;
        else if (input == "IDR(s)")
          return idr_s;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }
    };

    struct ReactionSolverType
    {
      enum Kind
      {
        ARKode,
        fixed_step
      };

      static const std::string pattern()
      {
        return "ARKode|fixed step";
      }

      static Kind
      parse(const std::string &input)
      {
        if (input == "ARKode")
          return ARKode;
        else if (input == "fixed step")
          return fixed_step;
        else
          AssertThrow(false, ExcNotImplemented());

        return Kind();
      }
    };

    using CompositionalFieldDescription DEAL_II_DEPRECATED = aspect::CompositionalFieldDescription;

    Parameters (ParameterHandler &prm,
                const MPI_Comm mpi_communicator);

    static
    void declare_parameters (ParameterHandler &prm, const unsigned int mpi_rank);

    void parse_parameters (ParameterHandler &prm,
                           const MPI_Comm mpi_communicator);

    void parse_geometry_dependent_parameters (ParameterHandler &prm,
                                              const GeometryModel::Interface<dim> &geometry_model);

    typename NonlinearSolver::Kind nonlinear_solver;
    typename NonlinearSolverFailureStrategy::Kind nonlinear_solver_failure_strategy;
    typename LinearSolverFailureStrategy::Kind linear_solver_failure_strategy;

    typename AdvectionStabilizationMethod::Kind advection_stabilization_method;
    double                         nonlinear_tolerance;
    bool                           resume_computation;
    double                         start_time;
    double                         end_time;
    double                         CFL_number;
    double                         maximum_time_step;
    double                         maximum_relative_increase_time_step;
    double                         maximum_first_time_step;
    bool                           use_artificial_viscosity_smoothing;
    bool                           use_conduction_timestep;
    bool                           convert_to_years;
    std::string                    output_directory;
    double                         surface_pressure;
    double                         adiabatic_surface_temperature;
    unsigned int                   timing_output_frequency;
    unsigned int                   max_nonlinear_iterations;
    unsigned int                   max_nonlinear_iterations_in_prerefinement;
    bool                           use_operator_splitting;
    std::string                    world_builder_file;
    unsigned int                   n_particle_managers;

    double                         temperature_solver_tolerance;
    double                         composition_solver_tolerance;

    // subsection: Advection solver parameters
    unsigned int                   advection_gmres_restart_length;

    // subsection: Stokes solver parameters
    bool                           use_direct_stokes_solver;
    bool                           use_bfbt;
    typename StokesSolverType::Kind stokes_solver_type;
    typename StokesGMGType::Kind stokes_gmg_type;
    typename StokesKrylovType::Kind stokes_krylov_type;
    unsigned int                    idr_s_parameter;

    double                         linear_stokes_solver_tolerance;
    unsigned int                   n_cheap_stokes_solver_steps;
    unsigned int                   n_expensive_stokes_solver_steps;
    double                         linear_solver_A_block_tolerance;
    bool                           use_full_A_block_preconditioner;
    bool                           force_nonsymmetric_A_block_solver;
    double                         linear_solver_S_block_tolerance;
    unsigned int                   stokes_gmres_restart_length;

    // subsection: AMG parameters
    std::string                    AMG_smoother_type;
    unsigned int                   AMG_smoother_sweeps;
    double                         AMG_aggregation_threshold;
    bool                           AMG_output_details;

    // subsection: Operator splitting parameters
    typename ReactionSolverType::Kind reaction_solver_type;
    double                         ARKode_relative_tolerance;
    double                         reaction_time_step;
    unsigned int                   reaction_steps_per_advection_step;

    // subsection: Diffusion solver parameters
    double                         diffusion_length_scale;

    typename Formulation::Kind formulation;

    typename Formulation::MassConservation::Kind formulation_mass_conservation;

    typename Formulation::TemperatureEquation::Kind formulation_temperature_equation;

    bool                           enable_elasticity;

    bool                           include_melt_transport;
    bool                           enable_additional_stokes_rhs;
    bool                           enable_prescribed_dilation;

    std::set<types::boundary_id> fixed_heat_flux_boundary_indicators;

    typename NullspaceRemoval::Kind nullspace_removal;
    unsigned int                   initial_global_refinement;
    unsigned int                   initial_adaptive_refinement;
    double                         refinement_fraction;
    double                         coarsening_fraction;
    bool                           adapt_by_fraction_of_cells;
    unsigned int                   min_grid_level;
    std::vector<double>            additional_refinement_times;
    unsigned int                   adaptive_refinement_interval;
    bool                           skip_solvers_on_initial_refinement;
    bool                           skip_setup_initial_conditions_on_initial_refinement;
    bool                           run_postprocessors_on_initial_refinement;
    bool                           run_postprocessors_on_nonlinear_iterations;
    unsigned int                   stabilization_alpha;
    std::vector<double>            stabilization_c_R;
    std::vector<double>            stabilization_beta;
    double                         stabilization_gamma;
    double                         discontinuous_penalty;
    bool                           use_limiter_for_discontinuous_temperature_solution;
    std::vector<bool>              use_limiter_for_discontinuous_composition_solution;
    double                         global_temperature_max_preset;
    double                         global_temperature_min_preset;
    std::vector<double>            global_composition_max_preset;
    std::vector<double>            global_composition_min_preset;

    std::vector<std::string>       compositional_fields_with_disabled_boundary_entropy_viscosity;
    int                            checkpoint_time_secs;
    int                            checkpoint_steps;
    std::vector<double>            additional_checkpoint_times;
    unsigned int                   n_checkpoints_to_keep;
    unsigned int                   n_additional_checkpoints_to_keep;
    unsigned int                   resume_checkpoint_id;
    double                         resume_time;
    unsigned int                   stokes_velocity_degree;
    bool                           use_locally_conservative_discretization;
    bool                           use_equal_order_interpolation_for_stokes;
    bool                           use_discontinuous_temperature_discretization;
    std::vector<bool>              use_discontinuous_composition_discretization;
    bool                           have_discontinuous_composition_discretization;
    unsigned int                   temperature_degree;
    std::vector<unsigned int>      composition_degrees;
    unsigned int                   max_composition_degree;
    std::string                    pressure_normalization;
    MaterialModel::MaterialAveraging::AveragingOperation material_averaging;

    typename AdvectionFieldMethod::Kind temperature_method;

    unsigned int                   n_compositional_fields;
    std::vector<std::string>       names_of_compositional_fields;
    std::vector<aspect::CompositionalFieldDescription>  composition_descriptions;
    unsigned int                   n_chemical_compositions;
    std::vector<unsigned int>      chemical_composition_indices;
    double                         minimum_composition_fraction;
    bool                           use_pressure_gradient_for_darcy_field;

    std::vector<typename AdvectionFieldMethod::Kind> compositional_field_methods;

    std::map<unsigned int, std::pair<std::string,unsigned int>> mapped_particle_properties;

    std::vector<unsigned int>      normalized_fields;
    bool                           mesh_deformation_enabled;
    bool                           volume_of_fluid_tracking_enabled;
  };

}
#endif