ProblemStat.hpp 21.3 KB
Newer Older
1
2
3
4
#pragma once

#include <list>
#include <map>
5
6
7
8
9
#include <memory>
#include <string>
#include <tuple>
#include <utility>
#include <vector>
10

11
12
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
13
14
#include <dune/common/shared_ptr.hh>

15
#include <dune/grid/common/grid.hh>
16

17
#include <amdis/AdaptInfo.hpp>
18
#include <amdis/AdaptiveGrid.hpp>
19
#include <amdis/BiLinearForm.hpp>
20
#include <amdis/CreatorInterface.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
21
#include <amdis/CreatorMap.hpp>
22
#include <amdis/DirichletBC.hpp>
23
#include <amdis/DOFVector.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
24
//#include <amdis/Estimator.hpp>
25
26
27
#include <amdis/Flag.hpp>
#include <amdis/Initfile.hpp>
#include <amdis/LinearAlgebra.hpp>
28
#include <amdis/LinearForm.hpp>
29
#include <amdis/Marker.hpp>
30
#include <amdis/MeshCreator.hpp>
31
#include <amdis/OperatorList.hpp>
32
#include <amdis/PeriodicBC.hpp>
33
34
35
#include <amdis/ProblemStatBase.hpp>
#include <amdis/ProblemStatTraits.hpp>
#include <amdis/StandardProblemIteration.hpp>
36
#include <amdis/common/SharedPtr.hpp>
37
#include <amdis/common/TupleUtility.hpp>
38
#include <amdis/common/TypeTraits.hpp>
39
#include <amdis/GridFunctions.hpp>
40
#include <amdis/gridfunctions/DiscreteFunction.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
41
#include <amdis/io/FileWriterBase.hpp>
42
#include <amdis/typetree/Traits.hpp>
43
#include <amdis/typetree/TreePath.hpp>
44

