[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
parameters_ode_solver.cpp
1 #include "parameters/parameters_ode_solver.h"
2 
3 namespace PHiLiP {
4 namespace Parameters {
5 
6 void ODESolverParam::declare_parameters (dealii::ParameterHandler &prm)
7 {
8  prm.enter_subsection("ODE solver");
9  {
10 
11  prm.declare_entry("ode_output", "verbose",
12  dealii::Patterns::Selection("quiet|verbose"),
13  "State whether output from ODE solver should be printed. "
14  "Choices are <quiet|verbose>.");
15 
16  prm.declare_entry("output_solution_every_x_steps", "-1",
17  dealii::Patterns::Integer(-1,dealii::Patterns::Integer::max_int_value),
18  "Outputs the solution every x steps in .vtk file");
19 
20  prm.declare_entry("output_solution_every_dt_time_intervals", "0.0",
21  dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
22  "Outputs the solution at time intervals of dt in .vtk file");
23 
24  prm.declare_entry("output_solution_at_fixed_times", "false",
25  dealii::Patterns::Bool(),
26  "Output solution at fixed times. False by default.");
27 
28  prm.declare_entry("output_solution_fixed_times_string", " ",
29  dealii::Patterns::FileName(dealii::Patterns::FileName::FileType::input),
30  "String of the times at which to output the velocity field. "
31  "Example: '0.0 1.0 2.0 3.0 ' or '0.0 1.0 2.0 3.0'");
32 
33  prm.declare_entry("output_solution_at_exact_fixed_times", "false",
34  dealii::Patterns::Bool(),
35  "Output solution at exact fixed times by decreasing the time step on the fly. False by default. "
36  "NOTE: Should be set to false if doing stability studies so that the time step is never influenced by solution file soutput times.");
37 
38  prm.declare_entry("ode_solver_type", "implicit",
39  dealii::Patterns::Selection(
40  " runge_kutta | "
41  " low_storage_runge_kutta | "
42  " implicit | "
43  " rrk_explicit | "
44  " pod_galerkin | "
45  " pod_petrov_galerkin | "
46  " hyper_reduced_petrov_galerkin | "
47  " pod_galerkin_runge_kutta "),
48  "Type of ODE solver to use."
49  "Choices are "
50  " <runge_kutta | "
51  " low_storage_runge_kutta | "
52  " implicit | "
53  " rrk_explicit | "
54  " pod_galerkin | "
55  " pod_petrov_galerkin | "
56  " hyper_reduced_petrov_galerkin | "
57  " pod_galerkin_runge_kutta>.");
58 
59  prm.declare_entry("nonlinear_max_iterations", "500000",
60  dealii::Patterns::Integer(0,dealii::Patterns::Integer::max_int_value),
61  "Maximum nonlinear solver iterations");
62  prm.declare_entry("nonlinear_steady_residual_tolerance", "1e-13",
63  //dealii::Patterns::Double(1e-16,dealii::Patterns::Double::max_double_value),
64  dealii::Patterns::Double(1e-300,dealii::Patterns::Double::max_double_value),
65  "Nonlinear solver residual tolerance");
66  prm.declare_entry("initial_time_step", "100.0",
67  dealii::Patterns::Double(1e-16,dealii::Patterns::Double::max_double_value),
68  "Time step used in ODE solver.");
69  prm.declare_entry("time_step_factor_residual", "0.0",
70  dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
71  "Multiplies initial time-step by time_step_factor_residual*(-log10(residual_norm_decrease)).");
72  prm.declare_entry("time_step_factor_residual_exp", "1.0",
73  dealii::Patterns::Double(0,dealii::Patterns::Double::max_double_value),
74  "Scales initial time step by pow(time_step_factor_residual*(-log10(residual_norm_decrease)),time_step_factor_residual_exp).");
75 
76  prm.declare_entry("print_iteration_modulo", "1",
77  dealii::Patterns::Integer(0,dealii::Patterns::Integer::max_int_value),
78  "Print every print_iteration_modulo iterations of "
79  "the nonlinear solver");
80  prm.declare_entry("output_final_steady_state_solution_to_file", "false",
81  dealii::Patterns::Bool(),
82  "Output final steady state solution to file if set to true");
83  prm.declare_entry("steady_state_final_solution_filename", "solution_snapshot",
84  dealii::Patterns::Anything(),
85  "Filename to use when outputting solution to a file.");
86  prm.declare_entry("output_ode_solver_steady_state_convergence_table","false",
87  dealii::Patterns::Bool(),
88  "Set as false by default. If true, writes the linear solver convergence data "
89  "for steady state to a file named 'ode_solver_steady_state_convergence_data_table.txt'.");
90 
91  prm.declare_entry("initial_time", "0.0",
92  dealii::Patterns::Double(0, dealii::Patterns::Double::max_double_value),
93  "Initial time at which we initialize the ODE solver with.");
94 
95  prm.declare_entry("initial_iteration", "0",
96  dealii::Patterns::Integer(0, dealii::Patterns::Integer::max_int_value),
97  "Initial iteration at which we initialize the ODE solver with.");
98 
99  prm.declare_entry("initial_desired_time_for_output_solution_every_dt_time_intervals", "0.0",
100  dealii::Patterns::Double(0, dealii::Patterns::Double::max_double_value),
101  "Initial desired time for outputting the solution every dt time intervals "
102  "at which we initialize the ODE solver with.");
103 
104  prm.declare_entry("runge_kutta_method", "ssprk3_ex",
105  dealii::Patterns::Selection(
106  " rk4_ex | "
107  " ssprk3_ex | "
108  " heun2_ex | "
109  " euler_ex | "
110  " euler_im | "
111  " dirk_2_im | "
112  " dirk_3_im | "
113  " RK3_2_5F_3SStarPlus | "
114  " RK4_3_5_3SStar | "
115  " RK4_3_9F_3SStarPlus |"
116  " RK5_4_10F_3SStarPlus "),
117  "Runge-kutta method to use. Methods with _ex are explicit, and with _im are implicit. [3S*] and [3S*+] methods are low-storage RK methods"
118  "Choices are "
119  " <rk4_ex | "
120  " ssprk3_ex | "
121  " heun2_ex | "
122  " euler_ex | "
123  " euler_im | "
124  " dirk_2_im | "
125  " dirk_3_im | "
126  " RK4_3_5_3SStar | "
127  " RK3_2_5F_3SStarPlus | "
128  " RK5_4_10F_3SStarPlus |"
129  " RK4_3_9F_3SStarPlus >.");
130  prm.enter_subsection("rrk root solver");
131  {
132  prm.declare_entry("rrk_root_solver_output", "quiet",
133  dealii::Patterns::Selection("quiet|verbose"),
134  "State whether output from rrk root solver should be printed. "
135  "Choices are <quiet|verbose>.");
136 
137  prm.declare_entry("relaxation_runge_kutta_root_tolerance", "5e-10",
138  dealii::Patterns::Double(),
139  "Tolerance for root-finding problem in entropy RRK ode solver."
140  "Defult 5E-10 is suitable in most cases.");
141  }
142  prm.leave_subsection();
143 
144  prm.enter_subsection("low-storage rk solver");
145  {
146  prm.declare_entry("atol", "0.001",
147  dealii::Patterns::Double(),
148  "Absolute Tolerance for automatic step size controller");
149 
150  prm.declare_entry("rtol", "0.001",
151  dealii::Patterns::Double(),
152  "Relative Tolerance for automatic step size controller");
153 
154  prm.declare_entry("beta1", "0.70",
155  dealii::Patterns::Double(),
156  "Beta Controller 1 for automatic step size controller");
157 
158  prm.declare_entry("beta2", "-0.23",
159  dealii::Patterns::Double(),
160  "Beta controller 2 for automatic step size controller");
161 
162  prm.declare_entry("beta3", "0.0",
163  dealii::Patterns::Double(),
164  "Beta controller 3 for automatic step size controller");
165  }
166  prm.leave_subsection();
167  }
168  prm.leave_subsection();
169 }
170 
171 void ODESolverParam::parse_parameters (dealii::ParameterHandler &prm)
172 {
173  prm.enter_subsection("ODE solver");
174  {
175  const std::string output_string = prm.get("ode_output");
176  if (output_string == "quiet") ode_output = OutputEnum::quiet;
177  else if (output_string == "verbose") ode_output = OutputEnum::verbose;
178 
179  output_solution_every_x_steps = prm.get_integer("output_solution_every_x_steps");
180  output_solution_every_dt_time_intervals = prm.get_double("output_solution_every_dt_time_intervals");
181  output_solution_at_fixed_times = prm.get_bool("output_solution_at_fixed_times");
182  output_solution_fixed_times_string = prm.get("output_solution_fixed_times_string");
184  output_solution_at_exact_fixed_times = prm.get_bool("output_solution_at_exact_fixed_times");
185 
186  // Assign ode_solver_type and the allocate AD matrix dRdW flag
187  const std::string solver_string = prm.get("ode_solver_type");
188  if (solver_string == "runge_kutta") { ode_solver_type = ODESolverEnum::runge_kutta_solver;
189  allocate_matrix_dRdW = false; }
190  else if (solver_string == "low_storage_runge_kutta") { ode_solver_type = ODESolverEnum::low_storage_runge_kutta_solver;
191  allocate_matrix_dRdW = false; }
192  else if (solver_string == "implicit") { ode_solver_type = ODESolverEnum::implicit_solver;
193  allocate_matrix_dRdW = true; }
194  else if (solver_string == "rrk_explicit") { ode_solver_type = ODESolverEnum::rrk_explicit_solver;
195  allocate_matrix_dRdW = false; }
196  else if (solver_string == "pod_galerkin") { ode_solver_type = ODESolverEnum::pod_galerkin_solver;
197  allocate_matrix_dRdW = true; }
198  else if (solver_string == "pod_petrov_galerkin") { ode_solver_type = ODESolverEnum::pod_petrov_galerkin_solver;
199  allocate_matrix_dRdW = true; }
200  else if (solver_string == "hyper_reduced_petrov_galerkin") { ode_solver_type = ODESolverEnum::hyper_reduced_petrov_galerkin_solver;
201  allocate_matrix_dRdW = true; }
202  else if (solver_string == "pod_galerkin_runge_kutta") { ode_solver_type = ODESolverEnum::pod_galerkin_runge_kutta_solver;
203  allocate_matrix_dRdW = true; }
204 
205  nonlinear_steady_residual_tolerance = prm.get_double("nonlinear_steady_residual_tolerance");
206  nonlinear_max_iterations = prm.get_integer("nonlinear_max_iterations");
207  initial_time_step = prm.get_double("initial_time_step");
208  time_step_factor_residual = prm.get_double("time_step_factor_residual");
209  time_step_factor_residual_exp = prm.get_double("time_step_factor_residual_exp");
210 
211  print_iteration_modulo = prm.get_integer("print_iteration_modulo");
212  output_final_steady_state_solution_to_file = prm.get_bool("output_final_steady_state_solution_to_file");
213  steady_state_final_solution_filename = prm.get("steady_state_final_solution_filename");
214  output_ode_solver_steady_state_convergence_table = prm.get_bool("output_ode_solver_steady_state_convergence_table");
215 
216  initial_time = prm.get_double("initial_time");
217  initial_iteration = prm.get_integer("initial_iteration");
218  initial_desired_time_for_output_solution_every_dt_time_intervals = prm.get_double("initial_desired_time_for_output_solution_every_dt_time_intervals");
219 
220  const std::string rk_method_string = prm.get("runge_kutta_method");
221  if (rk_method_string == "rk4_ex"){
222  runge_kutta_method = RKMethodEnum::rk4_ex;
223  n_rk_stages = 4;
224  rk_order = 4;
225  }
226  else if (rk_method_string == "ssprk3_ex"){
227  runge_kutta_method = RKMethodEnum::ssprk3_ex;
228  n_rk_stages = 3;
229  rk_order = 3;
230  }
231  else if (rk_method_string == "heun2_ex"){
232  runge_kutta_method = RKMethodEnum::heun2_ex;
233  n_rk_stages = 2;
234  rk_order = 2;
235  }
236  else if (rk_method_string == "euler_ex"){
237  runge_kutta_method = RKMethodEnum::euler_ex;
238  n_rk_stages = 1;
239  rk_order = 1;
240  }
241  else if (rk_method_string == "euler_im"){
242  runge_kutta_method = RKMethodEnum::euler_im;
243  n_rk_stages = 1;
244  rk_order = 1;
245  }
246  else if (rk_method_string == "dirk_2_im"){
247  runge_kutta_method = RKMethodEnum::dirk_2_im;
248  n_rk_stages = 2;
249  rk_order = 2;
250  }
251  else if (rk_method_string == "dirk_3_im"){
252  runge_kutta_method = RKMethodEnum::dirk_3_im;
253  n_rk_stages = 3;
254  rk_order = 3;
255  }
256  else if (rk_method_string == "RK3_2_5F_3SStarPlus"){
257  runge_kutta_method = RKMethodEnum::RK3_2_5F_3SStarPlus;
258  n_rk_stages = 5;
259  num_delta = 5;
260  rk_order = 3;
261  is_3Sstarplus = true;
262  }
263  else if (rk_method_string == "RK4_3_5_3SStar"){
264  runge_kutta_method = RKMethodEnum::RK4_3_5_3SStar;
265  n_rk_stages = 5;
266  num_delta = 7;
267  rk_order = 4;
268  is_3Sstarplus = false;
269  }
270  else if (rk_method_string == "RK4_3_9F_3SStarPlus"){
271  runge_kutta_method = RKMethodEnum::RK4_3_9F_3SStarPlus;
272  n_rk_stages = 9;
273  num_delta = 9;
274  rk_order = 4;
275  is_3Sstarplus = true;
276  }
277  else if (rk_method_string == "RK5_4_10F_3SStarPlus"){
278  runge_kutta_method = RKMethodEnum::RK5_4_10F_3SStarPlus;
279  n_rk_stages = 10;
280  num_delta = 10;
281  rk_order = 5;
282  is_3Sstarplus = true;
283  }
284 
285  prm.enter_subsection("rrk root solver");
286  {
287  const std::string output_string_rrk = prm.get("rrk_root_solver_output");
288  if (output_string_rrk == "verbose") rrk_root_solver_output = verbose;
289  else if (output_string_rrk == "quiet") rrk_root_solver_output = quiet;
290 
291  relaxation_runge_kutta_root_tolerance = prm.get_double("relaxation_runge_kutta_root_tolerance");
292  }
293  prm.leave_subsection();
294 
295  prm.enter_subsection("low-storage rk solver");
296  {
297  atol = prm.get_double("atol");
298  rtol = prm.get_double("rtol");
299  beta1 = prm.get_double("beta1");
300  beta2 = prm.get_double("beta2");
301  beta3 = prm.get_double("beta3");
302  }
303  prm.leave_subsection();
304 
305  }
306  prm.leave_subsection();
307 }
308 
309 } // Parameters namespace
310 } // PHiLiP namespace
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.
Definition: ADTypes.hpp:10
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.
unsigned int initial_iteration
Initial iteration at which we initialize the ODE solver with.
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.
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.
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.