DUNE-ACFEM (2.5.1)

parabolicfemscheme.hh
1#ifndef __DUNE_ACFEM_PARABOLIC_FEMSCHEME_HH__
2#define __DUNE_ACFEM_PARABOLIC_FEMSCHEME_HH__
3
4// iostream includes
5#include <iostream>
6
7// include grid part
8#include <dune/fem/gridpart/adaptiveleafgridpart.hh>
9
10// include discrete function space
11#include <dune/fem/space/lagrange.hh>
12
13// adaptation ...
14#include <dune/fem/function/adaptivefunction.hh>
15#include <dune/fem/space/common/adaptmanager.hh>
16
17// include discrete function
18#include <dune/fem/function/blockvectorfunction.hh>
19
20// Newton's iteration
21#include <dune/fem/solver/newtoninverseoperator.hh>
22
23// natural interpolation
24#include <dune/fem/space/common/interpolate.hh>
25
26// include norms
27#include <dune/fem/misc/l2norm.hh>
28#include <dune/fem/misc/h1norm.hh>
29
30// include parameter handling
31#include <dune/fem/io/parameter.hh>
32
33// include vtk output and checkpointing
34#include <dune/fem/io/file/dataoutput.hh>
35#include <dune/fem/io/file/datawriter.hh>
36
37// blah blah blah if only Dune were as hip as it ought to be
38
39// local includes
40#include "../algorithms/femschemeinterface.hh"
41#include "../models/modelinterface.hh"
42#include "../functions/gridfunctionwrapper.hh"
43#include "../operators/ellipticoperator.hh"
44#include "../operators/functionals/l2innerproductfunctional.hh"
45#include "../operators/functionals/l2boundaryfunctional.hh"
46#include "../operators/functionals/functionalexpression.hh"
47#include "../operators/constraints/dirichletconstraints.hh"
48#include "../operators/constraints/emptyblockconstraints.hh"
49#include "../estimators/parabolicestimator.hh"
50#include "../estimators/trueerrorestimator.hh"
51#include "../algorithms/marking.hh"
52#include "../models/timeview.hh"
53#include "../common/dataoutput.hh"
54
55// Select an appropriate solver, depending on ModelType and solver-family (ISTL, PESC ...)
56#include "../common/solverselector.hh"
57
58namespace Dune {
59
60 namespace ACFem {
61
89 template<class DiscreteFunction,
90 class TimeProvider,
91 class ImplicitModel,
92 class ExplicitModel,
93 class InitialValue,
94 template<class> class QuadratureTraits = DefaultQuadratureTraits>
96 : public TransientAdaptiveFemScheme
97 {
98 public:
100 typedef DiscreteFunction DiscreteFunctionType;
101
103 typedef typename DiscreteFunctionType::GridType GridType;
104
106 typedef typename DiscreteFunctionType::GridPartType GridPartType;
107
109 typedef typename DiscreteFunctionType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
110
112 typedef TimeProvider TimeProviderType;
113
115 typedef ImplicitModel ImplicitModelImplementationType;
117 typedef typename ImplicitModelType::OperatorPartsType ImplicitOperatorPartsType;
119 typedef ExplicitModel ExplicitModelImplementationType;
121 typedef typename ExplicitModelType::OperatorPartsType ExplicitOperatorPartsType;
123
124 typedef
125 decltype(std::declval<ImplicitModelType>() + std::declval<ExplicitModelType>())
126 SumModelType;
128
130 typedef InitialValue InitialValueType;
131 typedef
132 typename GridFunctionConverter<InitialValueType, GridPartType>::WrappedGridFunctionType
133 InitialValueFunctionType;
134
136 typedef typename ImplicitModelType::FunctionSpaceType FunctionSpaceType;
137
138 // or: (must coincide)
139 //typedef typenema DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType;
140
141 // choose type of discrete function, Matrix implementation and solver implementation
143 typedef typename SolverSelectorType::LinearOperatorType LinearOperatorType;
144 typedef typename SolverSelectorType::LinearInverseOperatorType LinearInverseOperatorType;
145
147 typedef
148 Fem::NewtonInverseOperator<LinearOperatorType, LinearInverseOperatorType>
150
153 typedef Dune::Fem::RestrictProlongDefault<DiscreteFunctionType> RestrictionProlongationType;
154
156 typedef Dune::Fem::AdaptationManager<GridType, RestrictionProlongationType> AdaptationManagerType;
157
160
163
166
168 typedef
169 decltype(std::declval<DirichletBoundaryFunctionType>()
170 /
171 std::declval<typename DirichletWeightFunctionType::GridFunctionType>())
173
176
179
181 typedef DifferentiableEllipticOperator<LinearOperatorType, ImplicitOperatorPartsType, ConstraintsOperatorType, QuadratureTraits> ImplicitOperatorType;
182
185 typedef EllipticOperator<ExplicitOperatorPartsType,
188
200 // Isn't ther a more simpler way for doing this?
201 typedef typename SumModelType::BulkForcesFunctionType BulkForcesFunctionType;
202 typedef
203 L2InnerProductFunctional<DiscreteFunctionSpaceType, BulkForcesFunctionType>
204 BulkForcesFunctional;
205 typedef typename SumModelType::NeumannBoundaryFunctionType BoundaryFluxFunctionType;
206 typedef
207 L2BoundaryFunctional<DiscreteFunctionSpaceType,
208 BoundaryFluxFunctionType,
210 BoundaryFluxFunctional;
212
214 typedef
217 ExplicitModelImplementationType>
219
221 typedef MarkingStrategy<EstimatorType> MarkingStrategyType;
222
224 typedef
227
229 typedef MarkingStrategy<InitialEstimatorType> InitialMarkingStrategyType;
230
232 typedef InitialValueFunctionType ExactSolutionFunctionType; // slight abuse
233
235 typedef std::tuple<DiscreteFunctionType *, const ExactSolutionFunctionType *> IOTupleType;
236
238 typedef DataOutput<GridType, IOTupleType> DataOutputType;
239
241 typedef Dune::Fem::CheckPointer<GridType> CheckPointerType;
242
266 const TimeProviderType& timeProvider,
267 const ImplicitModelType& implicitModel,
268 const ExplicitModelType& explicitModel,
269 const InitialValueType& initialValue,
270 const std::string name = "heat")
271 : grid_(solution.space().gridPart().grid()),
272 timeProvider_(timeProvider),
273 implicitModel_(implicitModel),
274 explicitModel_(explicitModel),
275 name_(name),
276 gridPart_(solution.space().gridPart()),
277 discreteSpace_(solution.space()),
278 initialValue_(wrapToGridFunction("Initial Value", initialValue, gridPart_, discreteSpace_.order()+1)),
279 solution_(solution),
280 oldSolution_("old solution", discreteSpace_),
281 residual_("residual", discreteSpace_),
282 update_("update", discreteSpace_),
283 estimator_(timeProvider_, implicitModel_, explicitModel_, oldSolution_),
284 markingStrategy_(estimator_, name_ + ".adaptation.space"),
285 initialEstimator_(initialValue_()),
286 initialMarkingStrategy_(initialEstimator_, name_ + ".adaptation" + ".initial"),
287 // restriction/prolongation operator
288 restrictProlong_(0),
289 // adaptation manager
290 adaptationManager_(0),
291 // create Dirichlet contraints
292 boundaryValues_(implicitModel_.dirichletBoundaryFunction(gridPart_)),
293 boundaryWeight_(implicitModel_.dirichletWeightFunction(gridPart_)),
294 dirichletConstraints_(discreteSpace_,
295 boundaryValues_ / boundaryWeight_.function(),
296 implicitModel_.dirichletIndicator()),
297 // the elliptic operator (implicit)
298 implicitOperator_(implicitModel_.operatorParts(), dirichletConstraints_),
299 explicitOperator_(explicitModel_.operatorParts()),
300 // create linear operator (domainSpace,rangeSpace)
301 linearOperator_("assembled elliptic operator", discreteSpace_, discreteSpace_),
302 solverEps_(Dune::Fem::Parameter::getValue<double>(name_ + ".solvereps", 1e-8)),
303 // the right-hand-side
304 sumModel_(implicitModel_ + explicitModel_),
305 bulkForces_(sumModel_.bulkForcesFunction(gridPart_)),
306 boundaryFlux_(sumModel_.neumannBoundaryFunction(gridPart_)),
307 // other stuff
308 sequence_(-1),
309 ioTuple_(&solution_, &initialValue_()),
310 dataOutput_(0),
311 checkPointer_(grid_, timeProvider_)
312 {
313 // Make the solution of the current timestep persistent for check-pointing.
314 Dune::Fem::persistenceManager << solution_;
315
316 oldSolution_.clear();
317 solution_.clear();
318 }
319
320 virtual ~ParabolicFemScheme() {
321 if (dataOutput_) {
322 delete dataOutput_;
323 }
324 if (adaptationManager_) {
325 delete adaptationManager_;
326 }
327 if (restrictProlong_) {
328 delete restrictProlong_;
329 }
330 }
331
332 virtual void initialize()
333 {
334 // initialValue_ already has to be a grid-function
335
336 // apply natural interpolation
337 interpolate(initialValue_(), oldSolution_);
338
339 // copy it over
340 solution_.assign(oldSolution_);
341
342 // std::cout << solution_ << std::endl;
343 }
344
345 virtual void restart(const std::string& name = "") {
346 std::string name_ = name;
347 if (name_ == "") {
348 name_ = Dune::Fem::CheckPointerParameters().checkPointPrefix();
349 }
350 checkPointer_.restoreData(grid_, name_);
351 next();
352 }
353
358 virtual void next()
359 {
360 oldSolution_.assign(solution_);
361 }
362
364 virtual void solve(bool forceMatrixAssembling = true)
365 {
366 // Need to copy the data again in case the mesh was adapted in
367 // between because the AdaptationManager only works on a single function ATM.
368 solution_.assign(oldSolution_);
369
370 if (ImplicitModelConstituents::hasDirichletBoundary) {
371 // set Dirichlet constraints, new time
372 dirichletConstraints_.rebuildValues(); // force update
373 dirichletConstraints_.constrain(solution_);
374 }
375
376 // If we have a RHS _function_ then also take that into
377 // account. Note that it would be comparatively cheap to
378 // omit the enclosing if-claus, because the default
379 // configuration results into an efficient no-op. There is
380 // also no point to optimize further into the case
381 // Forces-but-no-Neumann or Neumann-but-no-forces for the
382 // same reason.
383 auto rhsFunctional =
384 BulkForcesFunctional(discreteSpace_, bulkForces_)
385 +
386 BoundaryFluxFunctional(discreteSpace_, boundaryFlux_, implicitModel_.neumannIndicator())
387 +
388 sumModel_.forcesFunctional(discreteSpace_);
389
390 if (!isZero(rhsFunctional)) { // avoid copying thousands of zero values
391 rhsFunctional.coefficients(residual_);
392 } else {
393 residual_.clear();
394 }
395
396 // Compute the contribution from the time derivative
397 explicitOperator_(oldSolution_, update_);
398 residual_ -= update_;
399
400 // do no spoil the new Dirichlet-values
401 dirichletConstraints_.zeroConstrain(residual_);
402
403 if (ImplicitOperatorPartsType::isLinear) {
404 linearSolve(residual_, forceMatrixAssembling); // sove Lu = res
405 } else {
406 nonLinearSolve(residual_); // solve Lu = res
407 }
408 }
409
410 protected:
411
414 {
415 assert(!ImplicitOperatorPartsType::isLinear);
416
417 NonLinearInverseOperatorType newton(implicitOperator_);
418
419 newton(rhs, solution_);
420
421 assert(newton.converged());
422 }
423
428 virtual void linearSolve(DiscreteFunctionType& rhs, bool forceMatrixAssembling)
429 {
430 assert(ImplicitOperatorPartsType::isLinear);
431
432 if (sequence_ != discreteSpace_.sequence() || forceMatrixAssembling) {
433 // assemble linear operator (i.e. setup matrix) The jacobian
434 // incorporates the constraints defined by the constraint class.
435 implicitOperator_.jacobian(solution_, linearOperator_);
436
437 // update sequence number
438 sequence_ = discreteSpace_.sequence();
439 }
440
441 // inverse operator using linear operator
442 //
443 // Note: ATM, we use the difference quotient for the time
444 // derivative, i.e. 1/tau as factor in from of the
445 // mass-matrix. Consequently, the solver-tolerance should be
446 // scaled by the time-step size, at least when we use any
447 // preconditioner.
448 //double scaledEps_ = solverEps_ * timeProvider_.inverseDeltaT();
449 LinearInverseOperatorType solver(linearOperator_, solverEps_, solverEps_);
450
451 // Apply the affine linear operator to the start value. This
452 // computes the initial residual. In the linear case, Lagrange space,
453 // Dirichlet data g, this computes
454 //
455 // update_ = A u - g
456 //
457 // where it is allowed that A is affine-linear, without having
458 // to apply a non-linear solver for this trivial non-linearity
459 implicitOperator_(solution_, update_);
460 rhs -= update_; // rhs = rhs - A u ...
461
462 // solve system. Zero initial values for the update_ vector are
463 // ok, because the information about the previous solution already
464 // is contained in residual_.
465 update_.clear();
466 solver(rhs, update_);
467
468 // subtract the update, difference should be the solution.
469 solution_ += update_; // ... so u = u + invA(rhs - Au)
470 }
471
472 public:
473
475 virtual bool mark(const double tolerance)
476 {
477 return markingStrategy_.mark(tolerance);
478 }
479
481 virtual double estimate()
482 {
483 return estimator_.estimate(solution_);
484 }
485
487 virtual double timeEstimate()
488 {
489 return estimator_.timeEstimate();
490 }
491
492 virtual double initialEstimate()
493 {
494 return initialEstimator_.estimate(solution_);
495 }
496
497 virtual bool initialMarking(const double tolerance)
498 {
499 return initialMarkingStrategy_.mark(tolerance);
500 }
501
503 virtual void adapt()
504 {
505 // there can only one adaptation manager per grid, and the
506 // RestrictionProlongationType which determines which
507 // functions are restricted/prolongated is built into the
508 // AdaptationManagerType. In order to allow for override by
509 // derived classes, allocation that adaptationManager_
510 // dynamically (and thus not at all, if ThisType::adapt() is
511 // never called.
512
513 if (!adaptationManager_) {
514 // restriction/prolongation operator
515 restrictProlong_ = new RestrictionProlongationType(oldSolution_);
516
517 // adaptation manager
518 adaptationManager_ = new AdaptationManagerType(grid_, *restrictProlong_);
519 }
520
521 // adaptation and load balancing
522 adaptationManager_->adapt();
523 if (adaptationManager_->loadBalance()) {
524 // TODO: re-initialize stuff as needed.
525 }
526 }
527
529 virtual int output()
530 {
531 // write checkpoint in order to be able to restart after a crash
532 checkPointer_.write(timeProvider_);
533
534 if (!dataOutput_) {
535 // NOTE: this should allocated dynamically, otherwise a
536 // derived class has problems to define completely different
537 // IO-schemes Also DataOutputType likes to already generate
538 // some files during construction, so make sure this never happens
539 dataOutput_ = new DataOutputType(grid_, ioTuple_, timeProvider_, DataOutputParameters());
540 }
541
542 if (!dataOutput_->willWrite(timeProvider_)) {
543 return -1;
544 }
545
546 dataOutput_->write(timeProvider_);
547
548 return dataOutput_->writeStep() - 1;
549 }
550
552 virtual double residual() const
553 {
554 assert(0);
555 return -1;
556 }
557
565 virtual double error() const
566 {
567 // can also use H1-norm
568 typedef Dune::Fem::L2Norm<GridPartType> NormType;
569
570 // the DiscreteFunctionSpaceAdapter sets the order per default
571 // to 111 == \infty. We set therefore the order just high
572 // enough that we see the asymptotic error. If the polynomial
573 // order is D, then the error should decay with (h^D). In
574 // principle it should thus suffice to use a quadrature of
575 // degree D.
576 NormType norm(gridPart_, 2*discreteSpace_.order()+2);
577 return norm.distance(initialValue_(), solution_);
578 }
579
580 virtual size_t size() const
581 {
582 // NOTE: slaveDofs.size() is always one TOO LARGE.
583 size_t numberOfDofs = discreteSpace_.blockMapper().size() - discreteSpace_.slaveDofs().size() + 1;
584
585 numberOfDofs = grid_.comm().sum(numberOfDofs);
586
587 return numberOfDofs;
588 }
589
590 protected:
591 GridType& grid_; // hierarchical grid
592 const TimeProviderType& timeProvider_; // maybe change to TimeView later ...
593 const ImplicitModelType& implicitModel_; // new time-step
594 const ExplicitModelType& explicitModel_; // old time-step
595 const std::string name_;
596
597 GridPartType& gridPart_; // grid part(view), here the leaf grid the discrete space is build with
598
599 const DiscreteFunctionSpaceType& discreteSpace_; // discrete function space
600 ExpressionStorage<InitialValueFunctionType> initialValue_; // this is already the GRID function
601 DiscreteFunctionType& solution_; // the unknown
602 DiscreteFunctionType oldSolution_; // the solution from the previous timestep
603 DiscreteFunctionType residual_; // initial residual
604 DiscreteFunctionType update_; // distance to solution
605
606 EstimatorType estimator_; // residual error estimator
607 MarkingStrategyType markingStrategy_;
608
609 InitialEstimatorType initialEstimator_;
610 InitialMarkingStrategyType initialMarkingStrategy_;
611
612 RestrictionProlongationType *restrictProlong_; // local restriction/prolongation object
613 AdaptationManagerType *adaptationManager_; // adaptation manager handling adaptation
614
615 DirichletBoundaryFunctionType boundaryValues_;
616 DirichletWeightFunctionType boundaryWeight_;
617 ConstraintsOperatorType dirichletConstraints_; // dirichlet boundary constraints
618
619 ImplicitOperatorType implicitOperator_; // the affine-linear operator.
620 ExplicitOperatorType explicitOperator_;
621
622 LinearOperatorType linearOperator_; // the linear operator (i.e. jacobian of the operator)
623
624 const double solverEps_ ; // eps for linear solver
625
627
628 SumModelType sumModel_;
629 BulkForcesFunctionType bulkForces_;
630 BoundaryFluxFunctionType boundaryFlux_;
631
632 mutable int sequence_; // sequence number
633
634 IOTupleType ioTuple_; // tuple with pointers
635 DataOutputType *dataOutput_; // data output class
636 CheckPointerType checkPointer_;
637
638 };
639
641
642 } // ACFem
643
644} // Dune
645
646#endif // __PARABOLIC_FEMSCHEME_HH__
Potentially overwrite some parameters.
Definition: dataoutput.hh:37
void zeroConstrain(DiscreteFunctionType &w) const
Unconditionally set the values of all masked DoFs to zero.
Definition: bulkblockconstraints.hh:250
void constrain(DiscreteFunctionType &w) const
The solution operator; unconditionally install the given constraints into the argument.
Definition: bulkblockconstraints.hh:224
void rebuildValues()
Interpolate the Dirichlet values in any case, but leave slave-DoFs as is if nothing changed.
Definition: dirichletconstraints.hh:190
A class defining an elliptic operator.
Definition: ellipticoperator.hh:86
Interface class for second order elliptic models.
Definition: modelinterface.hh:192
TraitsType::DirichletBoundaryFunctionType DirichletBoundaryFunctionType
A BoundarySupportedFunction which must be sub-ordinate to the DirichletIndicatorType.
Definition: modelinterface.hh:241
TraitsType::DirichletIndicatorType DirichletIndicatorType
A BoundarySupportedFunction which must be sub-ordinate to the DirichletIndicatorType.
Definition: modelinterface.hh:243
TraitsType::DirichletWeightFunctionType DirichletWeightFunctionType
A BoundarySupportedFunction which must be sub-ordinate to the DirichletIndicatorType.
Definition: modelinterface.hh:242
TraitsType::NeumannIndicatorType NeumannIndicatorType
A function modelling "force" terms in the bulk-phase.
Definition: modelinterface.hh:228
NeumannIndicatorType neumannIndicator() const
Generate an object to identify parts of the boundary subject to Neumann boundary conditions.
Definition: modelinterface.hh:318
Residual estimator for the heat equation.
Definition: parabolicestimator.hh:60
RangeFieldType estimate(const DiscreteFunctionType &uh)
calculate estimator
Definition: parabolicestimator.hh:250
Basic parabolic fem-scheme class.
Definition: parabolicfemscheme.hh:97
DiscreteFunction DiscreteFunctionType
Type of the discrete solution function.
Definition: parabolicfemscheme.hh:100
InitialValue InitialValueType
Initial value.
Definition: parabolicfemscheme.hh:130
ImplicitModelType::DirichletBoundaryFunctionType DirichletBoundaryFunctionType
type of Dirichlet boundary values
Definition: parabolicfemscheme.hh:162
virtual void nonLinearSolve(DiscreteFunctionType &rhs)
Run the full Newton-scheme ...
Definition: parabolicfemscheme.hh:413
virtual double estimate()
calculate error estimator
Definition: parabolicfemscheme.hh:481
Dune::Fem::AdaptationManager< GridType, RestrictionProlongationType > AdaptationManagerType
type of adaptation manager handling adapation and DoF compression
Definition: parabolicfemscheme.hh:156
virtual int output()
data I/O
Definition: parabolicfemscheme.hh:529
std::tuple< DiscreteFunctionType *, const ExactSolutionFunctionType * > IOTupleType
type of input/output tuple
Definition: parabolicfemscheme.hh:235
MarkingStrategy< EstimatorType > MarkingStrategyType
type of marking strategy for space adaptivity
Definition: parabolicfemscheme.hh:221
virtual double error() const
Calculate L2/H1 error.
Definition: parabolicfemscheme.hh:565
DiscreteFunctionType::GridPartType GridPartType
type of the grid view
Definition: parabolicfemscheme.hh:106
ImplicitModelType::FunctionSpaceType FunctionSpaceType
type of function space (scalar functions, )
Definition: parabolicfemscheme.hh:136
MarkingStrategy< InitialEstimatorType > InitialMarkingStrategyType
Initial marking.
Definition: parabolicfemscheme.hh:229
ImplicitModelType::DirichletWeightFunctionType DirichletWeightFunctionType
type of Dirichlet weight function
Definition: parabolicfemscheme.hh:165
DifferentiableEllipticOperator< LinearOperatorType, ImplicitOperatorPartsType, ConstraintsOperatorType, QuadratureTraits > ImplicitOperatorType
define differential operator, implicit part
Definition: parabolicfemscheme.hh:181
DiscreteFunctionType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
type of the discrete function space
Definition: parabolicfemscheme.hh:109
TrueErrorEstimator< InitialValueFunctionType, Dune::Fem::L2Norm< GridPartType > > InitialEstimatorType
intial estimator
Definition: parabolicfemscheme.hh:226
ParabolicFemScheme(DiscreteFunctionType &solution, const TimeProviderType &timeProvider, const ImplicitModelType &implicitModel, const ExplicitModelType &explicitModel, const InitialValueType &initialValue, const std::string name="heat")
Constructor.
Definition: parabolicfemscheme.hh:265
DirichletConstraints< LinearOperatorType, EffectiveDirichletFunctionType, DirichletIndicatorType > ConstraintsOperatorType
type of Dirichlet constraints
Definition: parabolicfemscheme.hh:175
ImplicitModel ImplicitModelImplementationType
type of the mathematical model
Definition: parabolicfemscheme.hh:115
TimeProvider TimeProviderType
access to the simulation time and time-step
Definition: parabolicfemscheme.hh:112
Dune::Fem::RestrictProlongDefault< DiscreteFunctionType > RestrictionProlongationType
type of restriction/prolongation projection for adaptive simulations (use default here,...
Definition: parabolicfemscheme.hh:153
Fem::NewtonInverseOperator< LinearOperatorType, LinearInverseOperatorType > NonLinearInverseOperatorType
Non-linear solver.
Definition: parabolicfemscheme.hh:149
virtual void linearSolve(DiscreteFunctionType &rhs, bool forceMatrixAssembling)
Perform only one step of the Newton scheme for the affine-linear case.
Definition: parabolicfemscheme.hh:428
DataOutput< GridType, IOTupleType > DataOutputType
type of data writer (produces VTK data)
Definition: parabolicfemscheme.hh:238
InitialValueFunctionType ExactSolutionFunctionType
adapter to turn exact solution into a grid function (for visualization)
Definition: parabolicfemscheme.hh:232
DifferentiableEmptyBlockConstraints< LinearOperatorType > EmptyConstraintsType
empty constraints for the explicit operator (old solution is already constrained)
Definition: parabolicfemscheme.hh:178
virtual double timeEstimate()
return the most recent time estimate
Definition: parabolicfemscheme.hh:487
virtual void initialize()
initialize the solution
Definition: parabolicfemscheme.hh:332
EllipticOperator< ExplicitOperatorPartsType, DiscreteFunctionType, DiscreteFunctionType, EmptyConstraintsType, QuadratureTraits > ExplicitOperatorType
explicit part of differential operator.
Definition: parabolicfemscheme.hh:187
virtual void next()
Close the current time-step.
Definition: parabolicfemscheme.hh:358
virtual bool mark(const double tolerance)
mark elements for adaptation
Definition: parabolicfemscheme.hh:475
DiscreteFunctionType::GridType GridType
type of hierarchic grid
Definition: parabolicfemscheme.hh:103
virtual void adapt()
do the adaptation for a given marking
Definition: parabolicfemscheme.hh:503
ParabolicEulerEstimator< DiscreteFunctionType, TimeProviderType, ImplicitModelImplementationType, ExplicitModelImplementationType > EstimatorType
type of error estimator
Definition: parabolicfemscheme.hh:218
virtual void solve(bool forceMatrixAssembling=true)
solve the system
Definition: parabolicfemscheme.hh:364
virtual size_t size() const
return some measure about the number of DOFs in use
Definition: parabolicfemscheme.hh:580
Dune::Fem::CheckPointer< GridType > CheckPointerType
type of check-pointer (dumps unaltered simulation data)
Definition: parabolicfemscheme.hh:241
decltype(std::declval< DirichletBoundaryFunctionType >()/std::declval< typename DirichletWeightFunctionType::GridFunctionType >()) EffectiveDirichletFunctionType
type of effective Dirichlet values
Definition: parabolicfemscheme.hh:172
virtual double residual() const
calculate residual (in small l^2)
Definition: parabolicfemscheme.hh:552
ImplicitModelType::DirichletIndicatorType DirichletIndicatorType
types for various data
Definition: parabolicfemscheme.hh:159
"estimator" which return the "true" error with respect to some given function in some standard norm.
Definition: trueerrorestimator.hh:26
RangeFieldType estimate(const DiscreteFunctionType &uh)
calculate estimator
Definition: trueerrorestimator.hh:64
constexpr bool isZero(Expression &&)
Specialize to evaluate to true for zero expressions.
Definition: expressionoperations.hh:469
static auto wrapToGridFunction(const std::string &name, const FunctionImp &f, const GridPart &g, unsigned order=111)
Possibly wrap a function into a GridFunctionWrapper.
Definition: gridfunctionwrapper.hh:342
A helper class which identifies which components are provided by the given model.
Definition: modelinterface.hh:447
Select one appropriate (linear) solver depending on whether the model is symmetric and/or semidefinit...
Definition: solverselector.hh:90
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)