45
namespace AMDiS
46
{
47
48
49
50
  // forward declaration
  template <class Traits>
  class ProblemInstat;

51
  template <class Traits>
52
  class ProblemStat
Praetorius, Simon's avatar
Praetorius, Simon committed
53
      : public ProblemStatBase
54
      , public StandardProblemIterationAdaptor<ProblemStat<Traits>>
55
  {
56
    using Self = ProblemStat;
57

58
59
    friend class ProblemInstat<Traits>;

60
  public: // typedefs and static constants
61

62
    using GlobalBasis = typename Traits::GlobalBasis;
Praetorius, Simon's avatar
Praetorius, Simon committed
63
    using GridView    = typename GlobalBasis::GridView;
Praetorius, Simon's avatar
Praetorius, Simon committed
64
    using Grid        = AdaptiveGrid_t<typename GridView::Grid>;
Praetorius, Simon's avatar
Praetorius, Simon committed
65
    using Element     = typename GridView::template Codim<0>::Entity;
66
    using WorldVector = typename Element::Geometry::GlobalCoordinate;
67
    using WorldMatrix = FieldMatrix<typename WorldVector::field_type, WorldVector::dimension, WorldVector::dimension>;
68

69
    /// Dimension of the grid
70
    static constexpr int dim = Grid::dimension;
71

72
    /// Dimension of the world
73
    static constexpr int dow = Grid::dimensionworld;
74

Praetorius, Simon's avatar
Praetorius, Simon committed
75
76
77
78
79
  private:
    using LinAlgTraits = typename Traits::LinAlgTraits;
    using Mat = typename LinAlgTraits::template MatrixImpl<typename Traits::CoefficientType>;
    using Vec = typename LinAlgTraits::template VectorImpl<typename Traits::CoefficientType>;
    using PartitionSet = typename LinAlgTraits::PartitionSet;
80

Praetorius, Simon's avatar
Praetorius, Simon committed
81
82
83
84
85
  public:
    using LinearSolver = LinearSolverInterface<Mat,Vec>;
    using SystemMatrix = BiLinearForm<GlobalBasis, GlobalBasis, typename Traits::CoefficientType, LinAlgTraits>;
    using SystemVector = LinearForm<GlobalBasis, typename Traits::CoefficientType, LinAlgTraits>;
    using SolutionVector = DOFVector<GlobalBasis, typename Traits::CoefficientType, LinAlgTraits>;
86

87
  public:
88
89
    /**
     * \brief Constructor. Takes the name of the problem that is used to
90
     * access values corresponding to this problem in the parameter file.
91
     **/
92
    explicit ProblemStat(std::string const& name)
93
      : name_(name)
94
    {}
95

96
    /// Constructor taking additionally a grid that is used
97
    /// instead of the default created grid, \ref ProblemStat
98
99
    template <class Grid_>
    ProblemStat(std::string const& name, Grid_&& grid)
100
      : ProblemStat(name)
101
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
102
      adoptGrid(wrap_or_share(FWD(grid)));
103
104
    }

105
    /// \brief Constructor taking a grid and basis.
106
    /// Wraps both in shared pointers.
107
    template <class Grid_, class Basis_, class B_ = Underlying_t<Basis_>,
108
      REQUIRES(Concepts::GlobalBasis<B_>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
109
    ProblemStat(std::string const& name, Grid_&& grid, Basis_&& globalBasis)
110
      : ProblemStat(name, FWD(grid))
111
    {
112
      adoptGlobalBasis(wrap_or_share(FWD(globalBasis)));
113
    }
Praetorius, Simon's avatar
Praetorius, Simon committed
114

115
116
    /// \brief Constructor taking a grid and pre-basis factory to create a global basis
    /// on the fly.
117
    template <class Grid_, class PBF_, class GV_ = typename Underlying_t<Grid_>::LeafGridView,
118
      REQUIRES(Concepts::PreBasisFactory<PBF_, GV_, MultiIndex_t<PBF_>>)>
119
120
121
    ProblemStat(std::string const& name, Grid_&& grid, PBF_ const& preBasisFactory)
      : ProblemStat(name, FWD(grid))
    {
122
      adoptGlobalBasis(makeSharedPtr(GlobalBasis{grid_->leafGridView(), preBasisFactory}));
123
    }
124

Praetorius, Simon's avatar
Praetorius, Simon committed
125
    /// \brief Initialisation of the problem.
126
    /**
Praetorius, Simon's avatar
Praetorius, Simon committed
127
128
     * Parameters read in initialize()
     *   - `[GRID_NAME]->global refinements`:  nr of initial global refinements
129
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
130
    void initialize(Flag initFlag, Self* adoptProblem = nullptr, Flag adoptFlag = INIT_NOTHING);
131

Praetorius, Simon's avatar
Praetorius, Simon committed
132
133
134
135
136
137
138
139
    /// \brief Read the grid and solution from backup files and initialize the problem
    /**
     * Parameters read in restore() for problem with name 'PROB'
     *   - `[PROB]->restore->grid`:      name of the grid backup file
     *   - `[PROB]->restore->solution`:  name of the solution backup file
     **/
    void restore(Flag initFlag);

140

Praetorius, Simon's avatar
Praetorius, Simon committed
141
    /// Add an operator to \ref A.
142
    /** @{ */
Praetorius, Simon's avatar
Praetorius, Simon committed
143
144
145
146
147
148
149
    /// Operator evaluated on the whole element
    /**
     * Adds an operator to the list of element operators to be assembled in
     * quadrature points inside the element.
     *
     * \param op   A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator
     * \param row  TreePath identifying the sub-basis in the global basis tree
150
     *             corresponding to the row basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
151
     * \param col  TreePath identifying the sub-basis in the global basis tree
152
     *             corresponding to the column basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
153
154
155
156
157
158
159
     *
     * Example:
     * ```
     * auto op = makeOperator(tag::test_trial{}, 1.0/tau);
     * prob.addMatrixOperator(op, _0, _0);
     * ```
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
160
    template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
161
162
    void addMatrixOperator(Operator const& op, RowTreePath row = {}, ColTreePath col = {})
    {
163
164
165
166
167
168
169
      static constexpr bool isValidTreePath =
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, RowTreePath> &&
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, ColTreePath>;
      static_assert(isValidTreePath, "Invalid row and/or col treepath passed to addMatrixOperator!");

      if constexpr (isValidTreePath)
        systemMatrix_->addOperator(tag::element_operator<Element>{}, op, row, col);
170
    }
171

Praetorius, Simon's avatar
Praetorius, Simon committed
172
173
174
175
176
    /// Operator evaluated on the boundary of the domain with boundary index `b`
    /**
     * Adds an operator to the list of boundary operators to be assembled in
     * quadrature points on the boundary intersections.
     *
Praetorius, Simon's avatar
Praetorius, Simon committed
177
     * \param b    Boundary identifier where to assemble this operator. Can be
Praetorius, Simon's avatar
Praetorius, Simon committed
178
179
180
     *             constructed from an integer. \see BoundaryType
     * \param op   A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator
     * \param row  TreePath identifying the sub-basis in the global basis tree
181
     *             corresponding to the row basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
182
     * \param col  TreePath identifying the sub-basis in the global basis tree
183
     *             corresponding to the column basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
184
185
186
187
188
189
190
     *
     * Example:
     * ```
     * auto op = makeOperator(tag::test_trial{}, alpha);
     * prob.addMatrixOperator(BoundaryType{1}, op, _0, _0);
     * ```
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
191
    template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
192
193
    void addMatrixOperator(BoundaryType b, Operator const& op, RowTreePath row = {}, ColTreePath col = {})
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
194
      using I = typename GridView::Intersection;
195
196
197
198
199
200
201
      static constexpr bool isValidTreePath =
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, RowTreePath> &&
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, ColTreePath>;
      static_assert(isValidTreePath, "Invalid row and/or col treepath passed to addMatrixOperator!");

      if constexpr (isValidTreePath)
        systemMatrix_->addOperator(tag::boundary_operator<I>{*boundaryManager_,b}, op, row, col);
202
    }
203
    /** @} */
204
205


Praetorius, Simon's avatar
Praetorius, Simon committed
206
    /// Add an operator to \ref rhs.
207
    /** @{ */
Praetorius, Simon's avatar
Praetorius, Simon committed
208
209
210
211
212
213
214
    /// Operator evaluated on the whole element
    /**
     * Adds an operator to the list of element operators to be assembled in
     * quadrature points inside the element.
     *
     * \param op    A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator
     * \param path  TreePath identifying the sub-basis in the global basis tree
215
     *              corresponding to the row basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
216
217
218
219
220
221
222
     *
     * Example:
     * ```
     * auto op = makeOperator(tag::test{}, probInstat.getOldSolution(0) / tau);
     * prob.addVectorOperator(op, _0);
     * ```
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
223
    template <class Operator, class TreePath = RootTreePath>
224
225
    void addVectorOperator(Operator const& op, TreePath path = {})
    {
226
227
228
229
230
231
      static constexpr bool isValidTreePath =
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, TreePath>;
      static_assert(isValidTreePath, "Invalid treepath passed to addVectorOperator!");

      if constexpr (isValidTreePath)
        rhs_->addOperator(tag::element_operator<Element>{}, op, path);
232
    }
233

Praetorius, Simon's avatar
Praetorius, Simon committed
234
235
236
237
238
    /// Operator evaluated on the boundary of the domain with boundary index `b`
    /**
     * Adds an operator to the list of boundary operators to be assembled in
     * quadrature points on the boundary intersections.
     *
Praetorius, Simon's avatar
Praetorius, Simon committed
239
     * \param b     Boundary identifier where to assemble this operator. Can be
Praetorius, Simon's avatar
Praetorius, Simon committed
240
241
242
     *              constructed from an integer. \see BoundaryType
     * \param op    A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator
     * \param path  TreePath identifying the sub-basis in the global basis tree
243
     *              corresponding to the row basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
244
245
246
247
248
249
250
     *
     * Example:
     * ```
     * auto op = makeOperator(tag::test{}, [g](auto const& x) { return g(x); });
     * prob.addVectorOperator(BoundaryType{1}, op, _0);
     * ```
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
251
    template <class Operator, class TreePath = RootTreePath>
252
253
    void addVectorOperator(BoundaryType b, Operator const& op, TreePath path = {})
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
254
      using I = typename GridView::Intersection;
255
256
257
258
259
260
      static constexpr bool isValidTreePath =
        Concepts::ValidTreePath<typename GlobalBasis::LocalView::Tree, TreePath>;
      static_assert(isValidTreePath, "Invalid treepath passed to addVectorOperator!");

      if constexpr (isValidTreePath)
        rhs_->addOperator(tag::boundary_operator<I>{*boundaryManager_,b}, op, path);
261
    }
262
    /** @} */
263

264

Praetorius, Simon's avatar
Praetorius, Simon committed
265
266
267
268
269
270
271
272
273
274
    /// Add boundary conditions to the system
    /** @{ */
    /// Dirichlet boundary condition
    /**
     * Enforce Dirichlet boundary values for the solution vector on boundary
     * regions identified by the predicate.
     *
     * \param predicate  Functor `bool(WorldVector)` returning true for all
     *                   DOFs on the boundary that should be assigned a value.
     * \param row        TreePath identifying the sub-basis in the global basis tree
275
     *                   corresponding to the row basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
276
     * \param col        TreePath identifying the sub-basis in the global basis tree
277
     *                   corresponding to the column basis. \see makeTreePath()
Praetorius, Simon's avatar
Praetorius, Simon committed
278
279
280
281
282
283
284
285
286
     * \param values     Functor `Range(WorldVector)` or any \ref GridFunction
     *                   that is evaluated in the DOFs identified by the predicate.
     *
     * Example:
     * ```
     * prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8; }, 0, 0,
     *                     [](auto const& x) { return 0.0; });
     * ```
     **/
287
    template <class Predicate, class RowTreePath, class ColTreePath, class Values>
288
    void addDirichletBC(Predicate const& predicate,
289
                        RowTreePath row, ColTreePath col,
Praetorius, Simon's avatar
Praetorius, Simon committed
290
                        Values const& values);
291
292
293
294
295
296

    template <class RowTreePath, class ColTreePath, class Values>
    void addDirichletBC(BoundaryType id,
                        RowTreePath row, ColTreePath col,
                        Values const& values);

297
298
299
300
301
302
    template <class Identifier, class Values>
    void addDirichletBC(Identifier&& id, Values&& values)
    {
      addDirichletBC(FWD(id), RootTreePath{}, RootTreePath{}, FWD(values));
    }

303
304
305
    /// Add a periodic boundary conditions to the system, by specifying a face transformation
    /// y = A*x + b of coordinates. We assume, that A is orthonormal.
    void addPeriodicBC(BoundaryType id, WorldMatrix const& A, WorldVector const& b);
Praetorius, Simon's avatar
Praetorius, Simon committed
306
307
    /** @} */

308

309
  public:
310

311
312
313
314
315
    /// Implementation of \ref StandardProblemIteration::oneIteration.
    Flag oneIteration(AdaptInfo& adaptInfo, Flag toDo = FULL_ITERATION) override
    {
      return StandardProblemIteration::oneIteration(adaptInfo, toDo);
    }
316

317
    /// Implementation of \ref ProblemStatBase::buildAfterCoarse
318
319
320
321
    void buildAfterAdapt(AdaptInfo& adaptInfo,
                         Flag flag,
                         bool asmMatrix = true,
                         bool asmVector = true) override;
322

Praetorius, Simon's avatar
Praetorius, Simon committed
323
324
325
326
327
328
329
    /// \brief Assemble the linear system by calling \ref buildAfterAdapt with
    /// `asmMatrix` and `asmVector` set to true.
    void assemble(AdaptInfo& adaptInfo)
    {
      buildAfterAdapt(adaptInfo, Flag{0}, true, true);
    }

330
331
332
333
334
335
    /// Implementation of \ref ProblemStatBase::solve
    void solve(AdaptInfo& adaptInfo,
               bool createMatrixData = true,
               bool storeMatrixData = false) override;

    /// Implementation of \ref ProblemStatBase::estimate.
336
    void estimate(AdaptInfo& /*adaptInfo*/) override { /* do nothing. */ }
337
338
339
340
341
342
343
344
345
346
347
348
349

    /// Implementation of \ref ProblemStatBase::refineMesh.
    Flag adaptGrid(AdaptInfo& adaptInfo) override;

    /// Implementation of \ref ProblemStatBase::markElements.
    Flag markElements(AdaptInfo& adaptInfo) override;

    /// Uniform global grid coarsening by up to n level
    Flag globalCoarsen(int n) override;

    /// Uniform global refinement by n level
    Flag globalRefine(int n) override;

Praetorius, Simon's avatar
Praetorius, Simon committed
350
    /// Writes output files. If force=true write even if timestep out of write rhythm.
351
    void writeFiles(AdaptInfo& adaptInfo, bool force = false);
352

Praetorius, Simon's avatar
Praetorius, Simon committed
353

354
  public: // get-methods
355

Praetorius, Simon's avatar
Praetorius, Simon committed
356
    /// Implementation of \ref ProblemStatBase::name
357
    std::string const& name() const override { return name_; }
Praetorius, Simon's avatar
Praetorius, Simon committed
358

359
    /// Return the \ref grid_
360
361
    std::shared_ptr<Grid>       grid()       { return grid_; }
    std::shared_ptr<Grid const> grid() const { return grid_; }
Praetorius, Simon's avatar
Praetorius, Simon committed
362

363
    /// Return the gridView of the basis
364
    GridView gridView() const { return globalBasis_->gridView(); }
Praetorius, Simon's avatar
Praetorius, Simon committed
365

366
    /// Return the boundary manager to identify boundary segments
367
368
    std::shared_ptr<BoundaryManager<Grid>>       boundaryManager()       { return boundaryManager_; }
    std::shared_ptr<BoundaryManager<Grid> const> boundaryManager() const { return boundaryManager_; }
369

Praetorius, Simon's avatar
Praetorius, Simon committed
370
    /// Return the \ref globalBasis_
371
372
    std::shared_ptr<GlobalBasis>       globalBasis()       { return globalBasis_; }
    std::shared_ptr<GlobalBasis const> globalBasis() const { return globalBasis_; }
Praetorius, Simon's avatar
Praetorius, Simon committed
373
374

    /// Return a reference to the linear solver, \ref linearSolver
Praetorius, Simon's avatar
Praetorius, Simon committed
375
376
    std::shared_ptr<LinearSolver>       solver()       { return linearSolver_; }
    std::shared_ptr<LinearSolver const> solver() const { return linearSolver_; }
Praetorius, Simon's avatar
Praetorius, Simon committed
377

378
    /// Returns a reference to system-matrix, \ref systemMatrix_
379
380
    std::shared_ptr<SystemMatrix>       systemMatrix()       { return systemMatrix_; }
    std::shared_ptr<SystemMatrix const> systemMatrix() const { return systemMatrix_; }
381
382

    /// Returns a reference to the solution vector, \ref solution_
Praetorius, Simon's avatar
Praetorius, Simon committed
383
384
    std::shared_ptr<SolutionVector>       solutionVector()       { return solution_; }
    std::shared_ptr<SolutionVector const> solutionVector() const { return solution_; }
385
386

    /// Return a reference to the rhs system-vector, \ref rhs
387
388
    std::shared_ptr<SystemVector>       rhsVector()       { return rhs_; }
    std::shared_ptr<SystemVector const> rhsVector() const { return rhs_; }
389

390
391

    /// Return a mutable view to a solution component
392
393
    template <class... Indices>
    auto solution(Indices... ii)
394
    {
395
      assert(bool(solution_) && "You have to call initialize() before.");
396
      return valueOf(*solution_, ii...);
397
    }
398

399
    /// Return a const view to a solution component
400
401
    template <class... Indices>
    auto solution(Indices... ii) const
402
    {
403
      assert(bool(solution_) && "You have to call initialize() before.");
404
      return valueOf(*solution_, ii...);
405
    }
406
407


Praetorius, Simon's avatar
Praetorius, Simon committed
408
  public: // set-methods
409

410
    /// Set a new linear solver for the problem
Praetorius, Simon's avatar
Praetorius, Simon committed
411
412
    template <class Solver_>
    void setSolver(Solver_&& solver)
413
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
414
      linearSolver_ = wrap_or_share(FWD(solver));
Praetorius, Simon's avatar
Praetorius, Simon committed
415
416
    }

417

418
    /// Set the grid. Stores pointer and initializes feSpaces
Praetorius, Simon's avatar
Praetorius, Simon committed
419
    /// matrices and vectors, as well as markers and file-writers.
Praetorius, Simon's avatar
Praetorius, Simon committed
420
421
422
    /// If grid is given as reference, wrap it into a non-destroying shared_ptr
    template <class Grid_>
    void setGrid(Grid_&& grid)
423
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
424
      adoptGrid(wrap_or_share(FWD(grid)));
425
      createGlobalBasis();
426
      createMatricesAndVectors();
Praetorius, Simon's avatar
Praetorius, Simon committed
427
      createMarker();
428
429
430
      createFileWriter();
    }

