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 #ifdef LIBMESH_ENABLE_PERIODIC
62  periodic_boundaries(0),
63 #endif
64 
65  run_simulation(true), run_postprocess(false),
66 
67  fe_family(1, "LAGRANGE"), fe_order(1, 1),
68  extra_quadrature_order(0),
69 
70  analytic_jacobians(true), verify_analytic_jacobians(0.0),
71  numerical_jacobian_h(TOLERANCE),
72  print_solution_norms(false), print_solutions(false),
73  print_residual_norms(false), print_residuals(false),
74  print_jacobian_norms(false), print_jacobians(false),
75  print_element_solutions(false),
76  print_element_residuals(false),
77  print_element_jacobians(false),
78 
79  use_petsc_snes(false),
80  time_solver_quiet(true), solver_quiet(true), solver_verbose(false),
81  reuse_preconditioner(true),
82  require_residual_reduction(true),
83  min_step_length(1e-5),
84  max_linear_iterations(200000), max_nonlinear_iterations(20),
85  relative_step_tolerance(1.e-7), relative_residual_tolerance(1.e-10),
86  absolute_residual_tolerance(1.e-10),
87  initial_linear_tolerance(1.e-3), minimum_linear_tolerance(TOLERANCE*TOLERANCE),
88  linear_tolerance_multiplier(1.e-3),
89 
90  initial_sobolev_order(1),
91  initial_extra_quadrature(0),
92  refine_uniformly(false),
93  indicator_type("kelly"), patch_reuse(true), sobolev_order(1)
94 {
95 }
96 
98 {
99  for (std::map<subdomain_id_type, FunctionBase<Number> *>::iterator
100  i = initial_conditions.begin(); i != initial_conditions.end();
101  ++i)
102  delete i->second;
103 
104  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
105  i = dirichlet_conditions.begin(); i != dirichlet_conditions.end();
106  ++i)
107  delete i->second;
108 
109  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
110  i = neumann_conditions.begin(); i != neumann_conditions.end();
111  ++i)
112  delete i->second;
113 
114  for (std::map<int,
116  >::iterator
117  i = other_boundary_functions.begin(); i != other_boundary_functions.end();
118  ++i)
119  for (std::map<boundary_id_type, FunctionBase<Number> *>::iterator
120  j = i->second.begin(); j != i->second.end();
121  ++j)
122  delete j->second;
123 
124  for (std::map<int,
126  >::iterator
127  i = other_interior_functions.begin(); i != other_interior_functions.end();
128  ++i)
129  for (std::map<subdomain_id_type, FunctionBase<Number> *>::iterator
130  j = i->second.begin(); j != i->second.end();
131  ++j)
132  delete j->second;
133 }
134 
135 
136 UniquePtr<FunctionBase<Number>> new_function_base(const std::string & func_type,
137  const std::string & func_value)
138 {
139  if (func_type == "parsed")
141  else if (func_type == "zero")
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 (std::size_t 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 #ifdef LIBMESH_ENABLE_PERIODIC
274  GETPOT_REGISTER(periodic_boundaries);
275  const unsigned int n_periodic_bcs =
276  input.vector_variable_size("periodic_boundaries");
277 
278  if (n_periodic_bcs)
279  {
280  if (domaintype != "square" &&
281  domaintype != "cylinder" &&
282  domaintype != "file" &&
283  domaintype != "od2")
284  {
285  libMesh::out << "Periodic boundaries need rectilinear domains" << std::endl;;
286  libmesh_error();
287  }
288  for (unsigned int i=0; i != n_periodic_bcs; ++i)
289  {
290  unsigned int myboundary =
291  input("periodic_boundaries", -1, i);
292  unsigned int pairedboundary = 0;
293  RealVectorValue translation_vector;
294  if (dimension == 2)
295  switch (myboundary)
296  {
297  case 0:
298  pairedboundary = 2;
299  translation_vector = RealVectorValue(0., domain_edge_length);
300  break;
301  case 1:
302  pairedboundary = 3;
303  translation_vector = RealVectorValue(-domain_edge_width, 0);
304  break;
305  default:
306  libMesh::out << "Unrecognized periodic boundary id " <<
307  myboundary << std::endl;;
308  libmesh_error();
309  }
310  else if (dimension == 3)
311  switch (myboundary)
312  {
313  case 0:
314  pairedboundary = 5;
315  translation_vector = RealVectorValue(Real(0), Real(0), domain_edge_height);
316  break;
317  case 1:
318  pairedboundary = 3;
319  translation_vector = RealVectorValue(Real(0), domain_edge_length, Real(0));
320  break;
321  case 2:
322  pairedboundary = 4;
323  translation_vector = RealVectorValue(-domain_edge_width, Real(0), Real(0));
324  break;
325  default:
326  libMesh::out << "Unrecognized periodic boundary id " <<
327  myboundary << std::endl;;
328  libmesh_error();
329  }
330  periodic_boundaries.push_back(PeriodicBoundary(translation_vector));
331  periodic_boundaries[i].myboundary = myboundary;
332  periodic_boundaries[i].pairedboundary = pairedboundary;
333  }
334  }
335 #endif // LIBMESH_ENABLE_PERIODIC
336 
337  // Use std::string inputs so GetPot doesn't have to make a bunch
338  // of internal C string copies
339  std::string zero_string = "zero";
340  std::string empty_string = "";
341 
342  GETPOT_REGISTER(dirichlet_condition_types);
343  GETPOT_REGISTER(dirichlet_condition_values);
344  GETPOT_REGISTER(dirichlet_condition_boundaries);
345  GETPOT_REGISTER(dirichlet_condition_variables);
346 
347  const unsigned int n_dirichlet_conditions=
348  input.vector_variable_size("dirichlet_condition_types");
349 
350  if (n_dirichlet_conditions !=
351  input.vector_variable_size("dirichlet_condition_values"))
352  {
353  libMesh::out << "Error: " << n_dirichlet_conditions
354  << " Dirichlet condition types does not match "
355  << input.vector_variable_size("dirichlet_condition_values")
356  << " Dirichlet condition values." << std::endl;
357 
358  libmesh_error();
359  }
360 
361  if (n_dirichlet_conditions !=
362  input.vector_variable_size("dirichlet_condition_boundaries"))
363  {
364  libMesh::out << "Error: " << n_dirichlet_conditions
365  << " Dirichlet condition types does not match "
366  << input.vector_variable_size("dirichlet_condition_boundaries")
367  << " Dirichlet condition boundaries." << std::endl;
368 
369  libmesh_error();
370  }
371 
372  if (n_dirichlet_conditions !=
373  input.vector_variable_size("dirichlet_condition_variables"))
374  {
375  libMesh::out << "Error: " << n_dirichlet_conditions
376  << " Dirichlet condition types does not match "
377  << input.vector_variable_size("dirichlet_condition_variables")
378  << " Dirichlet condition variables sets." << std::endl;
379 
380  libmesh_error();
381  }
382 
383  for (unsigned int i=0; i != n_dirichlet_conditions; ++i)
384  {
385  const std::string func_type =
386  input("dirichlet_condition_types", zero_string, i);
387 
388  const std::string func_value =
389  input("dirichlet_condition_values", empty_string, i);
390 
391  const boundary_id_type func_boundary =
392  input("dirichlet_condition_boundaries", boundary_id_type(0), i);
393 
394  dirichlet_conditions[func_boundary] =
395  (new_function_base(func_type, func_value).release());
396 
397  const std::string variable_set =
398  input("dirichlet_condition_variables", empty_string, i);
399 
400  for (std::size_t i=0; i != variable_set.size(); ++i)
401  {
402  if (variable_set[i] == '1')
403  dirichlet_condition_variables[func_boundary].push_back(i);
404  else if (variable_set[i] != '0')
405  {
406  libMesh::out << "Unable to understand Dirichlet variable set"
407  << variable_set << std::endl;
408  libmesh_error();
409  }
410  }
411  }
412 
413  GETPOT_REGISTER(neumann_condition_types);
414  GETPOT_REGISTER(neumann_condition_values);
415  GETPOT_REGISTER(neumann_condition_boundaries);
416  GETPOT_REGISTER(neumann_condition_variables);
417 
418  const unsigned int n_neumann_conditions=
419  input.vector_variable_size("neumann_condition_types");
420 
421  if (n_neumann_conditions !=
422  input.vector_variable_size("neumann_condition_values"))
423  {
424  libMesh::out << "Error: " << n_neumann_conditions
425  << " Neumann condition types does not match "
426  << input.vector_variable_size("neumann_condition_values")
427  << " Neumann condition values." << std::endl;
428 
429  libmesh_error();
430  }
431 
432  if (n_neumann_conditions !=
433  input.vector_variable_size("neumann_condition_boundaries"))
434  {
435  libMesh::out << "Error: " << n_neumann_conditions
436  << " Neumann condition types does not match "
437  << input.vector_variable_size("neumann_condition_boundaries")
438  << " Neumann condition boundaries." << std::endl;
439 
440  libmesh_error();
441  }
442 
443  if (n_neumann_conditions !=
444  input.vector_variable_size("neumann_condition_variables"))
445  {
446  libMesh::out << "Error: " << n_neumann_conditions
447  << " Neumann condition types does not match "
448  << input.vector_variable_size("neumann_condition_variables")
449  << " Neumann condition variables sets." << std::endl;
450 
451  libmesh_error();
452  }
453 
454  for (unsigned int i=0; i != n_neumann_conditions; ++i)
455  {
456  const std::string func_type =
457  input("neumann_condition_types", zero_string, i);
458 
459  const std::string func_value =
460  input("neumann_condition_values", empty_string, i);
461 
462  const boundary_id_type func_boundary =
463  input("neumann_condition_boundaries", boundary_id_type(0), i);
464 
465  neumann_conditions[func_boundary] =
466  (new_function_base(func_type, func_value).release());
467 
468  const std::string variable_set =
469  input("neumann_condition_variables", empty_string, i);
470 
471  for (std::size_t i=0; i != variable_set.size(); ++i)
472  {
473  if (variable_set[i] == '1')
474  neumann_condition_variables[func_boundary].push_back(i);
475  else if (variable_set[i] != '0')
476  {
477  libMesh::out << "Unable to understand Neumann variable set"
478  << variable_set << std::endl;
479  libmesh_error();
480  }
481  }
482  }
483 
484  GETPOT_REGISTER(initial_condition_types);
485  GETPOT_REGISTER(initial_condition_values);
486  GETPOT_REGISTER(initial_condition_subdomains);
487 
488  const unsigned int n_initial_conditions=
489  input.vector_variable_size("initial_condition_types");
490 
491  if (n_initial_conditions !=
492  input.vector_variable_size("initial_condition_values"))
493  {
494  libMesh::out << "Error: " << n_initial_conditions
495  << " initial condition types does not match "
496  << input.vector_variable_size("initial_condition_values")
497  << " initial condition values." << std::endl;
498 
499  libmesh_error();
500  }
501 
502  if (n_initial_conditions !=
503  input.vector_variable_size("initial_condition_subdomains"))
504  {
505  libMesh::out << "Error: " << n_initial_conditions
506  << " initial condition types does not match "
507  << input.vector_variable_size("initial_condition_subdomains")
508  << " initial condition subdomains." << std::endl;
509 
510  libmesh_error();
511  }
512 
513  for (unsigned int i=0; i != n_initial_conditions; ++i)
514  {
515  const std::string func_type =
516  input("initial_condition_types", zero_string, i);
517 
518  const std::string func_value =
519  input("initial_condition_values", empty_string, i);
520 
521  const subdomain_id_type func_subdomain =
522  input("initial_condition_subdomains", subdomain_id_type(0), i);
523 
524  initial_conditions[func_subdomain] =
525  (new_function_base(func_type, func_value).release());
526  }
527 
528  GETPOT_REGISTER(other_interior_function_types);
529  GETPOT_REGISTER(other_interior_function_values);
530  GETPOT_REGISTER(other_interior_function_subdomains);
531  GETPOT_REGISTER(other_interior_function_ids);
532 
533  const unsigned int n_other_interior_functions =
534  input.vector_variable_size("other_interior_function_types");
535 
536  if (n_other_interior_functions !=
537  input.vector_variable_size("other_interior_function_values"))
538  {
539  libMesh::out << "Error: " << n_other_interior_functions
540  << " other interior function types does not match "
541  << input.vector_variable_size("other_interior_function_values")
542  << " other interior function values." << std::endl;
543 
544  libmesh_error();
545  }
546 
547  if (n_other_interior_functions !=
548  input.vector_variable_size("other_interior_function_subdomains"))
549  {
550  libMesh::out << "Error: " << n_other_interior_functions
551  << " other interior function types does not match "
552  << input.vector_variable_size("other_interior_function_subdomains")
553  << " other interior function subdomains." << std::endl;
554 
555  libmesh_error();
556  }
557 
558  if (n_other_interior_functions !=
559  input.vector_variable_size("other_interior_function_ids"))
560  {
561  libMesh::out << "Error: " << n_other_interior_functions
562  << " other interior function types does not match "
563  << input.vector_variable_size("other_interior_function_ids")
564  << " other interior function ids." << std::endl;
565 
566  libmesh_error();
567  }
568 
569  for (unsigned int i=0; i != n_other_interior_functions; ++i)
570  {
571  const std::string func_type =
572  input("other_interior_function_types", zero_string, i);
573 
574  const std::string func_value =
575  input("other_interior_function_values", empty_string, i);
576 
577  const subdomain_id_type func_subdomain =
578  input("other_interior_condition_subdomains", subdomain_id_type(0), i);
579 
580  const int func_id =
581  input("other_interior_condition_ids", int(0), i);
582 
583  other_interior_functions[func_id][func_subdomain] =
584  (new_function_base(func_type, func_value).release());
585  }
586 
587  GETPOT_REGISTER(other_boundary_function_types);
588  GETPOT_REGISTER(other_boundary_function_values);
589  GETPOT_REGISTER(other_boundary_function_boundaries);
590  GETPOT_REGISTER(other_boundary_function_ids);
591 
592  const unsigned int n_other_boundary_functions =
593  input.vector_variable_size("other_boundary_function_types");
594 
595  if (n_other_boundary_functions !=
596  input.vector_variable_size("other_boundary_function_values"))
597  {
598  libMesh::out << "Error: " << n_other_boundary_functions
599  << " other boundary function types does not match "
600  << input.vector_variable_size("other_boundary_function_values")
601  << " other boundary function values." << std::endl;
602 
603  libmesh_error();
604  }
605 
606  if (n_other_boundary_functions !=
607  input.vector_variable_size("other_boundary_function_boundaries"))
608  {
609  libMesh::out << "Error: " << n_other_boundary_functions
610  << " other boundary function types does not match "
611  << input.vector_variable_size("other_boundary_function_boundaries")
612  << " other boundary function boundaries." << std::endl;
613 
614  libmesh_error();
615  }
616 
617  if (n_other_boundary_functions !=
618  input.vector_variable_size("other_boundary_function_ids"))
619  {
620  libMesh::out << "Error: " << n_other_boundary_functions
621  << " other boundary function types does not match "
622  << input.vector_variable_size("other_boundary_function_ids")
623  << " other boundary function ids." << std::endl;
624 
625  libmesh_error();
626  }
627 
628  for (unsigned int i=0; i != n_other_boundary_functions; ++i)
629  {
630  const std::string func_type =
631  input("other_boundary_function_types", zero_string, i);
632 
633  const std::string func_value =
634  input("other_boundary_function_values", empty_string, i);
635 
636  const boundary_id_type func_boundary =
637  input("other_boundary_function_boundaries", boundary_id_type(0), i);
638 
639  const int func_id =
640  input("other_boundary_function_ids", int(0), i);
641 
642  other_boundary_functions[func_id][func_boundary] =
643  (new_function_base(func_type, func_value).release());
644  }
645 
646  GETPOT_INPUT(run_simulation);
647  GETPOT_INPUT(run_postprocess);
648 
649 
650  GETPOT_REGISTER(fe_family);
651  const unsigned int n_fe_family =
652  std::max(1u, input.vector_variable_size("fe_family"));
653  fe_family.resize(n_fe_family, "LAGRANGE");
654  for (unsigned int i=0; i != n_fe_family; ++i)
655  fe_family[i] = input("fe_family", fe_family[i].c_str(), i);
656  GETPOT_REGISTER(fe_order);
657  const unsigned int n_fe_order =
658  input.vector_variable_size("fe_order");
659  fe_order.resize(n_fe_order, 1);
660  for (unsigned int i=0; i != n_fe_order; ++i)
661  fe_order[i] = input("fe_order", (int)fe_order[i], i);
662  GETPOT_INPUT(extra_quadrature_order);
663 
664 
665  GETPOT_INPUT(analytic_jacobians);
666  GETPOT_INPUT(verify_analytic_jacobians);
667  GETPOT_INPUT(numerical_jacobian_h);
668  GETPOT_INPUT(print_solution_norms);
669  GETPOT_INPUT(print_solutions);
670  GETPOT_INPUT(print_residual_norms);
671  GETPOT_INPUT(print_residuals);
672  GETPOT_INPUT(print_jacobian_norms);
673  GETPOT_INPUT(print_jacobians);
674  GETPOT_INPUT(print_element_solutions);
675  GETPOT_INPUT(print_element_residuals);
676  GETPOT_INPUT(print_element_jacobians);
677 
678 
679  GETPOT_INPUT(use_petsc_snes);
680  GETPOT_INPUT(time_solver_quiet);
681  GETPOT_INPUT(solver_quiet);
682  GETPOT_INPUT(solver_verbose);
683  GETPOT_INPUT(reuse_preconditioner);
684  GETPOT_INPUT(require_residual_reduction);
685  GETPOT_INPUT(min_step_length);
686  GETPOT_INT_INPUT(max_linear_iterations);
687  GETPOT_INT_INPUT(max_nonlinear_iterations);
688  GETPOT_INPUT(relative_step_tolerance);
689  GETPOT_INPUT(relative_residual_tolerance);
690  GETPOT_INPUT(absolute_residual_tolerance);
691  GETPOT_INPUT(initial_linear_tolerance);
692  GETPOT_INPUT(minimum_linear_tolerance);
693  GETPOT_INPUT(linear_tolerance_multiplier);
694 
695 
696  GETPOT_INT_INPUT(initial_sobolev_order);
697  GETPOT_INT_INPUT(initial_extra_quadrature);
698  GETPOT_INPUT(refine_uniformly);
699  GETPOT_INPUT(indicator_type);
700  GETPOT_INPUT(patch_reuse);
701  GETPOT_INT_INPUT(sobolev_order);
702 
703  GETPOT_INPUT(system_config_file);
704 
705  std::vector<std::string> bad_variables =
706  input.unidentified_arguments(variable_names);
707 
708  if (this->comm().rank() == 0 && !bad_variables.empty())
709  {
710  libMesh::err << "ERROR: Unrecognized variables:" << std::endl;
711  for (std::size_t i = 0; i != bad_variables.size(); ++i)
712  libMesh::err << bad_variables[i] << std::endl;
713  libMesh::err << "Not found among recognized variables:" << std::endl;
714  for (std::size_t i = 0; i != variable_names.size(); ++i)
715  libMesh::err << variable_names[i] << std::endl;
716  libmesh_error();
717  }
718 }
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::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
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
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