AMDiS  0.1
The Adaptive Multi-Dimensional Simulation Toolbox
ProblemStat.inc.hpp
1 #pragma once
2 
3 #include <map>
4 #include <string>
5 #include <utility>
6 
7 #include <dune/common/hybridutilities.hh>
8 #include <dune/common/timer.hh>
9 #include <dune/functions/functionspacebases/subspacebasis.hh>
10 #include <dune/grid/common/capabilities.hh>
11 #include <dune/typetree/childextraction.hh>
12 
13 #include <amdis/AdaptInfo.hpp>
14 #include <amdis/BackupRestore.hpp>
15 #include <amdis/Assembler.hpp>
16 #include <amdis/GridFunctionOperator.hpp>
17 #include <amdis/io/FileWriterCreator.hpp>
18 #include <amdis/linearalgebra/SymmetryStructure.hpp>
19 
20 namespace AMDiS {
21 
22 template <class Traits>
24  Flag initFlag,
25  Self* adoptProblem,
26  Flag adoptFlag)
27 {
28  // create grids
29  if (!grid_) {
30  if (initFlag.isSet(CREATE_MESH) ||
31  (!adoptFlag.isSet(INIT_MESH) &&
32  (initFlag.isSet(INIT_SYSTEM) || initFlag.isSet(INIT_FE_SPACE)))) {
33  createGrid();
34  }
35 
36  if (adoptProblem &&
37  (adoptFlag.isSet(INIT_MESH) ||
38  adoptFlag.isSet(INIT_SYSTEM) ||
39  adoptFlag.isSet(INIT_FE_SPACE))) {
40  adoptGrid(adoptProblem->grid_, adoptProblem->boundaryManager_);
41  }
42  }
43 
44  if (!grid_)
45  warning("no grid created");
46 
47  // create fespace
48  if (!globalBasis_) {
49  if (initFlag.isSet(INIT_FE_SPACE) ||
50  (initFlag.isSet(INIT_SYSTEM) && !adoptFlag.isSet(INIT_FE_SPACE))) {
51  createGlobalBasis();
52  }
53 
54  if (adoptProblem &&
55  (adoptFlag.isSet(INIT_FE_SPACE) || adoptFlag.isSet(INIT_SYSTEM))) {
56  adoptGlobalBasis(adoptProblem->globalBasis_);
57  }
58  }
59 
60  if (!globalBasis_)
61  warning("no globalBasis created\n");
62 
63  // create system
64  if (initFlag.isSet(INIT_SYSTEM))
65  createMatricesAndVectors();
66 
67  if (adoptProblem && adoptFlag.isSet(INIT_SYSTEM)) {
68  systemMatrix_ = adoptProblem->systemMatrix_;
69  solution_ = adoptProblem->solution_;
70  rhs_ = adoptProblem->rhs_;
71  estimates_ = adoptProblem->estimates_;
72  }
73 
74 
75  // create solver
76  if (!linearSolver_) {
77  if (initFlag.isSet(INIT_SOLVER))
78  createSolver();
79 
80  if (adoptProblem && adoptFlag.isSet(INIT_SOLVER)) {
81  test_exit(!linearSolver_, "solver already created\n");
82  linearSolver_ = adoptProblem->linearSolver_;
83  }
84  }
85 
86  if (!linearSolver_) {
87  warning("no solver created\n");
88  }
89 
90  // create marker
91  if (initFlag.isSet(INIT_MARKER))
92  createMarker();
93 
94  if (adoptProblem && adoptFlag.isSet(INIT_MARKER))
95  marker_ = adoptProblem->marker_;
96 
97 
98  // create file writer
99  if (initFlag.isSet(INIT_FILEWRITER))
100  createFileWriter();
101 
102  solution_->resizeZero();
103 }
104 
105 
106 template <class Traits>
108 restore(Flag initFlag)
109 {
110  std::string grid_filename = Parameters::get<std::string>(name_ + "->restore->grid").value();
111  std::string solution_filename = Parameters::get<std::string>(name_ + "->restore->solution").value();
112  test_exit(filesystem::exists(grid_filename), "Restore file '{}' not found.", grid_filename);
113  test_exit(filesystem::exists(solution_filename), "Restore file '{}' not found.", solution_filename);
114 
115  // TODO(SP): implement BAckupRestore independent of wrapped grid
116  using HostGrid = typename Grid::HostGrid;
117 
118  // restore grid from file
119  if (Dune::Capabilities::hasBackupRestoreFacilities<HostGrid>::v)
120  adoptGrid(std::shared_ptr<HostGrid>(Dune::BackupRestoreFacility<HostGrid>::restore(grid_filename)));
121  else
122  adoptGrid(std::shared_ptr<HostGrid>(BackupRestoreByGridFactory<HostGrid>::restore(grid_filename)));
123 
124  // create fespace
125  if (initFlag.isSet(INIT_FE_SPACE) || initFlag.isSet(INIT_SYSTEM))
126  createGlobalBasis();
127 
128  // create system
129  if (initFlag.isSet(INIT_SYSTEM))
130  createMatricesAndVectors();
131 
132  // create solver
133  if (!linearSolver_ && initFlag.isSet(INIT_SOLVER))
134  createSolver();
135 
136  // create marker
137  if (initFlag.isSet(INIT_MARKER))
138  createMarker();
139 
140  // create file writer
141  if (initFlag.isSet(INIT_FILEWRITER))
142  createFileWriter();
143 
144  solution_->resize(sizeInfo(*globalBasis_));
145  solution_->restore(solution_filename);
146 }
147 
148 
149 template <class Traits>
151 {
152  Parameters::get(name_ + "->mesh", gridName_);
153 
154  MeshCreator<Grid> creator(gridName_);
155  grid_ = creator.create();
156 
157  boundaryManager_ = std::make_shared<BoundaryManager<Grid>>(grid_);
158  if (!creator.boundaryIds().empty())
159  boundaryManager_->setBoundaryIds(creator.boundaryIds());
160 
161  info(3,"Create grid:");
162  info(3,"#elements = {}" , grid_->size(0));
163  info(3,"#faces/edges = {}", grid_->size(1));
164  info(3,"#vertices = {}" , grid_->size(dim));
165  info(3,"overlap-size = {}", grid_->leafGridView().overlapSize(0));
166  info(3,"ghost-size = {}" , grid_->leafGridView().ghostSize(0));
167  info(3,"");
168 }
169 
170 
171 template <class T, class GV>
172 using HasCreate = decltype(T::create(std::declval<GV>()));
173 
174 
175 template <class Traits>
177 {
178  createGlobalBasisImpl(Dune::Std::is_detected<HasCreate,Traits,GridView>{});
179  initGlobalBasis();
180 }
181 
182 
183 template <class Traits>
185 {
186  assert( bool(grid_) );
187  static_assert(std::is_same_v<GridView, typename Grid::LeafGridView>, "");
188  auto basis = Traits::create(name_, grid_->leafGridView());
189  globalBasis_ = std::make_shared<GlobalBasis>(std::move(basis));
190 }
191 
192 
193 template <class Traits>
194 void ProblemStat<Traits>::createGlobalBasisImpl(std::false_type)
195 {
196  error_exit("Cannot create GlobalBasis from type. Pass a BasisCreator instead!");
197 }
198 
199 
200 template <class Traits>
202 
203 
204 template <class Traits>
206 {
207  systemMatrix_ = std::make_shared<SystemMatrix>(globalBasis_, globalBasis_);
208  std::string symmetryStr = "unknown";
209  Parameters::get(name_ + "->symmetry", symmetryStr);
210  systemMatrix_->setSymmetryStructure(symmetryStr);
211 
212  solution_ = std::make_shared<SolutionVector>(globalBasis_);
213  rhs_ = std::make_shared<SystemVector>(globalBasis_);
214 
215  auto localView = globalBasis_->localView();
216  for_each_node(localView.tree(), [&,this](auto&&, auto treePath) -> void
217  {
218  std::string i = to_string(treePath);
219  estimates_[i].resize(globalBasis_->gridView().indexSet().size(0));
220  for (std::size_t j = 0; j < estimates_[i].size(); j++)
221  estimates_[i][j] = 0.0; // TODO: Remove when estimate() is implemented
222  });
223 }
224 
225 
226 template <class Traits>
228 {
229  std::string solverName = "default";
230  Parameters::get(name_ + "->solver", solverName);
231 
232  auto solverCreator
233  = named(CreatorMap<LinearSolver>::getCreator(solverName, name_ + "->solver"));
234 
235  linearSolver_ = solverCreator->createWithString(name_ + "->solver");
236 }
237 
238 
239 template <class Traits>
241 {
242  marker_.clear();
243  auto localView = globalBasis_->localView();
244  for_each_node(localView.tree(), [&,this](auto&&, auto treePath) -> void
245  {
246  std::string componentName = name_ + "->marker[" + to_string(treePath) + "]";
247 
248  if (!Parameters::get<std::string>(componentName + "->strategy"))
249  return;
250 
251  std::string tp = to_string(treePath);
252  auto newMarker
253  = EstimatorMarker<Grid>::createMarker(componentName, tp, estimates_[tp], grid_);
254  assert(bool(newMarker));
255  this->addMarker(std::move(newMarker));
256  });
257 }
258 
259 
260 template <class Traits>
262 {
263  FileWriterCreator<SolutionVector> creator(solution_, boundaryManager_);
264 
265  filewriter_.clear();
266  auto localView = globalBasis_->localView();
267  for_each_node(localView.tree(), [&](auto const& /*node*/, auto treePath) -> void
268  {
269  std::string componentName = name_ + "->output[" + to_string(treePath) + "]";
270  auto format = Parameters::get<std::vector<std::string>>(componentName + "->format");
271 
272  if (!format && to_string(treePath).empty()) {
273  // alternative for root treepath
274  componentName = name_ + "->output";
275  format = Parameters::get<std::vector<std::string>>(componentName + "->format");
276  }
277 
278  if (!format)
279  return;
280 
281  for (std::string const& type : format.value()) {
282  auto writer = creator.create(type, componentName, treePath);
283  if (writer)
284  filewriter_.push_back(std::move(writer));
285  }
286  });
287 }
288 
289 
290 // Adds a Dirichlet boundary condition
291 template <class Traits>
292  template <class Predicate, class RowTreePath, class ColTreePath, class Values>
294 addDirichletBC(Predicate const& predicate, RowTreePath row, ColTreePath col, Values const& values)
295 {
296  static constexpr bool isValidPredicate = Concepts::Functor<Predicate, bool(WorldVector)>;
297  static_assert( Concepts::Functor<Predicate, bool(WorldVector)>,
298  "Function passed to addDirichletBC for `predicate` does not model the Functor<bool(WorldVector)> concept");
299 
300  static constexpr bool isValidTreePath =
301  Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, RowTreePath> &&
302  Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, ColTreePath>;
303  static_assert(isValidTreePath, "Invalid row and/or col treepath passed to addDirichletBC!");
304 
305  if constexpr (isValidPredicate && isValidTreePath) {
306  auto localView = globalBasis_->localView();
307  auto i = makeTreePath(row);
308  auto j = makeTreePath(col);
309  auto rowBasis = Dune::Functions::subspaceBasis(*globalBasis_, i);
310  auto colBasis = Dune::Functions::subspaceBasis(*globalBasis_, j);
311 
312  auto valueGridFct = makeGridFunction(values, this->gridView());
313 
314  auto bc = makeDirichletBC<SystemMatrix, SolutionVector, SystemVector>(
315  std::move(rowBasis), std::move(colBasis), {predicate}, valueGridFct);
316  boundaryConditions_[i][j].push_back(makeUniquePtr(std::move(bc)));
317  }
318 }
319 
320 
321 // Adds a Dirichlet boundary condition
322 template <class Traits>
323  template <class RowTreePath, class ColTreePath, class Values>
325 addDirichletBC(BoundaryType id, RowTreePath row, ColTreePath col, Values const& values)
326 {
327  static constexpr bool isValidTreePath =
328  Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, RowTreePath> &&
329  Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, ColTreePath>;
330  static_assert(isValidTreePath, "Invalid row and/or col treepath passed to addDirichletBC!");
331 
332  if constexpr (isValidTreePath) {
333  auto localView = globalBasis_->localView();
334  auto i = makeTreePath(row);
335  auto j = makeTreePath(col);
336  auto rowBasis = Dune::Functions::subspaceBasis(*globalBasis_, i);
337  auto colBasis = Dune::Functions::subspaceBasis(*globalBasis_, j);
338 
339  auto valueGridFct = makeGridFunction(values, this->gridView());
340 
341  auto bc = makeDirichletBC<SystemMatrix, SolutionVector, SystemVector>(
342  std::move(rowBasis), std::move(colBasis), {*boundaryManager_, id}, valueGridFct);
343  boundaryConditions_[i][j].push_back(makeUniquePtr(std::move(bc)));
344  }
345 }
346 
347 
348 template <class Traits>
350 addPeriodicBC(BoundaryType id, WorldMatrix const& matrix, WorldVector const& vector)
351 {
352  auto localView = globalBasis_->localView();
353  auto basis = Dune::Functions::subspaceBasis(*globalBasis_, makeTreePath());
354  auto bc = makePeriodicBC<SystemMatrix, SolutionVector, SystemVector>(
355  std::move(basis), {*boundaryManager_, id}, {matrix, vector});
356  boundaryConditions_[makeTreePath()][makeTreePath()].push_back(makeUniquePtr(std::move(bc)));
357 }
358 
359 
360 template <class Traits>
362 solve(AdaptInfo& /*adaptInfo*/, bool createMatrixData, bool storeMatrixData)
363 {
364  Dune::Timer t;
365 
366  SolverInfo solverInfo(name_ + "->solver");
367  solverInfo.setCreateMatrixData(createMatrixData);
368  solverInfo.setStoreMatrixData(storeMatrixData);
369 
370  solution_->resize();
371  linearSolver_->solve(*systemMatrix_, *solution_, *rhs_, solverInfo);
372 
373  if (solverInfo.info() > 0) {
374  msg("solution of discrete system needed {} seconds", t.elapsed());
375 
376  if (solverInfo.absResidual() >= 0.0) {
377  if (solverInfo.relResidual() >= 0.0)
378  msg("Residual norm: ||b-Ax|| = {}, ||b-Ax||/||b|| = {}",
379  solverInfo.absResidual(), solverInfo.relResidual());
380  else
381  msg("Residual norm: ||b-Ax|| = {}", solverInfo.absResidual());
382  }
383  }
384 
385  test_exit(!solverInfo.doBreak() && !solverInfo.error(), "Could not solver the linear system!");
386 }
387 
388 
389 template <class Traits>
392 {
393  Dune::Timer t;
394 
395  Flag markFlag = 0;
396  for (auto& currentMarker : marker_)
397  markFlag |= currentMarker.second->markGrid(adaptInfo);
398 
399  msg("markElements needed {} seconds", t.elapsed());
400 
401  return markFlag;
402 }
403 
404 
405 template <class Traits>
408 {
409  Dune::Timer t;
410  bool adapted = false;
411  // TODO(FM): Find a less expensive alternative to the loop adaption
412  for (int i = 0; i < n; ++i) {
413  // mark all entities for coarsening
414  for (const auto& element : elements(grid_->leafGridView()))
415  grid_->mark(-1, element);
416 
417  bool adaptedInLoop = grid_->preAdapt();
418  adaptedInLoop |= grid_->adapt();
419  grid_->postAdapt();
420  if (!adaptedInLoop)
421  break;
422  else
423  adapted = true;
424  }
425 
426  msg("globalCoarsen needed {} seconds", t.elapsed());
427  return adapted ? MESH_ADAPTED : Flag(0);
428 }
429 
430 
431 // grid has globalRefine(int, AdaptDataHandleInterface&)
432 template <class G>
433 using HasGlobalRefineADHI = decltype(
434  std::declval<G>().globalRefine(1,std::declval<typename G::ADHI&>()));
435 
436 template <class Traits>
439 {
440  Dune::Timer t;
441  if constexpr (Dune::Std::is_detected<HasGlobalRefineADHI, Grid>::value)
442  grid_->globalRefine(n, globalBasis_->globalRefineCallback());
443  else
444  grid_->globalRefine(n);
445 
446  msg("globalRefine needed {} seconds", t.elapsed());
447  return n > 0 ? MESH_ADAPTED : Flag(0);
448 }
449 
450 
451 template <class Traits>
453 adaptGrid(AdaptInfo& /*adaptInfo*/)
454 {
455  Dune::Timer t;
456 
457  bool adapted = grid_->preAdapt();
458  adapted |= grid_->adapt();
459  grid_->postAdapt();
460 
461  msg("adaptGrid needed {} seconds", t.elapsed());
462  return adapted ? MESH_ADAPTED : Flag(0);
463 }
464 
465 
466 template <class Traits>
468 buildAfterAdapt(AdaptInfo& /*adaptInfo*/, Flag /*flag*/, bool asmMatrix, bool asmVector)
469 {
470  Dune::Timer t;
471  Dune::Timer t2;
472 
473  auto localView = globalBasis_->localView();
474  for_each_node(localView.tree(), [&](auto&&, auto rowTp) -> void {
475  for_each_node(localView.tree(), [&](auto&&, auto colTp) -> void {
476  for (auto bc : boundaryConditions_[rowTp][colTp])
477  bc->init();
478  });
479  });
480 
481  t2.reset();
482 
483  // 1. init matrix and rhs vector and initialize dirichlet boundary conditions
484  systemMatrix_->init();
485  rhs_->init(sizeInfo(*globalBasis_), asmVector);
486 
487  // statistic about system size
488  if (Environment::mpiSize() > 1)
489  msg("{} local DOFs, {} global DOFs", rhs_->localSize(), rhs_->globalSize());
490  else
491  msg("{} local DOFs", rhs_->localSize());
492 
493  // 2. traverse grid and assemble operators on the elements
494  for (auto const& element : elements(gridView(), PartitionSet{})) {
495  localView.bind(element);
496 
497  if (asmMatrix)
498  systemMatrix_->assemble(localView, localView);
499  if (asmVector)
500  rhs_->assemble(localView);
501 
502  localView.unbind();
503  }
504 
505  // 3. finish matrix insertion and apply dirichlet boundary conditions
506  systemMatrix_->finish();
507  rhs_->finish();
508 
509  info(2," assemble operators needed {} seconds", t2.elapsed());
510  t2.reset();
511 
512  solution_->resize(sizeInfo(*globalBasis_));
513  for_each_node(localView.tree(), [&](auto&&, auto rowTp) -> void {
514  for_each_node(localView.tree(), [&](auto&&, auto colTp) -> void {
515  // finish boundary condition
516  for (auto bc : boundaryConditions_[rowTp][colTp])
517  bc->apply(*systemMatrix_, *solution_, *rhs_);
518  });
519  });
520 
521  info(2," assemble boundary conditions needed {} seconds", t2.elapsed());
522 
523  msg("fill-in of assembled matrix: {}", systemMatrix_->nnz());
524  msg("assemble needed {} seconds", t.elapsed());
525 }
526 
527 
528 template <class Traits>
530 writeFiles(AdaptInfo& adaptInfo, bool force)
531 {
532  Dune::Timer t;
533  for (auto writer : filewriter_)
534  writer->write(adaptInfo, force);
535  msg("writeFiles needed {} seconds", t.elapsed());
536 }
537 
538 } // end namespace AMDiS
constexpr bool isSet(Flag const &f) const
Checks whether all set bits of f.flags_ are set in flags_ too.
Definition: Flag.hpp:156
double relResidual() const
Returns relResidual_.
Definition: SolverInfo.hpp:50
void addPeriodicBC(BoundaryType id, WorldMatrix const &A, WorldVector const &b)
Definition: ProblemStat.inc.hpp:350
The Flag class encapsulates flags which represents simple information. Used e.g. while mesh traversal...
Definition: Flag.hpp:13
constexpr bool Functor
A Functor is a function F with signature Signature.
Definition: Concepts.hpp:134
void addDirichletBC(Predicate const &predicate, RowTreePath row, ColTreePath col, Values const &values)
Add boundary conditions to the system.
Definition: ProblemStat.inc.hpp:294
void error_exit(std::string const &str, Args &&... args)
print a message and exit
Definition: Output.hpp:142
constexpr bool Predicate
A predicate is a function that returns a boolean.
Definition: Concepts.hpp:142
decltype(auto) makeGridFunction(PreGridFct const &preGridFct, GridView const &gridView)
Generator for Gridfunctions from Expressions (PreGridfunctions)
Definition: GridFunction.hpp:154
A CreatorMap is used to construct objects, which types depends on key words determined at run time...
Definition: CreatorMap.hpp:29
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
Definition: BackupRestore.hpp:15
void setCreateMatrixData(bool b)
Sets createMatrixData_.
Definition: SolverInfo.hpp:104
A creator class for dune grids.
Definition: MeshCreator.hpp:52
Definition: ProblemStat.hpp:52
int info() const
Returns info.
Definition: SolverInfo.hpp:38
Creator class for filewriters depending on a given type name.
Definition: FileWriterCreator.hpp:25
void initialize(Flag initFlag, Self *adoptProblem=nullptr, Flag adoptFlag=INIT_NOTHING)
Initialisation of the problem.
Definition: ProblemStat.inc.hpp:23
int error() const
Returns error code in last run of an iterative solver.
Definition: SolverInfo.hpp:32
static int mpiSize()
Return the MPI_Size of the group created by Dune::MPIHelper.
Definition: Environment.hpp:74
void restore(Flag initFlag)
Read the grid and solution from backup files and initialize the problem.
Definition: ProblemStat.inc.hpp:108
void msg(std::string const &str, Args &&... args)
print a message
Definition: Output.hpp:98
static std::shared_ptr< Grid > create(std::string name)
Static create mthod. See create()
Definition: MeshCreator.hpp:71
static std::optional< T > get(std::string const &key)
Get parameter-values from parameter-tree.
Definition: Initfile.hpp:25
Flag markElements(AdaptInfo &adaptInfo) override
Implementation of ProblemStatBase::markElements.
Definition: ProblemStat.inc.hpp:391
CreatorInterfaceName< BaseClass > * named(CreatorInterface< BaseClass > *ptr)
cast a ptr of CreatorInterface to CreatorInterfaceName
Definition: CreatorInterface.hpp:64
Flag adaptGrid(AdaptInfo &adaptInfo) override
Implementation of ProblemStatBase::refineMesh.
Definition: ProblemStat.inc.hpp:453
void writeFiles(AdaptInfo &adaptInfo, bool force=false)
Writes output files. If force=true write even if timestep out of write rhythm.
Definition: ProblemStat.inc.hpp:530
void solve(AdaptInfo &adaptInfo, bool createMatrixData=true, bool storeMatrixData=false) override
Implementation of ProblemStatBase::solve.
Definition: ProblemStat.inc.hpp:362
Flag globalRefine(int n) override
Uniform global refinement by n level.
Definition: ProblemStat.inc.hpp:438
Holds adapt parameters and infos about the problem.
Definition: AdaptInfo.hpp:25
void info(int noInfoLevel, std::string const &str, Args &&... args)
prints a message, if Environment::infoLevel() >= noInfoLevel
Definition: Output.hpp:105
void for_each_node(Tree &&tree, PreFunc &&preFunc, LeafFunc &&leafFunc, PostFunc &&postFunc)
Traverse tree and visit each node.
Definition: Traversal.hpp:81
Definition: SolverInfo.hpp:11
void buildAfterAdapt(AdaptInfo &adaptInfo, Flag flag, bool asmMatrix=true, bool asmVector=true) override
Implementation of ProblemStatBase::buildAfterCoarse.
Definition: ProblemStat.inc.hpp:468
Flag globalCoarsen(int n) override
Uniform global grid coarsening by up to n level.
Definition: ProblemStat.inc.hpp:407
auto makeUniquePtr(Obj &&obj)
Create a unique_ptr by copy/move construction.
Definition: TypeTraits.hpp:107
void test_exit(bool condition, std::string const &str, Args &&... args)
test for condition and in case of failure print message and exit
Definition: Output.hpp:163
void setStoreMatrixData(bool b)
Sets storeMatrixData_.
Definition: SolverInfo.hpp:110
static std::unique_ptr< EstimatorMarker< Grid > > createMarker(std::string const &name, std::string const &component, Estimates const &est, std::shared_ptr< Grid > const &grid)
Creates a scalar marker depending on the strategy set in parameters.
Definition: Marker.inc.hpp:102
std::vector< int > const & boundaryIds() const
Return the filled vector of boundary ids. NOTE: not all creators support reading this.
Definition: MeshCreator.hpp:140
std::unique_ptr< FileWriterInterface > create(std::string type, std::string prefix, Indices... ii) const
Create a new FileWriter of type type
Definition: FileWriterCreator.hpp:47
double absResidual() const
Returns absResidual_.
Definition: SolverInfo.hpp:44