431
432
433
434
435

    /// Store the shared_ptr and the name of the marker in the problem
    /**
     * Note: multiple markers can be added but must have different names
     **/
Praetorius, Simon's avatar
Praetorius, Simon committed
436
437
    template <class Marker_>
    void addMarker(Marker_&& m)
Praetorius, Simon's avatar
Praetorius, Simon committed
438
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
439
440
      auto marker = wrap_or_share(FWD(m));
      auto it = marker_.emplace(marker->name(), marker);
Praetorius, Simon's avatar
Praetorius, Simon committed
441
      if (marker_.size() > 1)
442
        it.first->second->setMaximumMarking(true);
Praetorius, Simon's avatar
Praetorius, Simon committed
443
444
    }

445
446
447
448
449
450
451
452
453
454
455
    /// Remove a marker with the given name from the problem
    void removeMarker(std::string name)
    {
      std::size_t num = marker_.erase(name);
      test_warning(num == 1, "A marker with the given name '{}' does not exist.", name);
    }

    /// Remove a marker from the problem
    void removeMarker(Marker<Grid> const& marker)
    {
      removeMarker(marker.name());
Praetorius, Simon's avatar
Praetorius, Simon committed
456
457
    }

458

459
  protected: // initialization methods
460

Praetorius, Simon's avatar
Praetorius, Simon committed
461
462
463
464
465
466
    void createGlobalBasis();
    void createGrid();
    void createMatricesAndVectors();
    void createSolver();
    void createMarker();
    void createFileWriter();
