4 #ifndef DUNE_PDELAB_MULTISTEP_CACHE_HH
5 #define DUNE_PDELAB_MULTISTEP_CACHE_HH
12 #include <dune/common/exceptions.hh>
13 #include <dune/common/shared_ptr.hh>
45 template<
class Step =
int,
class Time =
double>
108 virtual bool isAffine(std::size_t order, Step step)
const
121 Step requested, Step available)
const
126 Step requested, Step available)
const
131 Step available)
const
155 virtual void preStep(Step step, Step stepsOfScheme_,
156 Time endTime_, Time dt_)
232 template<
class VectorU,
class VectorV,
class Matrix,
233 class Step = int,
class Time =
double>
238 typedef std::map<Step, shared_ptr<const Matrix> > MatrixMap;
239 typedef typename MatrixMap::const_iterator MatrixIterator;
240 typedef std::map<Step, shared_ptr<const VectorV> > ResidualMap;
241 typedef typename ResidualMap::const_iterator ResidualIterator;
242 typedef std::map<Step, shared_ptr<const VectorU> > UnknownMap;
243 typedef typename UnknownMap::const_iterator UnknownIterator;
246 std::vector<ResidualMap> residualValues;
249 std::vector<MatrixMap> jacobians;
250 std::vector<ResidualMap> zeroResiduals;
253 MatrixMap composedJacobians;
259 shared_ptr<Policy> policy;
267 shared_ptr<Policy>(
new Policy)) :
272 "MultiStepCache constructed with policy == NULL");
282 DUNE_THROW(
CacheError,
"MultiStepCache::setPolicy(): attempt to set "
305 shared_ptr<const VectorV>
307 if(order < residualValues.size()) {
308 ResidualIterator it = residualValues[order].find(step);
309 if(it != residualValues[order].end())
312 DUNE_THROW(
NotInCache,
"MultiStepCache::getResidualValue(): The "
313 "requested residual value "
314 "r_" << order <<
"(t_" << step <<
",u_" << step <<
") is "
328 const shared_ptr<const VectorV> &residualValue)
330 if(!policy->cacheResidualValue(order, step))
332 if(order >= residualValues.size())
333 residualValues.resize(order+1);
334 if(!residualValues[order].insert(std::make_pair(step, residualValue))
337 "r_" << order <<
"(t_" << step <<
", u_" << step <<
") "
338 "is already in the cache!");
360 shared_ptr<const Matrix>
362 if(order < jacobians.size()) {
363 MatrixIterator it = jacobians[order].find(step);
364 const MatrixIterator &end = jacobians[order].end();
369 for(it = jacobians[order].begin(); it != end; ++it)
370 if(policy->canReuseJacobian(order, step, it->first))
372 return jacobians[order][step] = it->second;
374 DUNE_THROW(
NotInCache,
"MultiStepCache::getJacobian(): The requested "
375 "Jacobian J(r_" << order <<
"|_t_" << step <<
") is not in "
389 const shared_ptr<const Matrix> &jacobian)
391 if(!policy->cacheJacobian(order, step))
393 if(order >= jacobians.size())
394 jacobians.resize(order+1);
395 if(!jacobians[order].insert(std::make_pair(step, jacobian)).second)
396 DUNE_THROW(
AlreadyInCache,
"MultiStepCache::setJacobian(): Jacobian "
397 "J(r_" << order <<
"|_t_" << step <<
") is already in "
420 shared_ptr<const VectorV>
422 if(order < zeroResiduals.size()) {
423 ResidualIterator it = zeroResiduals[order].find(step);
424 const ResidualIterator &end = zeroResiduals[order].end();
429 for(it = zeroResiduals[order].begin(); it != end; ++it)
430 if(policy->canReuseZeroResidual(order, step, it->first))
432 return zeroResiduals[order][step] = it->second;
434 DUNE_THROW(
NotInCache,
"MultiStepCache::getZeroResidual(): The "
435 "requested zero-residual "
436 "r_" << order <<
"(t_" << step <<
",0) is not in the "
450 const shared_ptr<const VectorV> &zeroResidual)
452 if(!policy->cacheZeroResidual(order, step))
454 if(order >= zeroResiduals.size())
455 zeroResiduals.resize(order+1);
456 if(!zeroResiduals[order].insert(std::make_pair(step, zeroResidual))
459 "r_" << order <<
"(t_" << step <<
", 0) is already in "
481 shared_ptr<const Matrix>
483 MatrixIterator it = composedJacobians.find(step);
484 const MatrixIterator &end = composedJacobians.end();
489 for(it = composedJacobians.begin(); it != end; ++it)
490 if(policy->canReuseComposedJacobian(step, it->first))
492 return composedJacobians[step] = it->second;
494 DUNE_THROW(
NotInCache,
"MultiStepCache::getComposedJacobian(): The "
495 "requested composed Jacobian for step " << step <<
" is "
507 const shared_ptr<const Matrix> &jacobian)
509 if(!policy->cacheComposedJacobian(step))
511 if(!composedJacobians.insert(std::make_pair(step, jacobian)).second)
512 DUNE_THROW(
AlreadyInCache,
"MultiStepCache::setComposedJacobian(): "
513 "Composed Jacobian for time step " << step <<
" is "
514 "already in the cache!");
531 shared_ptr<const VectorU>
533 UnknownIterator it = unknowns.find(step);
534 if(it != unknowns.end())
536 DUNE_THROW(
NotInCache,
"Unknowns u_" << step <<
" missing in the "
547 void setUnknowns(Step step,
const shared_ptr<const VectorU> &unknowns_) {
549 DUNE_THROW(
AlreadyInCache,
"Unknowns u_" << step <<
" are already "
551 unknowns[step] = unknowns_;
565 residualValues.clear();
567 zeroResiduals.clear();
568 composedJacobians.clear();
593 void preStep(Step step, Step stepsOfScheme, Time endTime, Time dt) {
594 policy->preStep(step, stepsOfScheme, endTime, dt);
597 for(std::size_t order = 0; order < residualValues.size(); ++order) {
598 const ResidualIterator end = residualValues[order].end();
599 ResidualIterator it = residualValues[order].begin();
601 if(policy->canEvictResidualValue(order, it->first))
602 residualValues[order].erase(it++);
608 for(std::size_t order = 0; order < jacobians.size(); ++order) {
609 const MatrixIterator end = jacobians[order].end();
610 MatrixIterator it = jacobians[order].begin();
612 if(policy->canEvictJacobian(order, it->first))
613 jacobians[order].erase(it++);
619 for(std::size_t order = 0; order < zeroResiduals.size(); ++order) {
620 const ResidualIterator end = zeroResiduals[order].end();
621 ResidualIterator it = zeroResiduals[order].begin();
623 if(policy->canEvictZeroResidual(order, it->first))
624 zeroResiduals[order].erase(it++);
631 const MatrixIterator end = composedJacobians.end();
632 MatrixIterator it = composedJacobians.begin();
634 if(policy->canEvictComposedJacobian(it->first))
635 composedJacobians.erase(it++);
642 const UnknownIterator end = unknowns.end();
643 UnknownIterator it = unknowns.begin();
645 if(policy->canEvictUnknowns(it->first))
646 unknowns.erase(it++);
666 #endif // DUNE_PDELAB_MULTISTEP_CACHE_HH
Step stepsOfScheme
The number of steps in the scheme, set by preStep().
Definition: cache.hh:54
void setJacobian(std::size_t order, Step step, const shared_ptr< const Matrix > &jacobian)
store a Jacobian in the cache
Definition: cache.hh:388
virtual bool canEvictUnknowns(Step step) const
determine whether a vector of unknowns can be removed
Definition: cache.hh:184
void setResidualValue(std::size_t order, Step step, const shared_ptr< const VectorV > &residualValue)
store a residual value in the cache
Definition: cache.hh:327
void setZeroResidual(std::size_t order, Step step, const shared_ptr< const VectorV > &zeroResidual)
store a zero-residual in the cache
Definition: cache.hh:449
virtual bool canEvictJacobian(std::size_t order, Step step) const
determine whether a Jacobian can be removed
Definition: cache.hh:177
void setUnknowns(Step step, const shared_ptr< const VectorU > &unknowns_)
store a vector of unknowns in the cache
Definition: cache.hh:547
virtual bool canReuseZeroResidual(std::size_t order, Step requested, Step available) const
Whether a zero-residual can be reused.
Definition: cache.hh:125
Time dt
Time step size.
Definition: cache.hh:75
virtual bool canReuseJacobian(std::size_t order, Step requested, Step available) const
Whether a jacobian can be reused.
Definition: cache.hh:120
virtual bool cacheResidualValue(std::size_t order, Step step) const
whether to cache residual evaluations
Definition: cache.hh:86
Step currentStep
The number of the step set via preStep().
Definition: cache.hh:61
shared_ptr< const Matrix > getJacobian(std::size_t order, Step step)
get a Jacobian from the cache
Definition: cache.hh:361
virtual bool hasPureLinearAlpha(std::size_t order, Step step) const
Whether the operator has purely linear alpha_*() methods.
Definition: cache.hh:116
void postStep()
Do some housekeeping after computing a time-step.
Definition: cache.hh:656
virtual bool canEvictZeroResidual(std::size_t order, Step step) const
determine whether a zero-residual can be removed
Definition: cache.hh:180
Policy class for the MultiStepCache.
Definition: cache.hh:46
virtual void postStep()
called after the new values for a time step have been computed
Definition: cache.hh:165
virtual bool cacheZeroResidual(std::size_t order, Step step) const
whether to cache zero-residuals (only relevant for affine operators)
Definition: cache.hh:94
virtual bool isComposedAffine(Step step) const
Whether the composed operator is affine.
Definition: cache.hh:112
virtual bool canEvictComposedJacobian(Step step) const
determine whether a composed Jacobian can be removed
Definition: cache.hh:187
shared_ptr< const Policy > getPolicy() const
Definition: cache.hh:277
void flushAll()
Flush all cached values.
Definition: cache.hh:564
virtual bool cacheComposedJacobian(Step step) const
whether to cache composed jacobians (only relevant for affine operators)
Definition: cache.hh:98
virtual bool canEvictResidualValue(std::size_t order, Step step) const
determine whether a residual value can be removed
Definition: cache.hh:173
virtual ~MultiStepCachePolicy()
any virtual base needs a virtual destructor
Definition: cache.hh:79
Exception thrown when a stored item is already in the cache.
Definition: cache.hh:28
Exception thrown when a requested item is not in the cache.
Definition: cache.hh:25
Time endTime
The end-time of the step set via preStep().
Definition: cache.hh:68
Cache for the CachedMultiStepGridOperatorSpace.
Definition: cache.hh:234
virtual bool isAffine(std::size_t order, Step step) const
Whether the operator is affine.
Definition: cache.hh:108
shared_ptr< const VectorU > getUnknowns(Step step) const
get vector of unknowns from the cache
Definition: cache.hh:532
shared_ptr< const VectorV > getResidualValue(std::size_t order, Step step) const
get residual value from the cache
Definition: cache.hh:306
Base exception for cache related errors.
Definition: cache.hh:22
void setPolicy(const shared_ptr< Policy > &policy_)
Definition: cache.hh:279
shared_ptr< Policy > getPolicy()
Definition: cache.hh:275
void preStep(Step step, Step stepsOfScheme, Time endTime, Time dt)
Do some housekeeping before computing a time-step.
Definition: cache.hh:593
virtual bool canReuseComposedJacobian(Step requested, Step available) const
Whether a composed jacobian can be reused.
Definition: cache.hh:130
shared_ptr< const Matrix > getComposedJacobian(Step step)
get a composed Jacobian from the cache
Definition: cache.hh:482
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
MultiStepCachePolicy< Step, Time > Policy
Definition: cache.hh:235
void setComposedJacobian(Step step, const shared_ptr< const Matrix > &jacobian)
store a composed Jacobian in the cache
Definition: cache.hh:506
virtual void preStep(Step step, Step stepsOfScheme_, Time endTime_, Time dt_)
method called before starting to compute a step
Definition: cache.hh:155
virtual bool cacheJacobian(std::size_t order, Step step) const
whether to cache jacobians (only relevant for affine operators)
Definition: cache.hh:90
MultiStepCache(const shared_ptr< Policy > &policy_=shared_ptr< Policy >(new Policy))
Definition: cache.hh:266
shared_ptr< const VectorV > getZeroResidual(std::size_t order, Step step) const
get zero-residual from the cache
Definition: cache.hh:421