1 #include "parameters/parameters_ode_solver.h"     2 #include <deal.II/base/mpi.h>     3 #include <deal.II/base/utilities.h>     4 #include <deal.II/base/conditional_ostream.h>    11     prm.enter_subsection(
"ODE solver");
    14         prm.declare_entry(
"ode_output", 
"verbose",
    15                           dealii::Patterns::Selection(
"quiet|verbose"),
    16                           "State whether output from ODE solver should be printed. "    17                           "Choices are <quiet|verbose>.");
    19         prm.declare_entry(
"output_solution_every_x_steps", 
"-1",
    20                           dealii::Patterns::Integer(-1,dealii::Patterns::Integer::max_int_value),
    21                           "Outputs the solution every x steps in .vtk file");
    23         prm.declare_entry(
"output_solution_every_dt_time_intervals", 
"0.0",
    24                           dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
    25                           "Outputs the solution at time intervals of dt in .vtk file");
    27         prm.declare_entry(
"output_solution_at_fixed_times", 
"false",
    28                           dealii::Patterns::Bool(),
    29                           "Output solution at fixed times. False by default.");
    31         prm.declare_entry(
"output_solution_fixed_times_string", 
" ",
    32                           dealii::Patterns::FileName(dealii::Patterns::FileName::FileType::input),
    33                           "String of the times at which to output the velocity field. "    34                           "Example: '0.0 1.0 2.0 3.0 ' or '0.0 1.0 2.0 3.0'");
    36         prm.declare_entry(
"output_solution_at_exact_fixed_times", 
"false",
    37                           dealii::Patterns::Bool(),
    38                           "Output solution at exact fixed times by decreasing the time step on the fly. False by default. "    39                           "NOTE: Should be set to false if doing stability studies so that the time step is never influenced by solution file soutput times.");
    41         prm.declare_entry(
"ode_solver_type", 
"implicit",
    42                           dealii::Patterns::Selection(
    44                           " low_storage_runge_kutta | "    48                           " pod_petrov_galerkin | "    49                           " hyper_reduced_petrov_galerkin | "    50                           " pod_galerkin_runge_kutta "),
    51                           "Type of ODE solver to use."    54                           " low_storage_runge_kutta | "    58                           " pod_petrov_galerkin | "    59                           " hyper_reduced_petrov_galerkin | "    60                           " pod_galerkin_runge_kutta>.");
    62         prm.declare_entry(
"nonlinear_max_iterations", 
"500000",
    63                           dealii::Patterns::Integer(0,dealii::Patterns::Integer::max_int_value),
    64                           "Maximum nonlinear solver iterations");
    65         prm.declare_entry(
"nonlinear_steady_residual_tolerance", 
"1e-13",
    67                           dealii::Patterns::Double(1e-300,dealii::Patterns::Double::max_double_value),
    68                           "Nonlinear solver residual tolerance");
    69         prm.declare_entry(
"initial_time_step", 
"100.0",
    70                           dealii::Patterns::Double(1e-16,dealii::Patterns::Double::max_double_value),
    71                           "Time step used in ODE solver.");
    72         prm.declare_entry(
"time_step_factor_residual", 
"0.0",
    73                           dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
    74                           "Multiplies initial time-step by time_step_factor_residual*(-log10(residual_norm_decrease)).");
    75         prm.declare_entry(
"time_step_factor_residual_exp", 
"1.0",
    76                           dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
    77                           "Scales initial time step by pow(time_step_factor_residual*(-log10(residual_norm_decrease)),time_step_factor_residual_exp).");
    79         prm.declare_entry(
"print_iteration_modulo", 
"1",
    80                           dealii::Patterns::Integer(0,dealii::Patterns::Integer::max_int_value),
    81                           "Print every print_iteration_modulo iterations of "    82                           "the nonlinear solver");
    83         prm.declare_entry(
"output_final_steady_state_solution_to_file", 
"false",
    84                           dealii::Patterns::Bool(),
    85                           "Output final steady state solution to file if set to true");
    86         prm.declare_entry(
"steady_state_final_solution_filename", 
"solution_snapshot",
    87                           dealii::Patterns::Anything(),
    88                           "Filename to use when outputting solution to a file.");
    89         prm.declare_entry(
"output_ode_solver_steady_state_convergence_table",
"false",
    90                           dealii::Patterns::Bool(),
    91                           "Set as false by default. If true, writes the linear solver convergence data "    92                           "for steady state to a file named 'ode_solver_steady_state_convergence_data_table.txt'.");
    94         prm.declare_entry(
"initial_time", 
"0.0",
    95                           dealii::Patterns::Double(0, dealii::Patterns::Double::max_double_value),
    96                           "Initial time at which we initialize the ODE solver with.");
    98         prm.declare_entry(
"initial_iteration", 
"0",
    99                           dealii::Patterns::Integer(0, dealii::Patterns::Integer::max_int_value),
   100                           "Initial iteration at which we initialize the ODE solver with.");
   102         prm.declare_entry(
"initial_desired_time_for_output_solution_every_dt_time_intervals", 
"0.0",
   103                            dealii::Patterns::Double(0, dealii::Patterns::Double::max_double_value),
   104                            "Initial desired time for outputting the solution every dt time intervals "   105                            "at which we initialize the ODE solver with.");
   107         prm.declare_entry(
"runge_kutta_method", 
"ssprk3_ex",
   108                           dealii::Patterns::Selection(
   116                           " RK3_2_5F_3SStarPlus | "   118                           " RK4_3_9F_3SStarPlus |"   119                           " RK5_4_10F_3SStarPlus "),
   120                           "Runge-kutta method to use. Methods with _ex are explicit, and with _im are implicit. [3S*] and [3S*+] methods are low-storage RK methods"   130                           " RK3_2_5F_3SStarPlus | "   131                           " RK5_4_10F_3SStarPlus |"   132                           " RK4_3_9F_3SStarPlus >.");
   133         prm.enter_subsection(
"rrk root solver");
   135             prm.declare_entry(
"rrk_root_solver_output", 
"quiet",
   136                               dealii::Patterns::Selection(
"quiet|verbose"),
   137                               "State whether output from rrk root solver should be printed. "   138                               "Choices are <quiet|verbose>.");
   140             prm.declare_entry(
"relaxation_runge_kutta_root_tolerance", 
"5e-10",
   141                               dealii::Patterns::Double(),
   142                               "Tolerance for root-finding problem in entropy RRK ode solver."   143                               "Defult 5E-10 is suitable in most cases.");
   144             prm.declare_entry(
"use_relaxation_runge_kutta",
"false",
   145                               dealii::Patterns::Bool(),
   146                               "Toggle using relaxation runge-kutta. "   147                               "Must use a RK ode solver."   151         prm.leave_subsection();
   153         prm.enter_subsection(
"low-storage rk solver");
   155             prm.declare_entry(
"atol", 
"0.001",
   156                           dealii::Patterns::Double(),
   157                           "Absolute Tolerance for automatic step size controller");
   159             prm.declare_entry(
"rtol", 
"0.001",
   160                           dealii::Patterns::Double(),
   161                           "Relative Tolerance for automatic step size controller");
   163             prm.declare_entry(
"beta1", 
"0.70",
   164                           dealii::Patterns::Double(),
   165                           "Beta Controller 1 for automatic step size controller");
   167             prm.declare_entry(
"beta2", 
"-0.23",
   168                           dealii::Patterns::Double(),
   169                           "Beta controller 2 for automatic step size controller");
   171             prm.declare_entry(
"beta3", 
"0.0",
   172                            dealii::Patterns::Double(),
   173                           "Beta controller 3 for automatic step size controller");
   175         prm.leave_subsection();
   177     prm.leave_subsection();
   182     prm.enter_subsection(
"ODE solver");
   184         const std::string output_string = prm.get(
"ode_output");
   185         if (output_string == 
"quiet")   
ode_output = OutputEnum::quiet;
   186         else if (output_string == 
"verbose") 
ode_output = OutputEnum::verbose;
   196         const std::string solver_string = prm.get(
"ode_solver_type");
   197         if (solver_string == 
"runge_kutta")              { 
ode_solver_type = ODESolverEnum::runge_kutta_solver;
   199         else if (solver_string == 
"low_storage_runge_kutta")  { 
ode_solver_type = ODESolverEnum::low_storage_runge_kutta_solver;
   201         else if (solver_string == 
"implicit")            { 
ode_solver_type = ODESolverEnum::implicit_solver;
   203         else if (solver_string == 
"rrk_explicit")        { 
ode_solver_type = ODESolverEnum::rrk_explicit_solver;
   205         else if (solver_string == 
"pod_galerkin")        { 
ode_solver_type = ODESolverEnum::pod_galerkin_solver;
   207         else if (solver_string == 
"pod_petrov_galerkin") { 
ode_solver_type = ODESolverEnum::pod_petrov_galerkin_solver;
   209         else if (solver_string == 
"hyper_reduced_petrov_galerkin") { 
ode_solver_type = ODESolverEnum::hyper_reduced_petrov_galerkin_solver;
   211         else if (solver_string == 
"pod_galerkin_runge_kutta") { 
ode_solver_type = ODESolverEnum::pod_galerkin_runge_kutta_solver;
   229         const std::string rk_method_string = prm.get(
"runge_kutta_method");
   230         if (rk_method_string == 
"rk4_ex"){
   235         else if (rk_method_string == 
"ssprk3_ex"){
   240         else if (rk_method_string == 
"heun2_ex"){
   245         else if (rk_method_string == 
"euler_ex"){
   250         else if (rk_method_string == 
"euler_im"){
   255         else if (rk_method_string == 
"dirk_2_im"){
   260         else if (rk_method_string == 
"dirk_3_im"){
   265         else if (rk_method_string == 
"RK3_2_5F_3SStarPlus"){
   272         else if (rk_method_string == 
"RK4_3_5_3SStar"){
   279         else if (rk_method_string == 
"RK4_3_9F_3SStarPlus"){
   286         else if (rk_method_string == 
"RK5_4_10F_3SStarPlus"){
   294         prm.enter_subsection(
"rrk root solver");
   296             const std::string output_string_rrk = prm.get(
"rrk_root_solver_output");
   305                 const int mpi_rank = dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
   306                 dealii::ConditionalOStream pcout(std::cout, mpi_rank==0);
   307                 pcout << 
"Warning: rrk_explicit_solver parameter is depreciated. " << std::endl
   308                       << 
"Backwards compatibility was verified upon implementation." <<std::endl;
   312         prm.leave_subsection();
   314         prm.enter_subsection(
"low-storage rk solver");
   316             atol = prm.get_double(
"atol");
   317             rtol = prm.get_double(
"rtol");
   318             beta1 = prm.get_double(
"beta1");
   319             beta2 = prm.get_double(
"beta2");
   320             beta3 = prm.get_double(
"beta3");
   322         prm.leave_subsection();
   325     prm.leave_subsection();
 bool is_3Sstarplus
True or false depending on what low-storage RK method is used. 
OutputEnum ode_output
verbose or quiet. 
double initial_time
Initial time at which we initialize the ODE solver with. 
double beta2
Second value for beta controller;. 
int output_solution_every_x_steps
Outputs the solution every x steps to .vtk file. 
int num_delta
Number of delta values in low-storage RK methods. 
bool output_solution_at_fixed_times
Flag for outputting solution at fixed times. 
double time_step_factor_residual_exp
Scales initial time step by pow(time_step_factor_residual*(-log10(residual_norm_decrease)),time_step_factor_residual_exp) 
Files for the baseline physics. 
OutputEnum rrk_root_solver_output
Do output for root solving routine. 
unsigned int print_iteration_modulo
If ode_output==verbose, print every print_iteration_modulo iterations. 
bool output_ode_solver_steady_state_convergence_table
unsigned int initial_iteration
Initial iteration at which we initialize the ODE solver with. 
double atol
Absolute tolerance. 
double output_solution_every_dt_time_intervals
Outputs the solution every dt time intervals to .vtk file. 
bool allocate_matrix_dRdW
Flag to signal that automatic differentiation (AD) matrix dRdW must be allocated. ...
int rk_order
Order of the RK method; assigned based on runge_kutta_method. 
std::string output_solution_fixed_times_string
String of fixed solution output times. 
double initial_time_step
Time step used in ODE solver. 
double nonlinear_steady_residual_tolerance
Tolerance to determine steady-state convergence. 
double rtol
Relative tolerance. 
double initial_desired_time_for_output_solution_every_dt_time_intervals
std::string steady_state_final_solution_filename
Filename to write final steady state solution. 
unsigned int number_of_fixed_times_to_output_solution
Number of fixed times to output the solution. 
bool output_final_steady_state_solution_to_file
Output final steady state solution to file. 
RKMethodEnum runge_kutta_method
Runge-kutta method. 
double time_step_factor_residual
Multiplies initial time-step by time_step_factor_residual*(-log10(residual_norm_decrease)) ...
void parse_parameters(dealii::ParameterHandler &prm)
Parses input file and sets the variables. 
bool output_solution_at_exact_fixed_times
Flag for outputting the solution at exact fixed times by decreasing the time step on the fly...
ODESolverEnum ode_solver_type
ODE solver type. 
double relaxation_runge_kutta_root_tolerance
Tolerance for RRK root solver, default value 5E-10. 
double beta3
Third value for beta controller;. 
int n_rk_stages
Number of stages for an RK method; assigned based on runge_kutta_method. 
bool use_relaxation_runge_kutta
Use relaxation runge-kutta. 
double beta1
First value for beta controller;. 
unsigned int nonlinear_max_iterations
Maximum number of iterations. 
static void declare_parameters(dealii::ParameterHandler &prm)
Declares the possible variables and sets the defaults.