467

468
    void adoptGlobalBasis(std::shared_ptr<GlobalBasis> globalBasis)
469
    {
470
      globalBasis_ = std::move(globalBasis);
471
      initGlobalBasis();
472
473
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
474
    void adoptGrid(std::shared_ptr<Grid> const& grid,
475
                   std::shared_ptr<BoundaryManager<Grid>> const& boundaryManager)
476
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
477
      grid_ = grid;
478
      boundaryManager_ = boundaryManager;
479
      Parameters::get(name_ + "->mesh", gridName_);
480
    }
481

Praetorius, Simon's avatar
Praetorius, Simon committed
482
483
484
485
486
487
    void adoptGrid(std::shared_ptr<Grid> const& grid)
    {
      adoptGrid(grid, std::make_shared<BoundaryManager<Grid>>(grid));
    }

    void adoptGrid(std::shared_ptr<typename Grid::HostGrid> const& hostGrid)
488
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
489
      auto grid = std::make_shared<Grid>(hostGrid);
490
491
492
      adoptGrid(grid, std::make_shared<BoundaryManager<Grid>>(grid));
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
493
  private:
494

Praetorius, Simon's avatar
Praetorius, Simon committed
495
496
    void createGlobalBasisImpl(std::true_type);
    void createGlobalBasisImpl(std::false_type);
