libMesh
femparameters.C
Go to the documentation of this file.
1 
2 #include "femparameters.h"
3 
4 #include "libmesh/parsed_function.h"
5 #include "libmesh/zero_function.h"
6 
7 using namespace libMesh;
8 
9 #define GETPOT_INPUT(A) { A = input(#A, A); \
10  const std::string stringval = input(#A, std::string()); \
11  variable_names.push_back(std::string(#A "=") + stringval); }
12 #define GETPOT_INT_INPUT(A) { A = input(#A, (int)A); \
13  const std::string stringval = input(#A, std::string()); \
14  variable_names.push_back(std::string(#A "=") + stringval); }
15 
16 #define GETPOT_FUNCTION_INPUT(A) { \
17  const std::string type = input(#A "_type", "zero"); \
18  const std::string value = input(#A "_value", ""); \
19  A = new_function_base(type, value); }
20 #define GETPOT_REGISTER(A) { \
21  std::string stringval = input(#A, std::string()); \
22  variable_names.push_back(std::string(#A "=") + stringval); }
23 
25  ParallelObject(comm_in),
26  initial_timestep(0), n_timesteps(100),
27  transient(true),
28  deltat_reductions(0),
29  timesolver_core("euler"),
30  end_time(std::numeric_limits<Real>::max()),
31  deltat(0.0001), timesolver_theta(0.5),
32  timesolver_maxgrowth(0.), timesolver_tolerance(0.),
33  timesolver_upper_tolerance(0.),
34  steadystate_tolerance(0.),
35  timesolver_norm(0,L2),
36 
37  dimension(2),
38  domaintype("square"), domainfile("mesh.xda"), elementtype("quad"),
39  elementorder(2),
40  domain_xmin(0.0), domain_ymin(0.0), domain_zmin(0.0),
41  domain_edge_width(1.0), domain_edge_length(1.0), domain_edge_height(1.0),
42  coarsegridx(1), coarsegridy(1), coarsegridz(1),
43  coarserefinements(0), extrarefinements(0),
44  mesh_redistribute_func("0"),
45 
46  nelem_target(8000), global_tolerance(0.0),
47  refine_fraction(0.3), coarsen_fraction(0.3), coarsen_threshold(10),
48  max_adaptivesteps(1),
49  initial_adaptivesteps(0),
50 
51  write_interval(10),
52  write_gmv_error(false), write_tecplot_error(false),
53  write_exodus_error(false),
54  output_xda(false), output_xdr(false),
55  output_bz2(true), output_gz(true),
56  output_gmv(false), output_tecplot(false),
57  output_exodus(false), output_nemesis(false),
58 
59  system_types(0),
60 
61  periodic_boundaries(0),
62 
63  run_simulation(true), run_postprocess(false),
64 
65  fe_family(1, "LAGRANGE"), fe_order(1, 1),
66  extra_quadrature_order(0),
67 
68  analytic_jacobians(true), verify_analytic_jacobians(0.0),
69  numerical_jacobian_h(TOLERANCE),
70  print_solution_norms(false), print_solutions(false),
71  print_residual_norms(false), print_residuals(false),
72  print_jacobian_norms(false), print_jacobians(false),
73  print_element_solutions(false),
74  print_element_residuals(false),
75  print_element_jacobians(false),
76 
77  use_petsc_snes(false),
78  time_solver_quiet(true), solver_quiet(true), solver_verbose(false),
79  reuse_preconditioner(true),
80  require_residual_reduction(true),
81  min_step_length(1e-5),
82  max_linear_iterations(200000), max_nonlinear_iterations(20),
83  relative_step_tolerance(1.e-7), relative_residual_tolerance(1.e-10),
84  absolute_residual_tolerance(1.e-10),
85  initial_linear_tolerance(1.e-3), minimum_linear_tolerance(TOLERANCE*TOLERANCE),
86  linear_tolerance_multiplier(1.e-3),
87 
88  initial_sobolev_order(1),
89  initial_extra_quadrature(0),
90  refine_uniformly(false),
91  indicator_type("kelly"), patch_reuse(true), sobolev_order(1)
92 {
93 }
94 
96 {
97  for (std::map<subdomain_id_type, FunctionBase<Number> *>::iterator
98  i = initial_conditions.begin(); i != initial_conditions.end();
99  ++i)
100  delete i->second;
101 
102  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
103  i = dirichlet_conditions.begin(); i != dirichlet_conditions.end();
104  ++i)
105  delete i->second;
106 
107  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
108  i = neumann_conditions.begin(); i != neumann_conditions.end();
109  ++i)
110  delete i->second;
111 
112  for (std::map<int,
114  >::iterator
115  i = other_boundary_functions.begin(); i != other_boundary_functions.end();
116  ++i)
117  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
118  j = i->second.begin(); j != i->second.end();
119  ++j)
120  delete j->second;
121 
122  for (std::map<int,
124  >::iterator
125  i = other_interior_functions.begin(); i != other_interior_functions.end();
126  ++i)
127  for (std::map<subdomain_id_type, FunctionBase<Number> *>::iterator
128  j = i->second.begin(); j != i->second.end();
129  ++j)
130  delete j->second;
131 }
132 
133 
134 AutoPtr<FunctionBase<Number>> new_function_base(const std::string& func_type,
135  const std::string& func_value)
136 {
137  if (func_type == "parsed")
139  (new ParsedFunction<Number>(func_value));
140  else if (func_type == "zero")
142  (new ZeroFunction<Number>);
143  else
144  libmesh_not_implemented();
145 
147 }
148 
149 
150 void FEMParameters::read(GetPot & input,
151  const std::vector<std::string> * other_variable_names)
152 {
153  std::vector<std::string> variable_names;
154  if (other_variable_names)
155  for (unsigned int i=0; i != other_variable_names->size(); ++i)
156  {
157  const std::string& name = (*other_variable_names)[i];
158  const std::string stringval = input(name, std::string());
159  variable_names.push_back(name + "=" + stringval);
160  }
161 
162  GETPOT_INT_INPUT(initial_timestep);
163  GETPOT_INT_INPUT(n_timesteps);
164  GETPOT_INPUT(transient);
165  GETPOT_INT_INPUT(deltat_reductions);
166  GETPOT_INPUT(timesolver_core);
167  GETPOT_INPUT(end_time);
168  GETPOT_INPUT(deltat);
169  GETPOT_INPUT(timesolver_theta);
170  GETPOT_INPUT(timesolver_maxgrowth);
171  GETPOT_INPUT(timesolver_tolerance);
172  GETPOT_INPUT(timesolver_upper_tolerance);
173  GETPOT_INPUT(steadystate_tolerance);
174 
175  GETPOT_REGISTER(timesolver_norm);
176  const unsigned int n_timesolver_norm = input.vector_variable_size("timesolver_norm");
177  timesolver_norm.resize(n_timesolver_norm, L2);
178  for (unsigned int i=0; i != n_timesolver_norm; ++i)
179  {
180  int current_norm = 0; // L2
181  if (timesolver_norm[i] == H1)
182  current_norm = 1;
183  if (timesolver_norm[i] == H2)
184  current_norm = 2;
185  current_norm = input("timesolver_norm", current_norm, i);
186  if (current_norm == 0)
187  timesolver_norm[i] = L2;
188  else if (current_norm == 1)
189  timesolver_norm[i] = H1;
190  else if (current_norm == 2)
191  timesolver_norm[i] = H2;
192  else
194  }
195 
196 
197  GETPOT_INT_INPUT(dimension);
198  GETPOT_INPUT(domaintype);
199  GETPOT_INPUT(domainfile);
200  GETPOT_INPUT(elementtype);
201  GETPOT_INPUT(elementorder);
202  GETPOT_INPUT(domain_xmin);
203  GETPOT_INPUT(domain_ymin);
204  GETPOT_INPUT(domain_zmin);
205  GETPOT_INPUT(domain_edge_width);
206  GETPOT_INPUT(domain_edge_length);
207  GETPOT_INPUT(domain_edge_height);
208  GETPOT_INT_INPUT(coarsegridx);
209  GETPOT_INT_INPUT(coarsegridy);
210  GETPOT_INT_INPUT(coarsegridz);
211  GETPOT_INT_INPUT(coarserefinements);
212  GETPOT_INT_INPUT(extrarefinements);
213  GETPOT_INPUT(mesh_redistribute_func);
214 
215 
216  GETPOT_INT_INPUT(nelem_target);
217  GETPOT_INPUT(global_tolerance);
218  GETPOT_INPUT(refine_fraction);
219  GETPOT_INPUT(coarsen_fraction);
220  GETPOT_INPUT(coarsen_threshold);
221  GETPOT_INT_INPUT(max_adaptivesteps);
222  GETPOT_INT_INPUT(initial_adaptivesteps);
223 
224 
225  GETPOT_INT_INPUT(write_interval);
226  GETPOT_INPUT(output_xda);
227  GETPOT_INPUT(output_xdr);
228  GETPOT_INPUT(output_gz);
229 #ifndef LIBMESH_HAVE_GZSTREAM
230  output_gz = false;
231 #endif
232  GETPOT_INPUT(output_bz2);
233 #ifndef LIBMESH_HAVE_BZ2
234  output_bz2 = false;
235 #endif
236  GETPOT_INPUT(output_gmv);
237  GETPOT_INPUT(write_gmv_error);
238 #ifndef LIBMESH_HAVE_GMV
239  output_gmv = false;
240  write_gmv_error = false;
241 #endif
242  GETPOT_INPUT(output_tecplot);
243  GETPOT_INPUT(write_tecplot_error);
244 #ifndef LIBMESH_HAVE_TECPLOT_API
245  output_tecplot = false;
246  write_tecplot_error = false;
247 #endif
248  GETPOT_INPUT(output_exodus);
249  GETPOT_INPUT(write_exodus_error);
250 #ifndef LIBMESH_HAVE_EXODUS_API
251  output_exodus = false;
252  write_exodus_error = false;
253 #endif
254  GETPOT_INPUT(output_nemesis);
255 #ifndef LIBMESH_HAVE_NEMESIS_API
256  output_nemesis = false;
257 #endif
258 
259 
260  GETPOT_REGISTER(system_types);
261  const unsigned int n_system_types =
262  input.vector_variable_size("system_types");
263  if (n_system_types)
264  {
265  system_types.resize(n_system_types, "");
266  for (unsigned int i=0; i != n_system_types; ++i)
267  {
268  system_types[i] = input("system_types", system_types[i], i);
269  }
270  }
271 
272 
273  GETPOT_REGISTER(periodic_boundaries);
274  const unsigned int n_periodic_bcs =
275  input.vector_variable_size("periodic_boundaries");
276 
277  if (n_periodic_bcs)
278  {
279  if (domaintype != "square" &&
280  domaintype != "cylinder" &&
281  domaintype != "file" &&
282  domaintype != "od2")
283  {
284  libMesh::out << "Periodic boundaries need rectilinear domains" << std::endl;;
285  libmesh_error();
286  }
287  for (unsigned int i=0; i != n_periodic_bcs; ++i)
288  {
289  unsigned int myboundary =
290  input("periodic_boundaries", -1, i);
291  unsigned int pairedboundary = 0;
292  RealVectorValue translation_vector;
293  if (dimension == 2)
294  switch (myboundary)
295  {
296  case 0:
297  pairedboundary = 2;
298  translation_vector = RealVectorValue(0., domain_edge_length);
299  break;
300  case 1:
301  pairedboundary = 3;
302  translation_vector = RealVectorValue(-domain_edge_width, 0);
303  break;
304  default:
305  libMesh::out << "Unrecognized periodic boundary id " <<
306  myboundary << std::endl;;
307  libmesh_error();
308  }
309  else if (dimension == 3)
310  switch (myboundary)
311  {
312  case 0:
313  pairedboundary = 5;
314  translation_vector = RealVectorValue(Real(0), Real(0), domain_edge_height);
315  break;
316  case 1:
317  pairedboundary = 3;
318  translation_vector = RealVectorValue(Real(0), domain_edge_length, Real(0));
319  break;
320  case 2:
321  pairedboundary = 4;
322  translation_vector = RealVectorValue(-domain_edge_width, Real(0), Real(0));
323  break;
324  default:
325  libMesh::out << "Unrecognized periodic boundary id " <<
326  myboundary << std::endl;;
327  libmesh_error();
328  }
329  periodic_boundaries.push_back(PeriodicBoundary(translation_vector));
330  periodic_boundaries[i].myboundary = myboundary;
331  periodic_boundaries[i].pairedboundary = pairedboundary;
332  }
333  }
334 
335  // Use std::string inputs so GetPot doesn't have to make a bunch
336  // of internal C string copies
337  std::string zero_string = "zero";
338  std::string empty_string = "";
339 
340  GETPOT_REGISTER(dirichlet_condition_types);
341  GETPOT_REGISTER(dirichlet_condition_values);
342  GETPOT_REGISTER(dirichlet_condition_boundaries);
343  GETPOT_REGISTER(dirichlet_condition_variables);
344 
345  const unsigned int n_dirichlet_conditions=
346  input.vector_variable_size("dirichlet_condition_types");
347 
348  if (n_dirichlet_conditions !=
349  input.vector_variable_size("dirichlet_condition_values"))
350  {
351  libMesh::out << "Error: " << n_dirichlet_conditions
352  << " Dirichlet condition types does not match "
353  << input.vector_variable_size("dirichlet_condition_values")
354  << " Dirichlet condition values." << std::endl;
355 
356  libmesh_error();
357  }
358 
359  if (n_dirichlet_conditions !=
360  input.vector_variable_size("dirichlet_condition_boundaries"))
361  {
362  libMesh::out << "Error: " << n_dirichlet_conditions
363  << " Dirichlet condition types does not match "
364  << input.vector_variable_size("dirichlet_condition_boundaries")
365  << " Dirichlet condition boundaries." << std::endl;
366 
367  libmesh_error();
368  }
369 
370  if (n_dirichlet_conditions !=
371  input.vector_variable_size("dirichlet_condition_variables"))
372  {
373  libMesh::out << "Error: " << n_dirichlet_conditions
374  << " Dirichlet condition types does not match "
375  << input.vector_variable_size("dirichlet_condition_variables")
376  << " Dirichlet condition variables sets." << std::endl;
377 
378  libmesh_error();
379  }
380 
381  for (unsigned int i=0; i != n_dirichlet_conditions; ++i)
382  {
383  const std::string func_type =
384  input("dirichlet_condition_types", zero_string, i);
385 
386  const std::string func_value =
387  input("dirichlet_condition_values", empty_string, i);
388 
389  const boundary_id_type func_boundary =
390  input("dirichlet_condition_boundaries", boundary_id_type(0), i);
391 
392  dirichlet_conditions[func_boundary] =
393  (new_function_base(func_type, func_value).release());
394 
395  const std::string variable_set =
396  input("dirichlet_condition_variables", empty_string, i);
397 
398  for (unsigned int i=0; i != variable_set.size(); ++i)
399  {
400  if (variable_set[i] == '1')
401  dirichlet_condition_variables[func_boundary].push_back(i);
402  else if (variable_set[i] != '0')
403  {
404  libMesh::out << "Unable to understand Dirichlet variable set"
405  << variable_set << std::endl;
406  libmesh_error();
407  }
408  }
409  }
410 
411  GETPOT_REGISTER(neumann_condition_types);
412  GETPOT_REGISTER(neumann_condition_values);
413  GETPOT_REGISTER(neumann_condition_boundaries);
414  GETPOT_REGISTER(neumann_condition_variables);
415 
416  const unsigned int n_neumann_conditions=
417  input.vector_variable_size("neumann_condition_types");
418 
419  if (n_neumann_conditions !=
420  input.vector_variable_size("neumann_condition_values"))
421  {
422  libMesh::out << "Error: " << n_neumann_conditions
423  << " Neumann condition types does not match "
424  << input.vector_variable_size("neumann_condition_values")
425  << " Neumann condition values." << std::endl;
426 
427  libmesh_error();
428  }
429 
430  if (n_neumann_conditions !=
431  input.vector_variable_size("neumann_condition_boundaries"))
432  {
433  libMesh::out << "Error: " << n_neumann_conditions
434  << " Neumann condition types does not match "
435  << input.vector_variable_size("neumann_condition_boundaries")
436  << " Neumann condition boundaries." << std::endl;
437 
438  libmesh_error();
439  }
440 
441  if (n_neumann_conditions !=
442  input.vector_variable_size("neumann_condition_variables"))
443  {
444  libMesh::out << "Error: " << n_neumann_conditions
445  << " Neumann condition types does not match "
446  << input.vector_variable_size("neumann_condition_variables")
447  << " Neumann condition variables sets." << std::endl;
448 
449  libmesh_error();
450  }
451 
452  for (unsigned int i=0; i != n_neumann_conditions; ++i)
453  {
454  const std::string func_type =
455  input("neumann_condition_types", zero_string, i);
456 
457  const std::string func_value =
458  input("neumann_condition_values", empty_string, i);
459 
460  const boundary_id_type func_boundary =
461  input("neumann_condition_boundaries", boundary_id_type(0), i);
462 
463  neumann_conditions[func_boundary] =
464  (new_function_base(func_type, func_value).release());
465 
466  const std::string variable_set =
467  input("neumann_condition_variables", empty_string, i);
468 
469  for (unsigned int i=0; i != variable_set.size(); ++i)
470  {
471  if (variable_set[i] == '1')
472  neumann_condition_variables[func_boundary].push_back(i);
473  else if (variable_set[i] != '0')
474  {
475  libMesh::out << "Unable to understand Neumann variable set"
476  << variable_set << std::endl;
477  libmesh_error();
478  }
479  }
480  }
481 
482  GETPOT_REGISTER(initial_condition_types);
483  GETPOT_REGISTER(initial_condition_values);
484  GETPOT_REGISTER(initial_condition_subdomains);
485 
486  const unsigned int n_initial_conditions=
487  input.vector_variable_size("initial_condition_types");
488 
489  if (n_initial_conditions !=
490  input.vector_variable_size("initial_condition_values"))
491  {
492  libMesh::out << "Error: " << n_initial_conditions
493  << " initial condition types does not match "
494  << input.vector_variable_size("initial_condition_values")
495  << " initial condition values." << std::endl;
496 
497  libmesh_error();
498  }
499 
500  if (n_initial_conditions !=
501  input.vector_variable_size("initial_condition_subdomains"))
502  {
503  libMesh::out << "Error: " << n_initial_conditions
504  << " initial condition types does not match "
505  << input.vector_variable_size("initial_condition_subdomains")
506  << " initial condition subdomains." << std::endl;
507 
508  libmesh_error();
509  }
510 
511  for (unsigned int i=0; i != n_initial_conditions; ++i)
512  {
513  const std::string func_type =
514  input("initial_condition_types", zero_string, i);
515 
516  const std::string func_value =
517  input("initial_condition_values", empty_string, i);
518 
519  const subdomain_id_type func_subdomain =
520  input("initial_condition_subdomains", subdomain_id_type(0), i);
521 
522  initial_conditions[func_subdomain] =
523  (new_function_base(func_type, func_value).release());
524  }
525 
526  GETPOT_REGISTER(other_interior_function_types);
527  GETPOT_REGISTER(other_interior_function_values);
528  GETPOT_REGISTER(other_interior_function_subdomains);
529  GETPOT_REGISTER(other_interior_function_ids);
530 
531  const unsigned int n_other_interior_functions =
532  input.vector_variable_size("other_interior_function_types");
533 
534  if (n_other_interior_functions !=
535  input.vector_variable_size("other_interior_function_values"))
536  {
537  libMesh::out << "Error: " << n_other_interior_functions
538  << " other interior function types does not match "
539  << input.vector_variable_size("other_interior_function_values")
540  << " other interior function values." << std::endl;
541 
542  libmesh_error();
543  }
544 
545  if (n_other_interior_functions !=
546  input.vector_variable_size("other_interior_function_subdomains"))
547  {
548  libMesh::out << "Error: " << n_other_interior_functions
549  << " other interior function types does not match "
550  << input.vector_variable_size("other_interior_function_subdomains")
551  << " other interior function subdomains." << std::endl;
552 
553  libmesh_error();
554  }
555 
556  if (n_other_interior_functions !=
557  input.vector_variable_size("other_interior_function_ids"))
558  {
559  libMesh::out << "Error: " << n_other_interior_functions
560  << " other interior function types does not match "
561  << input.vector_variable_size("other_interior_function_ids")
562  << " other interior function ids." << std::endl;
563 
564  libmesh_error();
565  }
566 
567  for (unsigned int i=0; i != n_other_interior_functions; ++i)
568  {
569  const std::string func_type =
570  input("other_interior_function_types", zero_string, i);
571 
572  const std::string func_value =
573  input("other_interior_function_values", empty_string, i);
574 
575  const subdomain_id_type func_subdomain =
576  input("other_interior_condition_subdomains", subdomain_id_type(0), i);
577 
578  const int func_id =
579  input("other_interior_condition_ids", int(0), i);
580 
581  other_interior_functions[func_id][func_subdomain] =
582  (new_function_base(func_type, func_value).release());
583  }
584 
585  GETPOT_REGISTER(other_boundary_function_types);
586  GETPOT_REGISTER(other_boundary_function_values);
587  GETPOT_REGISTER(other_boundary_function_boundaries);
588  GETPOT_REGISTER(other_boundary_function_ids);
589 
590  const unsigned int n_other_boundary_functions =
591  input.vector_variable_size("other_boundary_function_types");
592 
593  if (n_other_boundary_functions !=
594  input.vector_variable_size("other_boundary_function_values"))
595  {
596  libMesh::out << "Error: " << n_other_boundary_functions
597  << " other boundary function types does not match "
598  << input.vector_variable_size("other_boundary_function_values")
599  << " other boundary function values." << std::endl;
600 
601  libmesh_error();
602  }
603 
604  if (n_other_boundary_functions !=
605  input.vector_variable_size("other_boundary_function_boundaries"))
606  {
607  libMesh::out << "Error: " << n_other_boundary_functions
608  << " other boundary function types does not match "
609  << input.vector_variable_size("other_boundary_function_boundaries")
610  << " other boundary function boundaries." << std::endl;
611 
612  libmesh_error();
613  }
614 
615  if (n_other_boundary_functions !=
616  input.vector_variable_size("other_boundary_function_ids"))
617  {
618  libMesh::out << "Error: " << n_other_boundary_functions
619  << " other boundary function types does not match "
620  << input.vector_variable_size("other_boundary_function_ids")
621  << " other boundary function ids." << std::endl;
622 
623  libmesh_error();
624  }
625 
626  for (unsigned int i=0; i != n_other_boundary_functions; ++i)
627  {
628  const std::string func_type =
629  input("other_boundary_function_types", zero_string, i);
630 
631  const std::string func_value =
632  input("other_boundary_function_values", empty_string, i);
633 
634  const boundary_id_type func_boundary =
635  input("other_boundary_function_boundaries", boundary_id_type(0), i);
636 
637  const int func_id =
638  input("other_boundary_function_ids", int(0), i);
639 
640  other_boundary_functions[func_id][func_boundary] =
641  (new_function_base(func_type, func_value).release());
642  }
643 
644  GETPOT_INPUT(run_simulation);
645  GETPOT_INPUT(run_postprocess);
646 
647 
648  GETPOT_REGISTER(fe_family);
649  const unsigned int n_fe_family =
650  std::max(1u, input.vector_variable_size("fe_family"));
651  fe_family.resize(n_fe_family, "LAGRANGE");
652  for (unsigned int i=0; i != n_fe_family; ++i)
653  fe_family[i] = input("fe_family", fe_family[i].c_str(), i);
654  GETPOT_REGISTER(fe_order);
655  const unsigned int n_fe_order =
656  input.vector_variable_size("fe_order");
657  fe_order.resize(n_fe_order, 1);
658  for (unsigned int i=0; i != n_fe_order; ++i)
659  fe_order[i] = input("fe_order", (int)fe_order[i], i);
660  GETPOT_INPUT(extra_quadrature_order);
661 
662 
663  GETPOT_INPUT(analytic_jacobians);
664  GETPOT_INPUT(verify_analytic_jacobians);
665  GETPOT_INPUT(numerical_jacobian_h);
666  GETPOT_INPUT(print_solution_norms);
667  GETPOT_INPUT(print_solutions);
668  GETPOT_INPUT(print_residual_norms);
669  GETPOT_INPUT(print_residuals);
670  GETPOT_INPUT(print_jacobian_norms);
671  GETPOT_INPUT(print_jacobians);
672  GETPOT_INPUT(print_element_solutions);
673  GETPOT_INPUT(print_element_residuals);
674  GETPOT_INPUT(print_element_jacobians);
675 
676 
677  GETPOT_INPUT(use_petsc_snes);
678  GETPOT_INPUT(time_solver_quiet);
679  GETPOT_INPUT(solver_quiet);
680  GETPOT_INPUT(solver_verbose);
681  GETPOT_INPUT(reuse_preconditioner);
682  GETPOT_INPUT(require_residual_reduction);
683  GETPOT_INPUT(min_step_length);
684  GETPOT_INT_INPUT(max_linear_iterations);
685  GETPOT_INT_INPUT(max_nonlinear_iterations);
686  GETPOT_INPUT(relative_step_tolerance);
687  GETPOT_INPUT(relative_residual_tolerance);
688  GETPOT_INPUT(absolute_residual_tolerance);
689  GETPOT_INPUT(initial_linear_tolerance);
690  GETPOT_INPUT(minimum_linear_tolerance);
691  GETPOT_INPUT(linear_tolerance_multiplier);
692 
693 
694  GETPOT_INT_INPUT(initial_sobolev_order);
695  GETPOT_INT_INPUT(initial_extra_quadrature);
696  GETPOT_INPUT(refine_uniformly);
697  GETPOT_INPUT(indicator_type);
698  GETPOT_INPUT(patch_reuse);
699  GETPOT_INT_INPUT(sobolev_order);
700 
701  GETPOT_INPUT(system_config_file);
702 
703  std::vector<std::string> bad_variables =
704  input.unidentified_arguments(variable_names);
705 
706  if (this->comm().rank() == 0 && !bad_variables.empty())
707  {
708  libMesh::err << "ERROR: Unrecognized variables:" << std::endl;
709  for (unsigned int i = 0; i != bad_variables.size(); ++i)
710  libMesh::err << bad_variables[i] << std::endl;
711  libMesh::err << "Not found among recognized variables:" << std::endl;
712  for (unsigned int i = 0; i != variable_names.size(); ++i)
713  libMesh::err << variable_names[i] << std::endl;
714  libmesh_error();
715  }
716 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
unsigned int nelem_target
Definition: femparameters.h:56
OStreamProxy err
libMesh::Real timesolver_upper_tolerance
Definition: femparameters.h:38
libMesh::Real domain_ymin
Definition: femparameters.h:48
bool print_solution_norms
libMesh::Real deltat
Definition: femparameters.h:38
unsigned int initial_timestep
Definition: femparameters.h:34
bool print_element_residuals
libMesh::Real initial_linear_tolerance
unsigned int dimension
Definition: femparameters.h:45
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
unsigned int coarserefinements
Definition: femparameters.h:51
bool reuse_preconditioner
libMesh::Real timesolver_tolerance
Definition: femparameters.h:38
std::string domaintype
Definition: femparameters.h:46
ConstFunction that simply returns 0.
Definition: zero_function.h:35
A Function generated (via FParser) by parsing a mathematical expression.
libMesh::Real elementorder
Definition: femparameters.h:47
bool analytic_jacobians
VectorValue< Real > RealVectorValue
Useful typedefs to allow transparent switching between Real and Complex data types.
TestClass subdomain_id_type
Based on the 4-byte comment warning above, this probably doesn&#39;t work with exodusII at all...
Definition: id_types.h:43
libMesh::Real domain_edge_length
Definition: femparameters.h:49
std::string timesolver_core
Definition: femparameters.h:37
std::string system_config_file
bool require_residual_reduction
unsigned int coarsegridy
Definition: femparameters.h:50
libMesh::Real domain_zmin
Definition: femparameters.h:48
bool print_residual_norms
libMesh::Real timesolver_maxgrowth
Definition: femparameters.h:38
std::map< libMesh::boundary_id_type, libMesh::FunctionBase< libMesh::Number > * > neumann_conditions
Definition: femparameters.h:87
libMesh::Real relative_step_tolerance
std::string indicator_type
static const Real TOLERANCE
The libMesh namespace provides an interface to certain functionality in the library.
The definition of a periodic boundary.
libMesh::Real refine_fraction
Definition: femparameters.h:58
UniquePtr< FunctionBase< Number > > new_function_base(const std::string &func_type, const std::string &func_value)
libMesh::Real domain_edge_height
Definition: femparameters.h:49
long double max(long double a, double b)
unsigned int initial_extra_quadrature
std::vector< unsigned int > fe_order
std::string mesh_redistribute_func
Definition: femparameters.h:52
libMesh::Real coarsen_threshold
Definition: femparameters.h:58
bool write_gmv_error
Definition: femparameters.h:65
libMesh::Real end_time
Definition: femparameters.h:38
unsigned int max_linear_iterations
int8_t boundary_id_type
Definition: id_types.h:51
std::map< int, std::map< libMesh::subdomain_id_type, libMesh::FunctionBase< libMesh::Number > * > > other_interior_functions
Definition: femparameters.h:94
std::map< libMesh::subdomain_id_type, libMesh::FunctionBase< libMesh::Number > * > initial_conditions
Definition: femparameters.h:85
std::string domainfile
Definition: femparameters.h:46
libMesh::Real minimum_linear_tolerance
std::vector< libMesh::PeriodicBoundary > periodic_boundaries
Definition: femparameters.h:81
std::vector< std::string > fe_family
std::map< int, std::map< libMesh::boundary_id_type, libMesh::FunctionBase< libMesh::Number > * > > other_boundary_functions
Definition: femparameters.h:97
unsigned int coarsegridx
Definition: femparameters.h:50
libMesh::Real verify_analytic_jacobians
unsigned int write_interval
Definition: femparameters.h:64
unsigned int max_adaptivesteps
Definition: femparameters.h:59
std::map< libMesh::boundary_id_type, std::vector< unsigned int > > neumann_condition_variables
Definition: femparameters.h:90
int extra_quadrature_order
libMesh::Real global_tolerance
Definition: femparameters.h:57
std::string elementtype
Definition: femparameters.h:46
libMesh::Real absolute_residual_tolerance
This class forms the base class for all other classes that are expected to be implemented in parallel...
bool write_exodus_error
Definition: femparameters.h:65
libMesh::Real relative_residual_tolerance
bool print_element_solutions
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool print_jacobian_norms
unsigned int extrarefinements
Definition: femparameters.h:51
const Parallel::Communicator & comm() const
OStreamProxy out
libMesh::Real min_step_length
unsigned int n_timesteps
Definition: femparameters.h:34
libMesh::Real timesolver_theta
Definition: femparameters.h:38
void read(GetPot &input, const std::vector< std::string > *other_variable_names=libmesh_nullptr)
unsigned int coarsegridz
Definition: femparameters.h:50
bool print_element_jacobians
unsigned int max_nonlinear_iterations
std::map< libMesh::boundary_id_type, std::vector< unsigned int > > dirichlet_condition_variables
Definition: femparameters.h:90
unsigned int initial_adaptivesteps
Definition: femparameters.h:60
bool time_solver_quiet
unsigned int initial_sobolev_order
libMesh::Real numerical_jacobian_h
std::map< libMesh::boundary_id_type, libMesh::FunctionBase< libMesh::Number > * > dirichlet_conditions
Definition: femparameters.h:87
std::vector< libMesh::FEMNormType > timesolver_norm
Definition: femparameters.h:41
unsigned int sobolev_order
std::vector< std::string > system_types
Definition: femparameters.h:74
libMesh::Real linear_tolerance_multiplier
A simple smart pointer providing strict ownership semantics.
Definition: auto_ptr.h:210
FEMParameters(const libMesh::Parallel::Communicator &comm_in)
Definition: femparameters.C:24
bool write_tecplot_error
Definition: femparameters.h:65
libMesh::Real domain_edge_width
Definition: femparameters.h:49
libMesh::Real coarsen_fraction
Definition: femparameters.h:58
libMesh::Real steadystate_tolerance
Definition: femparameters.h:38
libMesh::Real domain_xmin
Definition: femparameters.h:48
unsigned int deltat_reductions
Definition: femparameters.h:36