dune-pdelab  2.0.0
cache.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=8 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_MULTISTEP_CACHE_HH
5 #define DUNE_PDELAB_MULTISTEP_CACHE_HH
6 
7 #include <cstddef>
8 #include <map>
9 #include <utility>
10 #include <vector>
11 
12 #include <dune/common/exceptions.hh>
13 #include <dune/common/shared_ptr.hh>
14 
15 namespace Dune {
16  namespace PDELab {
17 
20 
22  class CacheError : public Exception {};
23 
25  class NotInCache : public CacheError {};
26 
28  class AlreadyInCache : public CacheError {};
29 
31 
45  template<class Step = int, class Time = double>
47  protected:
49 
56 
63 
68  Time endTime;
70 
75  Time dt;
76 
77  public:
79  virtual ~MultiStepCachePolicy() {}
80 
83 
85 
86  virtual bool cacheResidualValue(std::size_t order, Step step) const
87  { return true; }
89 
90  virtual bool cacheJacobian(std::size_t order, Step step) const
91  { return true; }
93 
94  virtual bool cacheZeroResidual(std::size_t order, Step step) const
95  { return true; }
97 
98  virtual bool cacheComposedJacobian(Step step) const
99  { return true; }
100 
102 
105 
107 
108  virtual bool isAffine(std::size_t order, Step step) const
109  { return false; }
111 
112  virtual bool isComposedAffine(Step step) const
113  { return false; }
115 
116  virtual bool hasPureLinearAlpha(std::size_t order, Step step) const
117  { return false; }
119 
120  virtual bool canReuseJacobian(std::size_t order,
121  Step requested, Step available) const
122  { return false; }
124 
125  virtual bool canReuseZeroResidual(std::size_t order,
126  Step requested, Step available) const
127  { return false; }
129 
130  virtual bool canReuseComposedJacobian(Step requested,
131  Step available) const
132  { return false; }
133 
135 
138 
140 
155  virtual void preStep(Step step, Step stepsOfScheme_,
156  Time endTime_, Time dt_)
157  {
158  currentStep = step;
159  stepsOfScheme = stepsOfScheme_;
160  endTime = endTime_;
161  dt = dt_;
162  }
164 
165  virtual void postStep() { }
166 
168 
171 
173  virtual bool canEvictResidualValue(std::size_t order,
174  Step step) const
175  { return step + stepsOfScheme < currentStep; }
177  virtual bool canEvictJacobian(std::size_t order, Step step) const
178  { return step + stepsOfScheme < currentStep; }
180  virtual bool canEvictZeroResidual(std::size_t order,
181  Step step) const
182  { return step + stepsOfScheme < currentStep; }
184  virtual bool canEvictUnknowns(Step step) const
185  { return step + stepsOfScheme < currentStep; }
187  virtual bool canEvictComposedJacobian(Step step) const
188  { return step + stepsOfScheme < currentStep; }
189  };
190 
192 
232  template<class VectorU, class VectorV, class Matrix,
233  class Step = int, class Time = double>
234  struct MultiStepCache {
236 
237  private:
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;
244 
245  // non-linear caching across time steps
246  std::vector<ResidualMap> residualValues;
247 
248  // affine operators
249  std::vector<MatrixMap> jacobians;
250  std::vector<ResidualMap> zeroResiduals;
251 
252  // composed jacobians
253  MatrixMap composedJacobians;
254 
255  // old values of the unknowns
256  UnknownMap unknowns;
257 
258  // policy object
259  shared_ptr<Policy> policy;
260 
261  public:
262 
265 
266  MultiStepCache(const shared_ptr<Policy> &policy_ =
267  shared_ptr<Policy>(new Policy)) :
268  policy(policy_)
269  {
270  if(!policy)
271  DUNE_THROW(CacheError,
272  "MultiStepCache constructed with policy == NULL");
273  }
274 
275  shared_ptr<Policy> getPolicy()
276  { return policy; }
277  shared_ptr<const Policy> getPolicy() const
278  { return policy; }
279  void setPolicy(const shared_ptr<Policy>& policy_)
280  {
281  if(!policy_)
282  DUNE_THROW(CacheError, "MultiStepCache::setPolicy(): attempt to set "
283  "policy == NULL");
284  policy = policy_;
285  }
286 
288 
290 
293 
295 
305  shared_ptr<const VectorV>
306  getResidualValue(std::size_t order, Step step) const {
307  if(order < residualValues.size()) {
308  ResidualIterator it = residualValues[order].find(step);
309  if(it != residualValues[order].end())
310  return it->second;
311  }
312  DUNE_THROW(NotInCache, "MultiStepCache::getResidualValue(): The "
313  "requested residual value "
314  "r_" << order << "(t_" << step << ",u_" << step << ") is "
315  "not in the cache");
316  }
318 
327  void setResidualValue(std::size_t order, Step step,
328  const shared_ptr<const VectorV> &residualValue)
329  {
330  if(!policy->cacheResidualValue(order, step))
331  return;
332  if(order >= residualValues.size())
333  residualValues.resize(order+1);
334  if(!residualValues[order].insert(std::make_pair(step, residualValue))
335  .second)
336  DUNE_THROW(AlreadyInCache, "Residual value"
337  "r_" << order << "(t_" << step << ", u_" << step << ") "
338  "is already in the cache!");
339  }
340 
342 
345 
347 
360  shared_ptr<const Matrix>
361  getJacobian(std::size_t order, Step step) {
362  if(order < jacobians.size()) {
363  MatrixIterator it = jacobians[order].find(step);
364  const MatrixIterator &end = jacobians[order].end();
365  if(it != end)
366  return it->second;
367 
368  // try to copy from another step
369  for(it = jacobians[order].begin(); it != end; ++it)
370  if(policy->canReuseJacobian(order, step, it->first))
371  // assign and return value
372  return jacobians[order][step] = it->second;
373  }
374  DUNE_THROW(NotInCache, "MultiStepCache::getJacobian(): The requested "
375  "Jacobian J(r_" << order << "|_t_" << step << ") is not in "
376  "the cache");
377  }
379 
388  void setJacobian(std::size_t order, Step step,
389  const shared_ptr<const Matrix> &jacobian)
390  {
391  if(!policy->cacheJacobian(order, step))
392  return;
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 "
398  "the cache!");
399  }
400 
402 
405 
407 
420  shared_ptr<const VectorV>
421  getZeroResidual(std::size_t order, Step step) const {
422  if(order < zeroResiduals.size()) {
423  ResidualIterator it = zeroResiduals[order].find(step);
424  const ResidualIterator &end = zeroResiduals[order].end();
425  if(it != end)
426  return it->second;
427 
428  // try to copy from another step
429  for(it = zeroResiduals[order].begin(); it != end; ++it)
430  if(policy->canReuseZeroResidual(order, step, it->first))
431  // assign and return value
432  return zeroResiduals[order][step] = it->second;
433  }
434  DUNE_THROW(NotInCache, "MultiStepCache::getZeroResidual(): The "
435  "requested zero-residual "
436  "r_" << order << "(t_" << step << ",0) is not in the "
437  "cache");
438  }
440 
449  void setZeroResidual(std::size_t order, Step step,
450  const shared_ptr<const VectorV> &zeroResidual)
451  {
452  if(!policy->cacheZeroResidual(order, step))
453  return;
454  if(order >= zeroResiduals.size())
455  zeroResiduals.resize(order+1);
456  if(!zeroResiduals[order].insert(std::make_pair(step, zeroResidual))
457  .second)
458  DUNE_THROW(AlreadyInCache, "Zero-residual "
459  "r_" << order << "(t_" << step << ", 0) is already in "
460  "the cache!");
461  }
462 
464 
467 
469 
481  shared_ptr<const Matrix>
482  getComposedJacobian(Step step) {
483  MatrixIterator it = composedJacobians.find(step);
484  const MatrixIterator &end = composedJacobians.end();
485  if(it != end)
486  return it->second;
487 
488  // try to copy from another step
489  for(it = composedJacobians.begin(); it != end; ++it)
490  if(policy->canReuseComposedJacobian(step, it->first))
491  // assign and return value
492  return composedJacobians[step] = it->second;
493 
494  DUNE_THROW(NotInCache, "MultiStepCache::getComposedJacobian(): The "
495  "requested composed Jacobian for step " << step << " is "
496  "not in the cache");
497  }
499 
506  void setComposedJacobian(Step step,
507  const shared_ptr<const Matrix> &jacobian)
508  {
509  if(!policy->cacheComposedJacobian(step))
510  return;
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!");
515  }
516 
518 
521 
523 
531  shared_ptr<const VectorU>
532  getUnknowns(Step step) const {
533  UnknownIterator it = unknowns.find(step);
534  if(it != unknowns.end())
535  return it->second;
536  DUNE_THROW(NotInCache, "Unknowns u_" << step << " missing in the "
537  "cache!");
538  }
540 
547  void setUnknowns(Step step, const shared_ptr<const VectorU> &unknowns_) {
548  if(unknowns[step])
549  DUNE_THROW(AlreadyInCache, "Unknowns u_" << step << " are already "
550  "in the cache!");
551  unknowns[step] = unknowns_;
552  }
553 
555 
558 
560 
564  void flushAll() {
565  residualValues.clear();
566  jacobians.clear();
567  zeroResiduals.clear();
568  composedJacobians.clear();
569  unknowns.clear();
570  }
571 
573 
575 
577 
593  void preStep(Step step, Step stepsOfScheme, Time endTime, Time dt) {
594  policy->preStep(step, stepsOfScheme, endTime, dt);
595 
596  // residual values
597  for(std::size_t order = 0; order < residualValues.size(); ++order) {
598  const ResidualIterator end = residualValues[order].end();
599  ResidualIterator it = residualValues[order].begin();
600  while(it != end)
601  if(policy->canEvictResidualValue(order, it->first))
602  residualValues[order].erase(it++);
603  else
604  ++it;
605  }
606 
607  // Jacobians
608  for(std::size_t order = 0; order < jacobians.size(); ++order) {
609  const MatrixIterator end = jacobians[order].end();
610  MatrixIterator it = jacobians[order].begin();
611  while(it != end)
612  if(policy->canEvictJacobian(order, it->first))
613  jacobians[order].erase(it++);
614  else
615  ++it;
616  }
617 
618  // zero-residual
619  for(std::size_t order = 0; order < zeroResiduals.size(); ++order) {
620  const ResidualIterator end = zeroResiduals[order].end();
621  ResidualIterator it = zeroResiduals[order].begin();
622  while(it != end)
623  if(policy->canEvictZeroResidual(order, it->first))
624  zeroResiduals[order].erase(it++);
625  else
626  ++it;
627  }
628 
629  // composed Jacobians
630  {
631  const MatrixIterator end = composedJacobians.end();
632  MatrixIterator it = composedJacobians.begin();
633  while(it != end)
634  if(policy->canEvictComposedJacobian(it->first))
635  composedJacobians.erase(it++);
636  else
637  ++it;
638  }
639 
640  // unknowns
641  {
642  const UnknownIterator end = unknowns.end();
643  UnknownIterator it = unknowns.begin();
644  while(it != end)
645  if(policy->canEvictUnknowns(it->first))
646  unknowns.erase(it++);
647  else
648  ++it;
649  }
650  }
651 
653 
656  void postStep() { policy->postStep(); }
657 
659 
660  };
661 
663  } // namespace PDELab
664 } // namespace Dune
665 
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