497

498
    void initGlobalBasis();
499

500
  private:
501
    /// Name of this problem.
502
    std::string name_;
503

504
    /// Grid of this problem.
Praetorius, Simon's avatar
Praetorius, Simon committed
505
    std::shared_ptr<Grid> grid_;
506

507
    /// Name of the grid
508
    std::string gridName_ = "mesh";
509

510
511
512
    /// Management of boundary conditions
    std::shared_ptr<BoundaryManager<Grid>> boundaryManager_;

513
    /// FE space of this problem.
514
    std::shared_ptr<GlobalBasis> globalBasis_;
515

516
    /// A FileWriter object
517
    std::list<std::shared_ptr<FileWriterInterface>> filewriter_;
518

519
    /// Pointer to the adaptation markers
520
    std::map<std::string, std::shared_ptr<Marker<Grid>>> marker_;
Praetorius, Simon's avatar
Praetorius, Simon committed
521
522
523

    /// Pointer to the estimators for this problem
//    std::vector<Estimator*> estimator;
524

525
    /// An object of the linearSolver Interface
Praetorius, Simon's avatar
Praetorius, Simon committed
526
    std::shared_ptr<LinearSolver> linearSolver_;
527

Praetorius, Simon's avatar
Praetorius, Simon committed
528
    /// Matrix that is filled during assembling
