C++ Wynn's epsilon algorithm with caching

out.txt

//note this can work with complex numbers, vectors, or tensors as the series data... not just doubles
we3 3.16667
we4 3.13334

//if the map caching debug output were turned on, the output would be:
we3 caching 1.25
caching -0.750002
caching 3.16667
3.16667
we4 caching -1.75
caching 3.13334
3.13334

WynnEpsilon.cpp

//#include <iostream> //only to debug cache
#include <map>

/**
 * Functor that returns an estimate for the converged
 * value of a series of numbers or higher order math
 * objects like tensors or vectors.
 *
 * Wynn's epsilon method is used, and the intermediate
 * epsilon results are cached in a two level map. A 
 * reference to the series container supplied to
 * the constructor is kept. If the caller puts more
 * objects in the container and then calls the
 * functor again, the cached values are reused while
 * calculating the updated convergence estimate.
 */
template<class S>
struct WynnEpsilon{
 /**
  * ctor callers must supply a container that
  * has the value_type typedef and that supports
  * operator[]
  *
  * @tparam S container type
  * @param s container of math objects
  */
 WynnEpsilon(S const&s): s_(s) {} //ctor
 /**
  * Returns series convergence estimate.
  */
 auto operator()(){ //call functor through this
  int sz(s_.size()); //signed size
  int off(!(sz % 2)); //1 even, 0 odd
  return w(sz-2-off,off);
 }
private:
 //https://mathematica.stackexchange.com/questions/178536/
 ///Wynn's epsilon algorithm, with caching
 typename S::value_type w(int r,int n){
  switch(r){ //cases we don't cache in map m
   case -2: return 0.0; //always 0
   case -1: return s_[n]; //pulls from s
  }
  auto& mr(m_[r]); //reference to inner map at key r
  auto mrn(mr.find(n)); //inner map iterator for key n
  if(mrn!=mr.end()) //if key n exists
   return mrn->second; //return cached result
  //if not cached, calculate
  auto mrns(w(r-2,n+1)+1.0/(w(r-1,n+1)-w(r-1,n)));
  //std::cerr<<"caching "<<mrns<<"\n"; //debug caching
  mr.emplace(std::pair{n,mrns}); //then cache
  return mrns; //then return
 }
 S const& s_; //container reference
 std::map<int,std::map<int,typename S::value_type>> m_; //epsilon cache
};

#include<iostream>
#include<vector>

int main(int,char*[]){
 std::vector<double> v{4.,2.66667,3.46667};
 WynnEpsilon we{v};
 std::cerr<<"we3 "<<we()<<"\n";
 v.push_back(2.89524);
 std::cerr<<"we4 "<<we()<<"\n";
}

WynnEpsilonTest.cpp

//not needed for the other files to work
//just an example of how testing would be done
//using my Expect gist and an equal function

#include "Math/WynnEpsilon/WynnEpsilon.hpp"

#include "Math/Equal/Equal.hpp"
#include "Expect/Expect.hpp"

#include <vector>

int main(){
 Expect e("WynnEpsilonTest");
 {std::vector<double> v{4.,8./3,52./15};
 WynnEpsilon we{v,true};
 e(equal(we(),19./6),"Pi3Terms");
 v.push_back(304./105);
 e(equal(we(),47./15),"Pi4Terms");}
 {std::vector<double> v{1.,2.,9./4};
 WynnEpsilon we{v,false};
 e(equal(we(),7./3),"E3Terms");
 v.push_back(64./27);
 e(equal(we(),139./56),"E4Terms");}
 return e.result();
}