libMesh
adjoints_ex1.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // <h1>Adjoints Example 1 - Laplace Equation in the L-Shaped Domain with Adjoint based mesh refinement</h1>
21 // \author Roy Stogner
22 // \date 2011
23 //
24 // This example solves the Laplace equation on the classic "L-shaped"
25 // domain with adaptive mesh refinement. The exact
26 // solution is u(r,\theta) = r^{2/3} * \sin ( (2/3) * \theta). The kelly and
27 // adjoint residual error estimators are used to develop error indicators and
28 // guide mesh adaptation. Since we use the adjoint capabilities of libMesh in
29 // this example, we use the DiffSystem framework. This file (adjoints_ex1.C)
30 // contains the declaration of mesh and equation system objects, L-shaped.C
31 // contains the assembly of the system, element_qoi_derivative.C and
32 // side_qoi_derivative.C contain the RHS for the adjoint systems.
33 // Postprocessing to compute the QoIs is done in element_postprocess.C and
34 // side_postprocess.C.
35 
36 // The initial mesh contains three QUAD9 elements which represent the
37 // standard quadrants I, II, and III of the domain [-1,1]x[-1,1],
38 // i.e.
39 // Element 0: [-1,0]x[ 0,1]
40 // Element 1: [ 0,1]x[ 0,1]
41 // Element 2: [-1,0]x[-1,0]
42 // The mesh is provided in the standard libMesh ASCII format file
43 // named "lshaped.xda". In addition, an input file named "general.in"
44 // is provided which allows the user to set several parameters for
45 // the solution so that the problem can be re-run without a
46 // re-compile. The solution technique employed is to have a
47 // refinement loop with a linear (forward and adjoint) solve inside followed by a
48 // refinement of the grid and projection of the solution to the new grid
49 // In the final loop iteration, there is no additional
50 // refinement after the solve. In the input file "general.in", the variable
51 // "max_adaptivesteps" controls the number of refinement steps, and
52 // "refine_fraction" / "coarsen_fraction" determine the number of
53 // elements which will be refined / coarsened at each step.
54 
55 // C++ includes
56 #include <iostream>
57 #include <iomanip>
58 
59 // General libMesh includes
60 #include "libmesh/equation_systems.h"
61 #include "libmesh/linear_solver.h"
62 #include "libmesh/error_vector.h"
63 #include "libmesh/mesh.h"
64 #include "libmesh/mesh_refinement.h"
65 #include "libmesh/newton_solver.h"
66 #include "libmesh/numeric_vector.h"
67 #include "libmesh/steady_solver.h"
68 #include "libmesh/system_norm.h"
69 
70 // Error Estimator includes
71 #include "libmesh/kelly_error_estimator.h"
72 #include "libmesh/patch_recovery_error_estimator.h"
73 
74 // Adjoint Related includes
75 #include "libmesh/adjoint_residual_error_estimator.h"
76 #include "libmesh/qoi_set.h"
77 
78 // libMesh I/O includes
79 #include "libmesh/getpot.h"
80 #include "libmesh/gmv_io.h"
81 #include "libmesh/exodusII_io.h"
82 
83 // Local includes
84 #include "femparameters.h"
85 #include "L-shaped.h"
86 
87 // Bring in everything from the libMesh namespace
88 using namespace libMesh;
89 
90 // Local function declarations
91 
92 // Number output files, the files are give a prefix of primal or adjoint_i depending on
93 // whether the output is the primal solution or the dual solution for the ith QoI
94 
95 // Write gmv output
96 
98  unsigned int a_step, // The adaptive step count
99  std::string solution_type, // primal or adjoint solve
100  FEMParameters & param)
101 {
102  // Ignore parameters when there are no output formats available.
103  libmesh_ignore(es);
104  libmesh_ignore(a_step);
105  libmesh_ignore(solution_type);
106  libmesh_ignore(param);
107 
108 #ifdef LIBMESH_HAVE_GMV
109  if (param.output_gmv)
110  {
111  MeshBase & mesh = es.get_mesh();
112 
113  std::ostringstream file_name_gmv;
114  file_name_gmv << solution_type
115  << ".out.gmv."
116  << std::setw(2)
117  << std::setfill('0')
118  << std::right
119  << a_step;
120 
122  (file_name_gmv.str(), es);
123  }
124 #endif
125 
126 #ifdef LIBMESH_HAVE_EXODUS_API
127  if (param.output_exodus)
128  {
129  MeshBase & mesh = es.get_mesh();
130 
131  // We write out one file per adaptive step. The files are named in
132  // the following way:
133  // foo.e
134  // foo.e-s002
135  // foo.e-s003
136  // ...
137  // so that, if you open the first one with Paraview, it actually
138  // opens the entire sequence of adapted files.
139  std::ostringstream file_name_exodus;
140 
141  file_name_exodus << solution_type << ".e";
142  if (a_step > 0)
143  file_name_exodus << "-s"
144  << std::setw(3)
145  << std::setfill('0')
146  << std::right
147  << a_step + 1;
148 
149  // We write each adaptive step as a pseudo "time" step, where the
150  // time simply matches the (1-based) adaptive step we are on.
151  ExodusII_IO(mesh).write_timestep(file_name_exodus.str(),
152  es,
153  1,
154  /*time=*/a_step + 1);
155  }
156 #endif
157 }
158 
159 // Set the parameters for the nonlinear and linear solvers to be used during the simulation
160 
162  FEMParameters & param)
163 {
164  // Use analytical jacobians?
165  system.analytic_jacobians() = param.analytic_jacobians;
166 
167  // Verify analytic jacobians against numerical ones?
169 
170  // Use the prescribed FE type
171  system.fe_family() = param.fe_family[0];
172  system.fe_order() = param.fe_order[0];
173 
174  // More desperate debugging options
176  system.print_solutions = param.print_solutions;
178  system.print_residuals = param.print_residuals;
180  system.print_jacobians = param.print_jacobians;
181 
182  // No transient time solver
183  system.time_solver =
184  UniquePtr<TimeSolver>(new SteadySolver(system));
185 
186  // Nonlinear solver options
187  {
188  NewtonSolver * solver = new NewtonSolver(system);
189  system.time_solver->diff_solver() = UniquePtr<DiffSolver>(solver);
190 
191  solver->quiet = param.solver_quiet;
192  solver->verbose = param.solver_verbose;
194  solver->minsteplength = param.min_step_length;
199  if (system.time_solver->reduce_deltat_on_diffsolver_failure)
200  {
201  solver->continue_after_max_iterations = true;
202  solver->continue_after_backtrack_failure = true;
203  }
204 
205  // And the linear solver options
209  }
210 }
211 
212 // Build the mesh refinement object and set parameters for refining/coarsening etc
213 
214 #ifdef LIBMESH_ENABLE_AMR
215 
217  FEMParameters & param)
218 {
219  MeshRefinement * mesh_refinement = new MeshRefinement(mesh);
220  mesh_refinement->coarsen_by_parents() = true;
221  mesh_refinement->absolute_global_tolerance() = param.global_tolerance;
222  mesh_refinement->nelem_target() = param.nelem_target;
223  mesh_refinement->refine_fraction() = param.refine_fraction;
224  mesh_refinement->coarsen_fraction() = param.coarsen_fraction;
225  mesh_refinement->coarsen_threshold() = param.coarsen_threshold;
226 
227  return UniquePtr<MeshRefinement>(mesh_refinement);
228 }
229 
230 #endif // LIBMESH_ENABLE_AMR
231 
232 // This is where we declare the error estimators to be built and used for
233 // mesh refinement. The adjoint residual estimator needs two estimators.
234 // One for the forward component of the estimate and one for the adjoint
235 // weighting factor. Here we use the Patch Recovery indicator to estimate both the
236 // forward and adjoint weights. The H1 seminorm component of the error is used
237 // as dictated by the weak form the Laplace equation.
238 
240  QoISet & qois)
241 {
242  if (param.indicator_type == "kelly")
243  {
244  libMesh::out << "Using Kelly Error Estimator" << std::endl;
245 
247  }
248  else if (param.indicator_type == "adjoint_residual")
249  {
250  libMesh::out << "Using Adjoint Residual Error Estimator with Patch Recovery Weights" << std::endl;
251 
252  AdjointResidualErrorEstimator * adjoint_residual_estimator = new AdjointResidualErrorEstimator;
253 
254  adjoint_residual_estimator->qoi_set() = qois;
255 
256  adjoint_residual_estimator->error_plot_suffix = "error.gmv";
257 
259  adjoint_residual_estimator->primal_error_estimator().reset(p1);
260 
262  adjoint_residual_estimator->dual_error_estimator().reset(p2);
263 
264  adjoint_residual_estimator->primal_error_estimator()->error_norm.set_type(0, H1_SEMINORM);
265  p1->set_patch_reuse(param.patch_reuse);
266 
267  adjoint_residual_estimator->dual_error_estimator()->error_norm.set_type(0, H1_SEMINORM);
268  p2->set_patch_reuse(param.patch_reuse);
269 
270  return UniquePtr<ErrorEstimator>(adjoint_residual_estimator);
271  }
272  else
273  libmesh_error_msg("Unknown indicator_type = " << param.indicator_type);
274 }
275 
276 // The main program.
277 int main (int argc, char ** argv)
278 {
279  // Initialize libMesh.
280  LibMeshInit init (argc, argv);
281 
282  // This example requires a linear solver package.
283  libmesh_example_requires(libMesh::default_solver_package() != INVALID_SOLVER_PACKAGE,
284  "--enable-petsc, --enable-trilinos, or --enable-eigen");
285 
286  // Skip adaptive examples on a non-adaptive libMesh build
287 #ifndef LIBMESH_ENABLE_AMR
288  libmesh_example_requires(false, "--enable-amr");
289 #else
290 
291  libMesh::out << "Started " << argv[0] << std::endl;
292 
293  // Make sure the general input file exists, and parse it
294  {
295  std::ifstream i("general.in");
296  if (!i)
297  libmesh_error_msg('[' << init.comm().rank() << "] Can't find general.in; exiting early.");
298  }
299  GetPot infile("general.in");
300 
301  // Read in parameters from the input file
302  FEMParameters param(init.comm());
303  param.read(infile);
304 
305  // Skip this default-2D example if libMesh was compiled as 1D-only.
306  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");
307 
308  // Create a mesh, with dimension to be overridden later, distributed
309  // across the default MPI communicator.
310  Mesh mesh(init.comm());
311 
312  // And an object to refine it
313  UniquePtr<MeshRefinement> mesh_refinement =
314  build_mesh_refinement(mesh, param);
315 
316  // And an EquationSystems to run on it
317  EquationSystems equation_systems (mesh);
318 
319  libMesh::out << "Reading in and building the mesh" << std::endl;
320 
321  // Read in the mesh
322  mesh.read(param.domainfile.c_str());
323  // Make all the elements of the mesh second order so we can compute
324  // with a higher order basis
326 
327  // Create a mesh refinement object to do the initial uniform refinements
328  // on the coarse grid read in from lshaped.xda
329  MeshRefinement initial_uniform_refinements(mesh);
330  initial_uniform_refinements.uniformly_refine(param.coarserefinements);
331 
332  libMesh::out << "Building system" << std::endl;
333 
334  // Build the FEMSystem
335  LaplaceSystem & system = equation_systems.add_system<LaplaceSystem> ("LaplaceSystem");
336 
337  // Set its parameters
338  set_system_parameters(system, param);
339 
340  libMesh::out << "Initializing systems" << std::endl;
341 
342  equation_systems.init ();
343 
344  // Print information about the mesh and system to the screen.
345  mesh.print_info();
346  equation_systems.print_info();
347  LinearSolver<Number> * linear_solver = system.get_linear_solver();
348 
349  {
350  // Adaptively solve the timestep
351  unsigned int a_step = 0;
352  for (; a_step != param.max_adaptivesteps; ++a_step)
353  {
354  // We can't adapt to both a tolerance and a
355  // target mesh size
356  if (param.global_tolerance != 0.)
357  libmesh_assert_equal_to (param.nelem_target, 0);
358  // If we aren't adapting to a tolerance we need a
359  // target mesh size
360  else
361  libmesh_assert_greater (param.nelem_target, 0);
362 
363  linear_solver->reuse_preconditioner(false);
364 
365  // Solve the forward problem
366  system.solve();
367 
368  // Write out the computed primal solution
369  write_output(equation_systems, a_step, "primal", param);
370 
371  // Get a pointer to the primal solution vector
372  NumericVector<Number> & primal_solution = *system.solution;
373 
374  // Declare a QoISet object, we need this object to set weights for our QoI error contributions
375  QoISet qois;
376 
377  // Declare a qoi_indices vector, each index will correspond to a QoI
378  std::vector<unsigned int> qoi_indices;
379  qoi_indices.push_back(0);
380  qoi_indices.push_back(1);
381  qois.add_indices(qoi_indices);
382 
383  // Set weights for each index, these will weight the contribution of each QoI in the final error
384  // estimate to be used for flagging elements for refinement
385  qois.set_weight(0, 0.5);
386  qois.set_weight(1, 0.5);
387 
388  // Make sure we get the contributions to the adjoint RHS from the sides
389  system.assemble_qoi_sides = true;
390 
391  // We are about to solve the adjoint system, but before we do this we see the same preconditioner
392  // flag to reuse the preconditioner from the forward solver
393  linear_solver->reuse_preconditioner(param.reuse_preconditioner);
394 
395  // Solve the adjoint system. This takes the transpose of the stiffness matrix and then
396  // solves the resulting system
397  system.adjoint_solve();
398 
399  // Now that we have solved the adjoint, set the adjoint_already_solved boolean to true, so we dont solve unnecessarily in the error estimator
400  system.set_adjoint_already_solved(true);
401 
402  // Get a pointer to the solution vector of the adjoint problem for QoI 0
403  NumericVector<Number> & dual_solution_0 = system.get_adjoint_solution(0);
404 
405  // Swap the primal and dual solutions so we can write out the adjoint solution
406  primal_solution.swap(dual_solution_0);
407  write_output(equation_systems, a_step, "adjoint_0", param);
408 
409  // Swap back
410  primal_solution.swap(dual_solution_0);
411 
412  // Get a pointer to the solution vector of the adjoint problem for QoI 0
413  NumericVector<Number> & dual_solution_1 = system.get_adjoint_solution(1);
414 
415  // Swap again
416  primal_solution.swap(dual_solution_1);
417  write_output(equation_systems, a_step, "adjoint_1", param);
418 
419  // Swap back again
420  primal_solution.swap(dual_solution_1);
421 
422  libMesh::out << "Adaptive step "
423  << a_step
424  << ", we have "
425  << mesh.n_active_elem()
426  << " active elements and "
427  << equation_systems.n_active_dofs()
428  << " active dofs."
429  << std::endl;
430 
431  // Postprocess, compute the approximate QoIs and write them out to the console
432  libMesh::out << "Postprocessing: " << std::endl;
433  system.postprocess_sides = true;
434  system.postprocess();
435  Number QoI_0_computed = system.get_QoI_value("computed", 0);
436  Number QoI_0_exact = system.get_QoI_value("exact", 0);
437  Number QoI_1_computed = system.get_QoI_value("computed", 1);
438  Number QoI_1_exact = system.get_QoI_value("exact", 1);
439 
440  libMesh::out << "The relative error in QoI 0 is "
441  << std::setprecision(17)
442  << std::abs(QoI_0_computed - QoI_0_exact) / std::abs(QoI_0_exact)
443  << std::endl;
444 
445  libMesh::out << "The relative error in QoI 1 is "
446  << std::setprecision(17)
447  << std::abs(QoI_1_computed - QoI_1_exact) / std::abs(QoI_1_exact)
448  << std::endl
449  << std::endl;
450 
451  // Now we construct the data structures for the mesh refinement process
452  ErrorVector error;
453 
454  // Build an error estimator object
455  UniquePtr<ErrorEstimator> error_estimator =
456  build_error_estimator(param, qois);
457 
458  // Estimate the error in each element using the Adjoint Residual or Kelly error estimator
459  error_estimator->estimate_error(system, error);
460 
461  // We have to refine either based on reaching an error tolerance or
462  // a number of elements target, which should be verified above
463  // Otherwise we flag elements by error tolerance or nelem target
464 
465  // Uniform refinement
466  if (param.refine_uniformly)
467  {
468  mesh_refinement->uniformly_refine(1);
469  }
470  // Adaptively refine based on reaching an error tolerance
471  else if (param.global_tolerance >= 0. && param.nelem_target == 0.)
472  {
473  mesh_refinement->flag_elements_by_error_tolerance (error);
474 
475  mesh_refinement->refine_and_coarsen_elements();
476  }
477  // Adaptively refine based on reaching a target number of elements
478  else
479  {
480  if (mesh.n_active_elem() >= param.nelem_target)
481  {
482  libMesh::out << "We reached the target number of elements." << std::endl << std::endl;
483  break;
484  }
485 
486  mesh_refinement->flag_elements_by_nelem_target (error);
487 
488  mesh_refinement->refine_and_coarsen_elements();
489  }
490 
491  // Dont forget to reinit the system after each adaptive refinement !
492  equation_systems.reinit();
493 
494  libMesh::out << "Refined mesh to "
495  << mesh.n_active_elem()
496  << " active elements and "
497  << equation_systems.n_active_dofs()
498  << " active dofs."
499  << std::endl;
500  }
501 
502  // Do one last solve if necessary
503  if (a_step == param.max_adaptivesteps)
504  {
505  linear_solver->reuse_preconditioner(false);
506  system.solve();
507 
508  write_output(equation_systems, a_step, "primal", param);
509 
510  NumericVector<Number> & primal_solution = *system.solution;
511 
512  QoISet qois;
513  std::vector<unsigned int> qoi_indices;
514 
515  qoi_indices.push_back(0);
516  qoi_indices.push_back(1);
517  qois.add_indices(qoi_indices);
518 
519  qois.set_weight(0, 0.5);
520  qois.set_weight(1, 0.5);
521 
522  system.assemble_qoi_sides = true;
523  linear_solver->reuse_preconditioner(param.reuse_preconditioner);
524  system.adjoint_solve();
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  NumericVector<Number> & dual_solution_0 = system.get_adjoint_solution(0);
530 
531  primal_solution.swap(dual_solution_0);
532  write_output(equation_systems, a_step, "adjoint_0", param);
533 
534  primal_solution.swap(dual_solution_0);
535 
536  NumericVector<Number> & dual_solution_1 = system.get_adjoint_solution(1);
537 
538  primal_solution.swap(dual_solution_1);
539  write_output(equation_systems, a_step, "adjoint_1", param);
540 
541  primal_solution.swap(dual_solution_1);
542 
543  libMesh::out << "Adaptive step "
544  << a_step
545  << ", we have "
546  << mesh.n_active_elem()
547  << " active elements and "
548  << equation_systems.n_active_dofs()
549  << " active dofs."
550  << std::endl;
551 
552  libMesh::out << "Postprocessing: " << std::endl;
553  system.postprocess_sides = true;
554  system.postprocess();
555 
556  Number QoI_0_computed = system.get_QoI_value("computed", 0);
557  Number QoI_0_exact = system.get_QoI_value("exact", 0);
558  Number QoI_1_computed = system.get_QoI_value("computed", 1);
559  Number QoI_1_exact = system.get_QoI_value("exact", 1);
560 
561  libMesh::out << "The relative error in QoI 0 is "
562  << std::setprecision(17)
563  << std::abs(QoI_0_computed - QoI_0_exact) / std::abs(QoI_0_exact)
564  << std::endl;
565 
566  libMesh::out << "The relative error in QoI 1 is "
567  << std::setprecision(17)
568  << std::abs(QoI_1_computed - QoI_1_exact) / std::abs(QoI_1_exact)
569  << std::endl
570  << std::endl;
571 
572  // Hard coded asserts to ensure that the actual numbers we are getting are what they should be
573  libmesh_assert_less(std::abs(QoI_0_computed - QoI_0_exact)/std::abs(QoI_0_exact), 4.e-5);
574  libmesh_assert_less(std::abs(QoI_1_computed - QoI_1_exact)/std::abs(QoI_1_exact), 1.e-4);
575  }
576  }
577 
578  libMesh::err << '[' << mesh.processor_id()
579  << "] Completing output."
580  << std::endl;
581 
582 #endif // #ifndef LIBMESH_ENABLE_AMR
583 
584  // All done.
585  return 0;
586 }
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
unsigned int nelem_target
Definition: femparameters.h:56
OStreamProxy err
bool print_solution_norms
This class implements a goal oriented error indicator, by weighting residual-based estimates from the...
double abs(double a)
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
This is the EquationSystems class.
libMesh::Real initial_linear_tolerance
virtual dof_id_type n_active_elem() const =0
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
void write_output(EquationSystems &es, unsigned int a_step, std::string solution_type, FEMParameters &param)
Definition: adjoints_ex1.C:97
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
void set_weight(unsigned int, Real)
Set the weight for this index.
Definition: qoi_set.h:229
void set_adjoint_already_solved(bool setting)
Setter for the adjoint_already_solved boolean.
Definition: system.h:378
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
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
bool print_jacobian_norms
Set print_jacobian_norms to true to print |J| whenever it is assembled.
Definition: diff_system.h:329
Real & absolute_global_tolerance()
If absolute_global_tolerance is set to a nonzero value, methods like flag_elements_by_global_toleranc...
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices=QoISet()) libmesh_override
This function sets the _is_adjoint boolean member of TimeSolver to true and then calls the adjoint_so...
Definition: diff_system.C:164
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Definition: exodusII_io.h:52
Number & get_QoI_value(std::string type, unsigned int QoI_index)
Definition: L-shaped.h:32
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
QoISet & qoi_set()
Access to the QoISet (default: weight all QoIs equally) to use when computing errors.
Real & coarsen_fraction()
The coarsen_fraction sets either a desired target or a desired maximum number of elements to flag for...
MeshBase & mesh
bool print_residual_norms
Real & refine_fraction()
The refine_fraction sets either a desired target or a desired maximum number of elements to flag for ...
libMesh::Real relative_step_tolerance
std::string indicator_type
bool & analytic_jacobians()
Definition: L-shaped.h:26
bool postprocess_sides
If postprocess_sides is true (it is false by default), the postprocessing loop will loop over all sid...
Definition: diff_system.h:302
UniquePtr< ErrorEstimator > & primal_error_estimator()
Access to the "subestimator" (default: PatchRecovery) to use on the primal/forward solution...
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
The libMesh namespace provides an interface to certain functionality in the library.
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
This class implements writing meshes in the GMV format.
Definition: gmv_io.h:46
libMesh::Real refine_fraction
Definition: femparameters.h:58
void set_system_parameters(LaplaceSystem &system, FEMParameters &param)
Definition: adjoints_ex1.C:161
virtual LinearSolver< Number > * get_linear_solver() const libmesh_override
Definition: diff_system.C:175
bool print_jacobians
Set print_jacobians to true to print J whenever it is assembled.
Definition: diff_system.h:334
bool & coarsen_by_parents()
If coarsen_by_parents is true, complete groups of sibling elements (elements with the same parent) wi...
This class defines a solver which uses the default libMesh linear solver in a quasiNewton method to h...
Definition: newton_solver.h:46
This is the MeshBase class.
Definition: mesh_base.h:68
std::vector< unsigned int > fe_order
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
PetscDiffSolver & solver
SolverPackage default_solver_package()
Definition: libmesh.C:995
std::string error_plot_suffix
To aid in investigating error estimator behavior, set this string to a suffix with which to plot (pre...
libMesh::Real coarsen_threshold
Definition: femparameters.h:58
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
virtual void solve() libmesh_override
Invokes the solver associated with the system.
Definition: fem_system.C:1056
int main(int argc, char **argv)
Definition: adjoints_ex1.C:277
This is the MeshRefinement class.
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
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 print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
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
virtual void reinit()
Reinitialize all the systems.
std::vector< std::string > fe_family
libMesh::Real verify_analytic_jacobians
Real & coarsen_threshold()
The coarsen_threshold provides hysteresis in AMR/C strategies.
bool print_residuals
Set print_residuals to true to print F whenever it is assembled.
Definition: diff_system.h:324
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
virtual System & add_system(const std::string &system_type, const std::string &name)
Add the system of type system_type named name to the systems array.
This class implements a TimeSolver which does a single solve of the steady state problem.
Definition: steady_solver.h:47
const Parallel::Communicator & comm() const
Definition: libmesh.h:81
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
libMesh::Real global_tolerance
Definition: femparameters.h:57
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
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 &)
unsigned int & fe_order()
Definition: L-shaped.h:25
This class implements the Kelly error indicator which is based on the flux jumps between elements...
bool verbose
The DiffSolver may print a lot more to libMesh::out if verbose is set to true; default is false...
Definition: diff_solver.h:168
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
virtual void reuse_preconditioner(bool)
Set the same_preconditioner flag, which indicates if we reuse the same preconditioner for subsequent ...
bool print_jacobian_norms
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
Real minsteplength
If the quasi-Newton step length must be reduced to below this factor to give a residual reduction...
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
This class implements the Patch Recovery error indicator.
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
libMesh::Real min_step_length
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
UniquePtr< ErrorEstimator > build_error_estimator(FEMParameters &param, QoISet &qois)
Definition: adjoints_ex1.C:239
const MeshBase & get_mesh() const
unsigned int max_nonlinear_iterations
unsigned int rank() const
Definition: parallel.h:724
UniquePtr< MeshRefinement > build_mesh_refinement(MeshBase &mesh, FEMParameters &param)
Definition: adjoints_ex1.C:216
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
virtual void init()
Initialize all the systems.
std::size_t n_active_dofs() const
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
UniquePtr< ErrorEstimator > & dual_error_estimator()
Access to the "subestimator" (default: PatchRecovery) to use on the dual/adjoint solution.
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.
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448
libMesh::Real linear_tolerance_multiplier
virtual void postprocess()
Runs a postprocessing loop over all elements, and if postprocess_sides is true over all sides...
Definition: L-shaped.C:168
Real relative_residual_tolerance
Definition: diff_solver.h:192
libMesh::Real coarsen_fraction
Definition: femparameters.h:58
processor_id_type processor_id() const
void uniformly_refine(unsigned int n=1)
Uniformly refines the mesh n times.