529
    std::shared_ptr<SystemMatrix> systemMatrix_;
530

Praetorius, Simon's avatar
Praetorius, Simon committed
531
    /// Vector with the solution components
Praetorius, Simon's avatar
Praetorius, Simon committed
532
    std::shared_ptr<SolutionVector> solution_;
533

Praetorius, Simon's avatar
Praetorius, Simon committed
534
    /// Vector (load-vector) corresponding to the right-hand side
535
    /// of the equation, filled during assembling
536
    std::shared_ptr<SystemVector> rhs_;
537

Praetorius, Simon's avatar
Praetorius, Simon committed
538
539
540
541
    /// A vector with the local element error estimates
    /// for each node in the basis tree, indexed by [to_string(treePath)][element index]
    std::map<std::string, std::vector<double>> estimates_;

542
  private: // some internal data-structures
543
    BoundaryConditions<SystemMatrix, SolutionVector, SystemVector, GlobalBasis, GlobalBasis> boundaryConditions_;
544
  };
545

Praetorius, Simon's avatar
Praetorius, Simon committed
546

547
548
549
550
551
552
  namespace Impl
  {
    template <class Grid, class B, class = void>
    struct DeducedProblemTraits;

    template <class Grid, class PB>
553
    struct DeducedProblemTraits<Grid,GlobalBasis<PB>,void>
554
    {
555
      using type = DefaultProblemTraits<GlobalBasis<PB>>;
556
557
558
559
    };

    template <class G, class PBF>
    struct DeducedProblemTraits<G,PBF,
560
      std::enable_if_t<Concepts::PreBasisFactory<PBF, typename G::LeafGridView, MultiIndex_t<PBF>>>>
561
562
563
    {
      using Grid = AdaptiveGrid_t<G>;
      using GridView = typename Grid::LeafGridView;
564
      using Basis = decltype(GlobalBasis{std::declval<GridView>(),std::declval<PBF>()});
565
566
567
568
569
570
571
572
573

      using type = DefaultProblemTraits<Basis>;
    };

    template <class Grid, class Basis>
    using DeducedProblemTraits_t = typename DeducedProblemTraits<Grid,Basis>::type;
  }


574
  // Deduction guide
575
576
577
  template <class Grid, class Basis>
  ProblemStat(std::string name, Grid&& grid, Basis&& globalBasis)
    -> ProblemStat<Impl::DeducedProblemTraits_t<Underlying_t<Grid>,Underlying_t<Basis>>>;
578

579

580
581
582
583
  // mark templates as explicitly instantiated in cpp file
  extern template class ProblemStat<YaspGridBasis<2,1>>;
  extern template class ProblemStat<YaspGridBasis<2,2>>;

584
} // end namespace AMDiS
585
586

#include "ProblemStat.inc.hpp"