libMesh
Functions
adjoints_ex2.C File Reference

Go to the source code of this file.

Functions

void write_output (EquationSystems &es, unsigned int a_step, std::string solution_type, FEMParameters &param)
 
void set_system_parameters (LaplaceSystem &system, FEMParameters &param)
 
UniquePtr< MeshRefinementbuild_mesh_refinement (MeshBase &mesh, FEMParameters &param)
 
UniquePtr< ErrorEstimatorbuild_error_estimator (FEMParameters &param)
 
int main (int argc, char **argv)
 

Function Documentation

UniquePtr<ErrorEstimator> build_error_estimator ( FEMParameters param)

Definition at line 243 of file adjoints_ex2.C.

References libMesh::AdjointResidualErrorEstimator::dual_error_estimator(), libMesh::AdjointResidualErrorEstimator::error_plot_suffix, libMesh::H1_SEMINORM, FEMParameters::indicator_type, libMesh::out, and libMesh::AdjointResidualErrorEstimator::primal_error_estimator().

Referenced by main().

244 {
245  if (param.indicator_type == "kelly")
246  {
247  libMesh::out << "Using Kelly Error Estimator" << std::endl;
248 
250  }
251  else if (param.indicator_type == "adjoint_residual")
252  {
253  libMesh::out << "Using Adjoint Residual Error Estimator with Patch Recovery Weights" << std::endl << std::endl;
254 
255  AdjointResidualErrorEstimator * adjoint_residual_estimator = new AdjointResidualErrorEstimator;
256 
257  adjoint_residual_estimator->error_plot_suffix = "error.gmv";
258 
260  adjoint_residual_estimator->primal_error_estimator().reset(p1);
261 
263  adjoint_residual_estimator->dual_error_estimator().reset(p2);
264 
265  adjoint_residual_estimator->primal_error_estimator()->error_norm.set_type(0, H1_SEMINORM);
266 
267  adjoint_residual_estimator->dual_error_estimator()->error_norm.set_type(0, H1_SEMINORM);
268 
269  return UniquePtr<ErrorEstimator>(adjoint_residual_estimator);
270  }
271  else
272  libmesh_error_msg("Unknown indicator_type = " << param.indicator_type);
273 }
This class implements a goal oriented error indicator, by weighting residual-based estimates from the...
std::string indicator_type
UniquePtr< ErrorEstimator > & primal_error_estimator()
Access to the "subestimator" (default: PatchRecovery) to use on the primal/forward solution...
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
std::string error_plot_suffix
To aid in investigating error estimator behavior, set this string to a suffix with which to plot (pre...
This class implements the Kelly error indicator which is based on the flux jumps between elements...
This class implements the Patch Recovery error indicator.
OStreamProxy out
UniquePtr< ErrorEstimator > & dual_error_estimator()
Access to the "subestimator" (default: PatchRecovery) to use on the dual/adjoint solution.
UniquePtr<MeshRefinement> build_mesh_refinement ( MeshBase mesh,
FEMParameters param 
)

Definition at line 220 of file adjoints_ex2.C.

References libMesh::MeshRefinement::absolute_global_tolerance(), libMesh::MeshRefinement::coarsen_by_parents(), FEMParameters::coarsen_fraction, libMesh::MeshRefinement::coarsen_fraction(), FEMParameters::coarsen_threshold, libMesh::MeshRefinement::coarsen_threshold(), FEMParameters::global_tolerance, FEMParameters::nelem_target, libMesh::MeshRefinement::nelem_target(), FEMParameters::refine_fraction, and libMesh::MeshRefinement::refine_fraction().

Referenced by main().

222 {
223  MeshRefinement * mesh_refinement = new MeshRefinement(mesh);
224  mesh_refinement->coarsen_by_parents() = true;
225  mesh_refinement->absolute_global_tolerance() = param.global_tolerance;
226  mesh_refinement->nelem_target() = param.nelem_target;
227  mesh_refinement->refine_fraction() = param.refine_fraction;
228  mesh_refinement->coarsen_fraction() = param.coarsen_fraction;
229  mesh_refinement->coarsen_threshold() = param.coarsen_threshold;
230 
231  return UniquePtr<MeshRefinement>(mesh_refinement);
232 }
unsigned int nelem_target
Definition: femparameters.h:56
Real & absolute_global_tolerance()
If absolute_global_tolerance is set to a nonzero value, methods like flag_elements_by_global_toleranc...
Real & coarsen_fraction()
The coarsen_fraction sets either a desired target or a desired maximum number of elements to flag for...
Real & refine_fraction()
The refine_fraction sets either a desired target or a desired maximum number of elements to flag for ...
libMesh::Real refine_fraction
Definition: femparameters.h:58
bool & coarsen_by_parents()
If coarsen_by_parents is true, complete groups of sibling elements (elements with the same parent) wi...
dof_id_type & nelem_target()
If nelem_target is set to a nonzero value, methods like flag_elements_by_nelem_target() will attempt ...
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
libMesh::Real coarsen_threshold
Definition: femparameters.h:58
This is the MeshRefinement class.
Real & coarsen_threshold()
The coarsen_threshold provides hysteresis in AMR/C strategies.
libMesh::Real global_tolerance
Definition: femparameters.h:57
libMesh::Real coarsen_fraction
Definition: femparameters.h:58
int main ( int  argc,
char **  argv 
)

Definition at line 276 of file adjoints_ex2.C.

References std::abs(), libMesh::QoISet::add_indices(), libMesh::EquationSystems::add_system(), libMesh::ImplicitSystem::adjoint_qoi_parameter_sensitivity(), libMesh::MeshBase::all_second_order(), libMesh::DifferentiableQoI::assemble_qoi_sides, libMesh::DifferentiableSystem::attach_qoi(), build_error_estimator(), build_mesh_refinement(), libMesh::LibMeshInit::comm(), libMesh::default_solver_package(), libMesh::err, libMesh::System::get_adjoint_solution(), LaplaceSystem::get_parameter_vector(), libMesh::TriangleWrapper::init(), libMesh::EquationSystems::init(), libMesh::INVALID_SOLVER_PACKAGE, mesh, libMesh::EquationSystems::n_active_dofs(), libMesh::MeshBase::n_active_elem(), libMesh::out, libMesh::EquationSystems::print_info(), libMesh::MeshBase::print_info(), libMesh::ParallelObject::processor_id(), libMesh::Parallel::Communicator::rank(), FEMParameters::read(), libMesh::MeshBase::read(), libMesh::EquationSystems::reinit(), libMesh::System::set_adjoint_already_solved(), set_system_parameters(), libMesh::System::solution, libMesh::FEMSystem::solve(), libMesh::NumericVector< T >::swap(), libMesh::MeshRefinement::uniformly_refine(), and write_output().

277 {
278  // Initialize libMesh.
279  LibMeshInit init (argc, argv);
280 
281  // This example requires a linear solver package.
282  libmesh_example_requires(libMesh::default_solver_package() != INVALID_SOLVER_PACKAGE,
283  "--enable-petsc, --enable-trilinos, or --enable-eigen");
284 
285  // Skip adaptive examples on a non-adaptive libMesh build
286 #ifndef LIBMESH_ENABLE_AMR
287  libmesh_example_requires(false, "--enable-amr");
288 #else
289 
290  libMesh::out << "Started " << argv[0] << std::endl;
291 
292  // Make sure the general input file exists, and parse it
293  {
294  std::ifstream i("general.in");
295  if (!i)
296  libmesh_error_msg('[' << init.comm().rank() << "] Can't find general.in; exiting early.");
297  }
298  GetPot infile("general.in");
299 
300  // Read in parameters from the input file
301  FEMParameters param(init.comm());
302  param.read(infile);
303 
304  // Skip this default-2D example if libMesh was compiled as 1D-only.
305  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");
306 
307  // Create a mesh, with dimension to be overridden later, distributed
308  // across the default MPI communicator.
309  Mesh mesh(init.comm());
310 
311  // And an object to refine it
312  UniquePtr<MeshRefinement> mesh_refinement =
313  build_mesh_refinement(mesh, param);
314 
315  // And an EquationSystems to run on it
316  EquationSystems equation_systems (mesh);
317 
318  libMesh::out << "Reading in and building the mesh" << std::endl;
319 
320  // Read in the mesh
321  mesh.read(param.domainfile.c_str());
322  // Make all the elements of the mesh second order so we can compute
323  // with a higher order basis
325 
326  // Create a mesh refinement object to do the initial uniform refinements
327  // on the coarse grid read in from lshaped.xda
328  MeshRefinement initial_uniform_refinements(mesh);
329  initial_uniform_refinements.uniformly_refine(param.coarserefinements);
330 
331  libMesh::out << "Building system" << std::endl;
332 
333  // Build the FEMSystem
334  LaplaceSystem & system = equation_systems.add_system<LaplaceSystem> ("LaplaceSystem");
335 
336  QoISet qois;
337 
338  std::vector<unsigned int> qoi_indices;
339  qoi_indices.push_back(0);
340  qois.add_indices(qoi_indices);
341 
342  qois.set_weight(0, 0.5);
343 
344  // Put some scope here to test that the cloning is working right
345  {
346  LaplaceQoI qoi;
347  system.attach_qoi(&qoi);
348  }
349 
350  // Set its parameters
351  set_system_parameters(system, param);
352 
353  libMesh::out << "Initializing systems" << std::endl;
354 
355  equation_systems.init ();
356 
357  // Print information about the mesh and system to the screen.
358  mesh.print_info();
359  equation_systems.print_info();
360 
361  {
362  // Adaptively solve the timestep
363  unsigned int a_step = 0;
364  for (; a_step != param.max_adaptivesteps; ++a_step)
365  {
366  // We can't adapt to both a tolerance and a
367  // target mesh size
368  if (param.global_tolerance != 0.)
369  libmesh_assert_equal_to (param.nelem_target, 0);
370  // If we aren't adapting to a tolerance we need a
371  // target mesh size
372  else
373  libmesh_assert_greater (param.nelem_target, 0);
374 
375  // Solve the forward problem
376  system.solve();
377 
378  // Write out the computed primal solution
379  write_output(equation_systems, a_step, "primal", param);
380 
381  // Get a pointer to the primal solution vector
382  NumericVector<Number> & primal_solution = *system.solution;
383 
384  // A SensitivityData object to hold the qois and parameters
385  SensitivityData sensitivities(qois, system, system.get_parameter_vector());
386 
387  // Make sure we get the contributions to the adjoint RHS from the sides
388  system.assemble_qoi_sides = true;
389 
390  // Here we solve the adjoint problem inside the adjoint_qoi_parameter_sensitivity
391  // function, so we have to set the adjoint_already_solved boolean to false
392  system.set_adjoint_already_solved(false);
393 
394  // Compute the sensitivities
395  system.adjoint_qoi_parameter_sensitivity(qois, system.get_parameter_vector(), sensitivities);
396 
397  // Now that we have solved the adjoint, set the adjoint_already_solved boolean to true, so we dont solve unnecessarily in the error estimator
398  system.set_adjoint_already_solved(true);
399 
400  GetPot infile_l_shaped("l-shaped.in");
401 
402  Number sensitivity_QoI_0_0_computed = sensitivities[0][0];
403  Number sensitivity_QoI_0_0_exact = infile_l_shaped("sensitivity_0_0", 0.0);
404  Number sensitivity_QoI_0_1_computed = sensitivities[0][1];
405  Number sensitivity_QoI_0_1_exact = infile_l_shaped("sensitivity_0_1", 0.0);
406 
407  libMesh::out << "Adaptive step "
408  << a_step
409  << ", we have "
410  << mesh.n_active_elem()
411  << " active elements and "
412  << equation_systems.n_active_dofs()
413  << " active dofs."
414  << std::endl;
415 
416  libMesh::out << "Sensitivity of QoI one to Parameter one is "
417  << sensitivity_QoI_0_0_computed
418  << std::endl;
419  libMesh::out << "Sensitivity of QoI one to Parameter two is "
420  << sensitivity_QoI_0_1_computed
421  << std::endl;
422 
423  libMesh::out << "The relative error in sensitivity QoI_0_0 is "
424  << std::setprecision(17)
425  << std::abs(sensitivity_QoI_0_0_computed - sensitivity_QoI_0_0_exact) / std::abs(sensitivity_QoI_0_0_exact)
426  << std::endl;
427 
428  libMesh::out << "The relative error in sensitivity QoI_0_1 is "
429  << std::setprecision(17)
430  << std::abs(sensitivity_QoI_0_1_computed - sensitivity_QoI_0_1_exact) / std::abs(sensitivity_QoI_0_1_exact)
431  << std::endl
432  << std::endl;
433 
434  // Get a pointer to the solution vector of the adjoint problem for QoI 0
435  NumericVector<Number> & dual_solution_0 = system.get_adjoint_solution(0);
436 
437  // Swap the primal and dual solutions so we can write out the adjoint solution
438  primal_solution.swap(dual_solution_0);
439  write_output(equation_systems, a_step, "adjoint_0", param);
440 
441  // Swap back
442  primal_solution.swap(dual_solution_0);
443 
444  // We have to refine either based on reaching an error tolerance or
445  // a number of elements target, which should be verified above
446  // Otherwise we flag elements by error tolerance or nelem target
447 
448  // Uniform refinement
449  if (param.refine_uniformly)
450  {
451  libMesh::out << "Refining Uniformly" << std::endl << std::endl;
452 
453  mesh_refinement->uniformly_refine(1);
454  }
455  // Adaptively refine based on reaching an error tolerance
456  else if (param.global_tolerance >= 0. && param.nelem_target == 0.)
457  {
458  // Now we construct the data structures for the mesh refinement process
459  ErrorVector error;
460 
461  // Build an error estimator object
462  UniquePtr<ErrorEstimator> error_estimator =
463  build_error_estimator(param);
464 
465  // Estimate the error in each element using the Adjoint Residual or Kelly error estimator
466  error_estimator->estimate_error(system, error);
467 
468  mesh_refinement->flag_elements_by_error_tolerance (error);
469 
470  mesh_refinement->refine_and_coarsen_elements();
471  }
472  // Adaptively refine based on reaching a target number of elements
473  else
474  {
475  // Now we construct the data structures for the mesh refinement process
476  ErrorVector error;
477 
478  // Build an error estimator object
479  UniquePtr<ErrorEstimator> error_estimator =
480  build_error_estimator(param);
481 
482  // Estimate the error in each element using the Adjoint Residual or Kelly error estimator
483  error_estimator->estimate_error(system, error);
484 
485  if (mesh.n_active_elem() >= param.nelem_target)
486  {
487  libMesh::out<<"We reached the target number of elements."<<std::endl <<std::endl;
488  break;
489  }
490 
491  mesh_refinement->flag_elements_by_nelem_target (error);
492 
493  mesh_refinement->refine_and_coarsen_elements();
494  }
495 
496  // Dont forget to reinit the system after each adaptive refinement !
497  equation_systems.reinit();
498 
499  libMesh::out << "Refined mesh to "
500  << mesh.n_active_elem()
501  << " active elements and "
502  << equation_systems.n_active_dofs()
503  << " active dofs."
504  << std::endl;
505  }
506 
507  // Do one last solve if necessary
508  if (a_step == param.max_adaptivesteps)
509  {
510  system.solve();
511 
512  write_output(equation_systems, a_step, "primal", param);
513 
514  NumericVector<Number> & primal_solution = *system.solution;
515 
516  SensitivityData sensitivities(qois, system, system.get_parameter_vector());
517 
518  system.assemble_qoi_sides = true;
519 
520  // Here we solve the adjoint problem inside the adjoint_qoi_parameter_sensitivity
521  // function, so we have to set the adjoint_already_solved boolean to false
522  system.set_adjoint_already_solved(false);
523 
524  system.adjoint_qoi_parameter_sensitivity(qois, system.get_parameter_vector(), sensitivities);
525 
526  // Now that we have solved the adjoint, set the adjoint_already_solved boolean to true, so we dont solve unnecessarily in the error estimator
527  system.set_adjoint_already_solved(true);
528 
529  GetPot infile_l_shaped("l-shaped.in");
530 
531  Number sensitivity_QoI_0_0_computed = sensitivities[0][0];
532  Number sensitivity_QoI_0_0_exact = infile_l_shaped("sensitivity_0_0", 0.0);
533  Number sensitivity_QoI_0_1_computed = sensitivities[0][1];
534  Number sensitivity_QoI_0_1_exact = infile_l_shaped("sensitivity_0_1", 0.0);
535 
536  libMesh::out << "Adaptive step "
537  << a_step
538  << ", we have "
539  << mesh.n_active_elem()
540  << " active elements and "
541  << equation_systems.n_active_dofs()
542  << " active dofs."
543  << std::endl;
544 
545  libMesh::out << "Sensitivity of QoI one to Parameter one is "
546  << sensitivity_QoI_0_0_computed
547  << std::endl;
548 
549  libMesh::out << "Sensitivity of QoI one to Parameter two is "
550  << sensitivity_QoI_0_1_computed
551  << std::endl;
552 
553  libMesh::out << "The error in sensitivity QoI_0_0 is "
554  << std::setprecision(17)
555  << std::abs(sensitivity_QoI_0_0_computed - sensitivity_QoI_0_0_exact)/sensitivity_QoI_0_0_exact
556  << std::endl;
557 
558  libMesh::out << "The error in sensitivity QoI_0_1 is "
559  << std::setprecision(17)
560  << std::abs(sensitivity_QoI_0_1_computed - sensitivity_QoI_0_1_exact)/sensitivity_QoI_0_1_exact
561  << std::endl
562  << std::endl;
563 
564  // Hard coded asserts to ensure that the actual numbers we are getting are what they should be
565  libmesh_assert_less(std::abs((sensitivity_QoI_0_0_computed - sensitivity_QoI_0_0_exact)/sensitivity_QoI_0_0_exact), 2.e-4);
566  libmesh_assert_less(std::abs((sensitivity_QoI_0_1_computed - sensitivity_QoI_0_1_exact)/sensitivity_QoI_0_1_exact), 2.e-4);
567 
568  NumericVector<Number> & dual_solution_0 = system.get_adjoint_solution(0);
569 
570  primal_solution.swap(dual_solution_0);
571  write_output(equation_systems, a_step, "adjoint_0", param);
572 
573  primal_solution.swap(dual_solution_0);
574  }
575  }
576 
577  libMesh::err << '[' << system.processor_id()
578  << "] Completing output."
579  << std::endl;
580 
581 #endif // #ifndef LIBMESH_ENABLE_AMR
582 
583  // All done.
584  return 0;
585 }
OStreamProxy err
double abs(double a)
This is the EquationSystems class.
virtual dof_id_type n_active_elem() const =0
void set_adjoint_already_solved(bool setting)
Setter for the adjoint_already_solved boolean.
Definition: system.h:378
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
UniquePtr< ErrorEstimator > build_error_estimator(FEMParameters &param)
Definition: adjoints_ex2.C:243
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
UniquePtr< MeshRefinement > build_mesh_refinement(MeshBase &mesh, FEMParameters &param)
Definition: adjoints_ex2.C:220
virtual void adjoint_qoi_parameter_sensitivity(const QoISet &qoi_indices, const ParameterVector &parameters, SensitivityData &sensitivities) libmesh_override
Solves for the derivative of each of the system&#39;s quantities of interest q in qoi[qoi_indices] with r...
MeshBase & mesh
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
SolverPackage default_solver_package()
Definition: libmesh.C:995
virtual void solve() libmesh_override
Invokes the solver associated with the system.
Definition: fem_system.C:1056
This is the MeshRefinement class.
Data structure for holding completed parameter sensitivity calculations.
virtual void all_second_order(const bool full_ordered=true)=0
Converts a (conforming, non-refined) mesh with linear elements into a mesh with second-order elements...
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
void attach_qoi(DifferentiableQoI *qoi_in)
Attach external QoI object.
Definition: diff_system.h:212
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
bool assemble_qoi_sides
If assemble_qoi_sides is true (it is false by default), the assembly loop for a quantity of interest ...
Definition: diff_qoi.h:82
void write_output(EquationSystems &es, unsigned int a_step, std::string solution_type, FEMParameters &param)
Definition: adjoints_ex2.C:105
OStreamProxy out
void add_indices(const std::vector< unsigned int > &indices)
Add this indices to the set to be calculated.
Definition: qoi_set.C:46
void read(GetPot &input, const std::vector< std::string > *other_variable_names=libmesh_nullptr)
NumericVector< Number > & get_adjoint_solution(unsigned int i=0)
Definition: system.C:989
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
virtual void read(const std::string &name, void *mesh_data=libmesh_nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false)=0
Interfaces for reading/writing a mesh to/from a file.
ParameterVector & get_parameter_vector()
Definition: L-shaped.h:35
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448
void set_system_parameters(LaplaceSystem &system, FEMParameters &param)
Definition: adjoints_ex2.C:167
processor_id_type processor_id() const
void set_system_parameters ( LaplaceSystem system,
FEMParameters param 
)

Definition at line 167 of file adjoints_ex2.C.

References LaplaceSystem::analytic_jacobians(), FEMParameters::analytic_jacobians, libMesh::DiffSolver::continue_after_backtrack_failure, libMesh::DiffSolver::continue_after_max_iterations, LaplaceSystem::fe_family(), FEMParameters::fe_family, LaplaceSystem::fe_order(), FEMParameters::fe_order, FEMParameters::initial_linear_tolerance, libMesh::DiffSolver::initial_linear_tolerance, FEMParameters::linear_tolerance_multiplier, libMesh::NewtonSolver::linear_tolerance_multiplier, FEMParameters::max_linear_iterations, libMesh::DiffSolver::max_linear_iterations, FEMParameters::max_nonlinear_iterations, libMesh::DiffSolver::max_nonlinear_iterations, FEMParameters::min_step_length, FEMParameters::minimum_linear_tolerance, libMesh::DiffSolver::minimum_linear_tolerance, libMesh::NewtonSolver::minsteplength, FEMParameters::print_jacobian_norms, libMesh::DifferentiableSystem::print_jacobian_norms, FEMParameters::print_jacobians, libMesh::DifferentiableSystem::print_jacobians, FEMParameters::print_residual_norms, libMesh::DifferentiableSystem::print_residual_norms, FEMParameters::print_residuals, libMesh::DifferentiableSystem::print_residuals, FEMParameters::print_solution_norms, libMesh::DifferentiableSystem::print_solution_norms, FEMParameters::print_solutions, libMesh::DifferentiableSystem::print_solutions, libMesh::DiffSolver::quiet, FEMParameters::relative_residual_tolerance, libMesh::DiffSolver::relative_residual_tolerance, FEMParameters::relative_step_tolerance, libMesh::DiffSolver::relative_step_tolerance, libMesh::NewtonSolver::require_residual_reduction, FEMParameters::require_residual_reduction, libMesh::solver, FEMParameters::solver_quiet, libMesh::DifferentiableSystem::time_solver, FEMParameters::verify_analytic_jacobians, and libMesh::FEMSystem::verify_analytic_jacobians.

Referenced by main().

168 {
169  // Use analytical jacobians?
170  system.analytic_jacobians() = param.analytic_jacobians;
171 
172  // Verify analytic jacobians against numerical ones?
174 
175  // Use the prescribed FE type
176  system.fe_family() = param.fe_family[0];
177  system.fe_order() = param.fe_order[0];
178 
179  // More desperate debugging options
181  system.print_solutions = param.print_solutions;
183  system.print_residuals = param.print_residuals;
185  system.print_jacobians = param.print_jacobians;
186 
187  // No transient time solver
188  system.time_solver =
189  UniquePtr<TimeSolver>(new SteadySolver(system));
190 
191  // Nonlinear solver options
192  {
193  NewtonSolver * solver = new NewtonSolver(system);
194  system.time_solver->diff_solver() = UniquePtr<DiffSolver>(solver);
195 
196  solver->quiet = param.solver_quiet;
198  solver->minsteplength = param.min_step_length;
203  if (system.time_solver->reduce_deltat_on_diffsolver_failure)
204  {
205  solver->continue_after_max_iterations = true;
206  solver->continue_after_backtrack_failure = true;
207  }
208 
209  // And the linear solver options
213  }
214 }
bool continue_after_max_iterations
Defaults to true, telling the DiffSolver to continue rather than exit when a solve has reached its ma...
Definition: diff_solver.h:174
bool print_solution_norms
UniquePtr< TimeSolver > time_solver
A pointer to the solver object we&#39;re going to use.
Definition: diff_system.h:221
Real minimum_linear_tolerance
The tolerance for linear solves is kept above this minimum.
Definition: diff_solver.h:215
libMesh::Real initial_linear_tolerance
Real verify_analytic_jacobians
If verify_analytic_jacobian is equal to zero (as it is by default), no numeric jacobians will be calc...
Definition: fem_system.h:203
bool analytic_jacobians
bool quiet
The DiffSolver should not print anything to libMesh::out unless quiet is set to false; default is tru...
Definition: diff_solver.h:162
bool print_jacobian_norms
Set print_jacobian_norms to true to print |J| whenever it is assembled.
Definition: diff_system.h:329
std::string & fe_family()
Definition: L-shaped.h:24
unsigned int max_nonlinear_iterations
The DiffSolver should exit in failure if max_nonlinear_iterations is exceeded and continue_after_max_...
Definition: diff_solver.h:156
bool require_residual_reduction
bool print_residual_norms
libMesh::Real relative_step_tolerance
bool & analytic_jacobians()
Definition: L-shaped.h:26
Real initial_linear_tolerance
Any required linear solves will at first be done with this tolerance; the DiffSolver may tighten the ...
Definition: diff_solver.h:210
bool print_jacobians
Set print_jacobians to true to print J whenever it is assembled.
Definition: diff_system.h:334
This class defines a solver which uses the default libMesh linear solver in a quasiNewton method to h...
Definition: newton_solver.h:46
std::vector< unsigned int > fe_order
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
PetscDiffSolver & solver
bool print_solution_norms
Set print_residual_norms to true to print |U| whenever it is used in an assembly() call...
Definition: diff_system.h:308
unsigned int max_linear_iterations
bool print_residual_norms
Set print_residual_norms to true to print |F| whenever it is assembled.
Definition: diff_system.h:319
unsigned int max_linear_iterations
Each linear solver step should exit after max_linear_iterations is exceeded.
Definition: diff_solver.h:148
libMesh::Real minimum_linear_tolerance
std::vector< std::string > fe_family
libMesh::Real verify_analytic_jacobians
bool print_residuals
Set print_residuals to true to print F whenever it is assembled.
Definition: diff_system.h:324
This class implements a TimeSolver which does a single solve of the steady state problem.
Definition: steady_solver.h:47
Real linear_tolerance_multiplier
The tolerance for linear solves is kept below this multiplier (which defaults to 1e-3) times the norm...
bool print_solutions
Set print_solutions to true to print U whenever it is used in an assembly() call. ...
Definition: diff_system.h:314
libMesh::Real relative_residual_tolerance
unsigned int & fe_order()
Definition: L-shaped.h:25
bool continue_after_backtrack_failure
Defaults to false, telling the DiffSolver to throw an error when the backtracking scheme fails to fin...
Definition: diff_solver.h:180
bool print_jacobian_norms
Real minsteplength
If the quasi-Newton step length must be reduced to below this factor to give a residual reduction...
libMesh::Real min_step_length
unsigned int max_nonlinear_iterations
bool require_residual_reduction
If this is set to true, the solver is forced to test the residual after each Newton step...
Definition: newton_solver.h:98
libMesh::Real linear_tolerance_multiplier
Real relative_residual_tolerance
Definition: diff_solver.h:192
void write_output ( EquationSystems es,
unsigned int  a_step,
std::string  solution_type,
FEMParameters param 
)

Definition at line 105 of file adjoints_ex2.C.

References libMesh::EquationSystems::get_mesh(), libMesh::libmesh_ignore(), mesh, FEMParameters::output_exodus, FEMParameters::output_gmv, libMesh::MeshOutput< MT >::write_equation_systems(), and libMesh::ExodusII_IO::write_timestep().

Referenced by main().

109 {
110  // Ignore parameters when there are no output formats available.
111  libmesh_ignore(es);
112  libmesh_ignore(a_step);
113  libmesh_ignore(solution_type);
114  libmesh_ignore(param);
115 
116 #ifdef LIBMESH_HAVE_GMV
117  if (param.output_gmv)
118  {
119  MeshBase & mesh = es.get_mesh();
120 
121  std::ostringstream file_name_gmv;
122  file_name_gmv << solution_type
123  << ".out.gmv."
124  << std::setw(2)
125  << std::setfill('0')
126  << std::right
127  << a_step;
128 
129  GMVIO(mesh).write_equation_systems(file_name_gmv.str(), es);
130  }
131 #endif
132 
133 #ifdef LIBMESH_HAVE_EXODUS_API
134  if (param.output_exodus)
135  {
136  MeshBase & mesh = es.get_mesh();
137 
138  // We write out one file per adaptive step. The files are named in
139  // the following way:
140  // foo.e
141  // foo.e-s002
142  // foo.e-s003
143  // ...
144  // so that, if you open the first one with Paraview, it actually
145  // opens the entire sequence of adapted files.
146  std::ostringstream file_name_exodus;
147 
148  file_name_exodus << solution_type << ".e";
149  if (a_step > 0)
150  file_name_exodus << "-s"
151  << std::setw(3)
152  << std::setfill('0')
153  << std::right
154  << a_step + 1;
155 
156  // We write each adaptive step as a pseudo "time" step, where the
157  // time simply matches the (1-based) adaptive step we are on.
158  ExodusII_IO(mesh).write_timestep(file_name_exodus.str(),
159  es,
160  1,
161  /*time=*/a_step + 1);
162  }
163 #endif
164 }
void write_timestep(const std::string &fname, const EquationSystems &es, const int timestep, const Real time)
Writes out the solution at a specific timestep.
Definition: exodusII_io.C:773
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
MeshBase & mesh
This class implements writing meshes in the GMV format.
Definition: gmv_io.h:46
This is the MeshBase class.
Definition: mesh_base.h:68
virtual void write_equation_systems(const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
This method implements writing a mesh with data to a specified file where the data is taken from the ...
Definition: mesh_output.C:31
void libmesh_ignore(const T &)
const MeshBase & get_mesh() const