bvp::interpolator< RBF > Class Template Reference

#include <interpolator.hpp>

Inheritance diagram for bvp::interpolator< RBF >:

Inheritance graph
[legend]
Collaboration diagram for bvp::interpolator< RBF >:

Collaboration graph
[legend]

List of all members.

Public Member Functions

Constructors
Constructors that take interpolation data perform the interpolation as part of the initialisation.

 interpolator ()
 Does not initialise the interpolator.
 interpolator (shared_ptr< linear_BVP2 > bvp)
 Interpolate given a BVP.
 interpolator (const map< point, double > &Xi)
 Interpolate given some data points and the value at those points.
 interpolator (shared_ptr< domain > Omega, const map< point, double > &Xi)
 Interpolate given a domain, data points and the value at those points.
Interpolation
Interpolate again either given new data or a different BVP.

void interpolate (const map< point, double > &Xi)
void interpolate (shared_ptr< linear_BVP2 > bvp)
Evaluations and derivatives
double operator() (const point &p) const
 Evaluation.
double at (const point &p) const
 Evaluation.
double d (const point &p, size_t k) const
 First derivative.
double d2 (const point &p, size_t k1, size_t k2) const
 Second derivatives.
Partial redefinitions
These functions allow for partial redefinition of the BVP as equired for the additive Schwartz domain decomposition method, a nd for other methods.

void set_f (const realfunc &f)
void set_g (const realfunc &g)
void set_f (const map< point, double > &f)
void set_g (const map< point, double > &g)
Linear arithmetic operators
These functions return a new interpolator.

interpolator< RBF > operator+ (const interpolator< RBF > &u) const
 Needs two operators on the same domain.
interpolator< RBF > operator- (const interpolator< RBF > &u) const
 Needs two operators on the same domain.
interpolator< RBF > operator* (double a) const
interpolator< RBF > operator/ (double a) const

Private Types

typedef std::pair
< linalg::point, std::vector
< size_t > > 
diff_data

Private Member Functions

void computecoeffs ()
void init (shared_ptr< linear_BVP2 > bvp)
void not_initted (int line, string file) const
size_t hash_value (const std::vector< RBF > &rbfs_in)

Private Attributes

shared_ptr< linear_BVP2thebvp
size_t n
size_t m
matrix M
bool initted
linalg::vector coeffs
std::vector< RBF > rbfs
size_t rbfs_hash
map< diff_data, double > remtable

template<typename RBF>
class bvp::interpolator< RBF >


Member Typedef Documentation

template<typename RBF>
typedef std::pair<linalg::point, std::vector<size_t> > bvp::interpolator< RBF >::diff_data [private]


Constructor & Destructor Documentation

template<typename RBF>
bvp::interpolator< RBF >::interpolator (  )  [inline]

Does not initialise the interpolator.

00020   {
00021     initted = false;
00022   }

template<typename RBF>
bvp::interpolator< RBF >::interpolator ( shared_ptr< linear_BVP2 bvp  )  [inline]

Interpolate given a BVP.

00026   {
00027     //Workaround because gdb can't break inside constructors. :-(
00028     init(bvp);
00029   }

Here is the call graph for this function:

template<typename RBF>
bvp::interpolator< RBF >::interpolator ( const map< point, double > &  Xi  )  [inline]

Interpolate given some data points and the value at those points.

00081   {
00082     interpolate(Xi);
00083   }

Here is the call graph for this function:

template<typename RBF>
bvp::interpolator< RBF >::interpolator ( shared_ptr< domain Omega,
const map< point, double > &  Xi 
) [inline]

Interpolate given a domain, data points and the value at those points.

Must provide domain information. The values of Xi must match points on the given domain Omega.

00034   {
00035     //Check that Xi matches the domain in size
00036     size_t Omega_size = Omega -> get_interior().size() 
00037                        + Omega -> get_boundary().size();
00038     if(Xi.size() != Omega_size)
00039     {
00040       badArgument exc;
00041       if(Xi.size() < Omega_size)
00042         exc.reason = "Did not provide enough interpolation data for every point in the given domain.";
00043       else
00044         exc.reason = "Provided more interpolation data than points in the given domain.";
00045       exc.line = __LINE__;
00046       exc.file = __FILE__;
00047       throw exc;
00048     }
00049 
00050     //Create a trivial bvp taking into account the given domain.
00051     shared_ptr<Id_op> Id(new Id_op);
00052     shared_ptr<dirichlet_op> D(new dirichlet_op);
00053     map<point, double> f, g;
00054 
00055     //Extract f and g from given information
00056     for(map<point, double>::const_iterator I = Xi.begin();
00057         I != Xi.end(); I++)
00058     {
00059       if(utils::contains(Omega -> get_interior(), I -> first))
00060         f[I -> first] = I -> second;
00061       else if(utils::contains(Omega -> get_boundary(), I -> first))
00062         g[I -> first] = I -> second;
00063       else
00064       {
00065         badArgument exc;
00066         exc.reason = "The interpolation data contains points not in the given domain.";
00067         exc.line = __LINE__;
00068         exc.file = __FILE__;
00069         throw exc;
00070       }
00071     } 
00072     
00073     shared_ptr<linear_BVP2> new_bvp(new linear_BVP2(Omega, Id, D, f,g));
00074     
00075     interpolate(new_bvp);
00076 
00077   }

Here is the call graph for this function:


Member Function Documentation

template<typename RBF>
void bvp::interpolator< RBF >::interpolate ( const map< point, double > &  Xi  )  [inline]

00087   {
00088     
00089     if(Xi.empty()){//Dude, wtf?
00090       badArgument exc;
00091       exc.reason = "Cannot interpolate if no data is given.";
00092       exc.line = __LINE__;
00093       exc.file = __FILE__;
00094       throw exc;
00095     }
00096 
00097     //Create a trivial bvp.
00098     shared_ptr<Id_op> Id(new Id_op);
00099     shared_ptr<dirichlet_op> D(new dirichlet_op);
00100     set<point> intr;
00101     set<point> bdry; //empty
00102     map<point, point> nrml; //empty
00103     map<point, double> g; //empty
00104     
00105     bool dim_set = false;
00106     size_t dimension = 0;
00107     for(map<point,double>::const_iterator i = Xi.begin(); i != Xi.end(); i++){
00108       if(!dim_set){
00109         dimension = (i -> first).size();
00110         dim_set = true;
00111       }
00112       else if(dimension != (i -> first).size()){
00113         badArgument exc;
00114         exc.reason = "Inconformant dimensions in interpolation data.";
00115         exc.line = __LINE__;
00116         exc.file = __FILE__;
00117         throw exc;
00118       }
00119       intr.insert( i->first);
00120     }
00121     shared_ptr<domain> Omega(new domain(dimension, intr,bdry,nrml));
00122     shared_ptr<linear_BVP2> bvp(new linear_BVP2(Omega, Id, D, Xi, g));
00123     
00124     init(bvp);
00125   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::interpolate ( shared_ptr< linear_BVP2 bvp  )  [inline]

00129   {
00130     init(bvp);
00131   }

Here is the call graph for this function:

template<typename RBF>
double bvp::interpolator< RBF >::operator() ( const point &  p  )  const [inline]

Evaluation.

Reimplemented from bvp::realfunc.

00208                                                           {
00209     return at(p);
00210   }

Here is the call graph for this function:

template<typename RBF>
double bvp::interpolator< RBF >::at ( const point &  p  )  const [inline, virtual]

Evaluation.

Reimplemented from bvp::realfunc.

00214                                                   {
00215     if(!initted){
00216       not_initted(__LINE__, __FILE__);
00217     }
00218     std::vector<size_t> alpha; //empty vector
00219     diff_data loc = std::make_pair(p,alpha);
00220     if(remtable.find(loc) != remtable.end())
00221       return remtable[loc];
00222     
00223     double result = 0;
00224     for(size_t i = 1; i <= coeffs.size(); i++)
00225       result += coeffs(i)*rbfs[i-1].at(p);
00226     
00227     remtable[loc] = result;
00228     return result;
00229   }

Here is the call graph for this function:

template<typename RBF>
double bvp::interpolator< RBF >::d ( const point &  p,
size_t  k 
) const [inline, virtual]

First derivative.

Reimplemented from bvp::realfunc.

00232                                                            {
00233     if(!initted){
00234       not_initted(__LINE__, __FILE__);
00235     }
00236     std::vector<size_t> alpha(k); alpha[k-1]++;
00237     diff_data loc = std::make_pair(p,alpha);
00238     if(remtable.find(loc) != remtable.end())
00239       return remtable[loc];
00240 
00241     double result = 0;
00242     for(size_t i = 1; i <= coeffs.size(); i++)
00243       result += coeffs(i)*rbfs[i-1].d(p,k);
00244 
00245     remtable[loc] = result;
00246     return result;
00247   }

Here is the call graph for this function:

template<typename RBF>
double bvp::interpolator< RBF >::d2 ( const point &  p,
size_t  k1,
size_t  k2 
) const [inline, virtual]

Second derivatives.

Reimplemented from bvp::realfunc.

00250                                                                         {
00251     if(!initted){
00252       not_initted(__LINE__, __FILE__);
00253     }
00254     std::vector<size_t> alpha(k1>k2?k1:k2); alpha[k1-1]++; alpha[k2-1]++;
00255     diff_data loc = std::make_pair(p,alpha);
00256     if(remtable.find(loc) != remtable.end())
00257       return remtable[loc];
00258 
00259     double result = 0;
00260     for(size_t i = 1; i <= coeffs.size(); i++)
00261       result += coeffs(i)*rbfs[i-1].d2(p,k1,k2);
00262     return result;
00263   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::set_f ( const realfunc f  )  [inline]

00332                                                 {
00333     if(!initted){
00334       not_initted(__LINE__, __FILE__);
00335     }
00336     thebvp -> set_f(f);
00337     computecoeffs();
00338   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::set_g ( const realfunc g  )  [inline]

00341                                                 {
00342     if(!initted){
00343       not_initted(__LINE__, __FILE__);
00344     }
00345     thebvp -> set_g(g);
00346     computecoeffs();
00347   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::set_f ( const map< point, double > &  f  )  [inline]

00350                                                           {
00351     if(!initted){
00352       not_initted(__LINE__, __FILE__);
00353     }
00354     thebvp -> set_f(f);
00355     computecoeffs();
00356   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::set_g ( const map< point, double > &  g  )  [inline]

00359                                                           {
00360     if(!initted){
00361       not_initted(__LINE__, __FILE__);
00362     }
00363     thebvp -> set_g(g);
00364     computecoeffs();
00365   }

Here is the call graph for this function:

template<typename RBF>
interpolator< RBF > bvp::interpolator< RBF >::operator+ ( const interpolator< RBF > &  u  )  const [inline]

Needs two operators on the same domain.

00270   {    
00271     if(this -> rbfs_hash != u.rbfs_hash){
00272       badArgument exc;
00273       exc.reason = 
00274         "Cannot add interpolators on different domains (interior and boundary must match)."; 
00275       exc.line = __LINE__;
00276       exc.file = __FILE__;
00277       throw exc;        
00278     }
00279     
00280     interpolator<RBF> out = *this;
00281     out.coeffs = (this -> coeffs) + u.coeffs;
00282    
00283     out.remtable.clear();
00284     
00285     return out;
00286   }

template<typename RBF>
interpolator< RBF > bvp::interpolator< RBF >::operator- ( const interpolator< RBF > &  u  )  const [inline]

Needs two operators on the same domain.

00291   {    
00292     if(this -> rbfs_hash != u.rbfs_hash){
00293       badArgument exc;
00294       exc.reason = 
00295         "Cannot subtract interpolators on different domains (interior and boundary must match)."; 
00296       exc.line = __LINE__;
00297       exc.file = __FILE__;
00298       throw exc;        
00299     }
00300 
00301     interpolator<RBF> out = *this;
00302     out.coeffs = (this -> coeffs) - u.coeffs;
00303     
00304     out.remtable.clear();
00305     
00306     return out;
00307   }

template<typename RBF>
interpolator< RBF > bvp::interpolator< RBF >::operator* ( double  a  )  const [inline]

00312   {
00313     interpolator<RBF> u = *this;
00314     u.coeffs = (this -> coeffs)*a;
00315     u.remtable.clear();
00316     return u;
00317   }

template<typename RBF>
interpolator< RBF > bvp::interpolator< RBF >::operator/ ( double  a  )  const [inline]

00322   {
00323     interpolator<RBF> u = *this;
00324     u.coeffs = (this -> coeffs)*(1/a);
00325     u.remtable.clear();
00326     return u;
00327   }

template<typename RBF>
void bvp::interpolator< RBF >::computecoeffs (  )  [inline, private]

00378                                        {
00379     using namespace std;
00380     linalg::vector rhs(n+m);
00381 
00382     map<point, double>::const_iterator I;
00383 
00384     I = (thebvp -> get_f()).begin();
00385     for(size_t i = 1; i <= n; i++){
00386       rhs(i) = I->second;
00387       I++;
00388     }
00389     I = (thebvp -> get_g()).begin();   
00390     for(size_t i = n+1; i <= n+m; i++){
00391       rhs(i) = I->second;
00392       I++;
00393     }
00394 
00395     coeffs = M.inv(rhs);
00396     remtable.clear();
00397   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::init ( shared_ptr< linear_BVP2 bvp  )  [inline, private]

00135   {
00136     thebvp = bvp;
00137 
00138     using namespace linalg;
00139     using std::set;
00140     
00141     shared_ptr<const domain> Omega = bvp -> get_domain();
00142     set<point> interior = Omega -> get_interior();
00143     set<point> boundary = Omega -> get_boundary();
00144     map<point, vector> normals = Omega -> get_normals();
00145     n = interior.size();
00146     m = boundary.size();
00147     
00148     vector temp(n+m);
00149     coeffs = temp;
00150     rbfs.reserve(n+m);
00151 
00152     
00153     RBF::set_dimension(Omega -> get_dimension());
00154     
00155     set<point>::iterator I;
00156     //Define all the rbfs...
00157     for(I = interior.begin(); I != interior.end(); I++){
00158       RBF r(*I);
00159       rbfs.push_back(r);
00160     }
00161     for(I = boundary.begin(); I != boundary.end(); I++){
00162       RBF r(*I);
00163       rbfs.push_back(r);
00164     }
00165     
00166     //Now define the matrix to be inverted...
00167     matrix Mtemp(n+m,n+m);
00168     M = Mtemp;
00169     shared_ptr<const linear_diff_op2> L = thebvp -> get_linear_diff_op2(); 
00170 
00171     shared_ptr<const bdry_diff_op> B = thebvp -> get_bdry_diff_op();
00172 
00173     I = interior.begin();
00174     for(size_t i = 1; i <= n; i++){
00175       for(size_t j = 1; j <= n+m; j++)
00176         M(i,j) = L -> at(rbfs[j-1], *I);
00177       I++;
00178     }
00179     
00180     map<point, vector>::iterator J;
00181     J = normals.begin();
00182     for(size_t i = n+1; i <= n+m; i++){
00183       for(size_t j = 1; j <= n+m; j++)
00184         M(i,j) = B -> at(rbfs[j-1], J->first, J->second);
00185       J++;
00186     }
00187    
00188     computecoeffs();
00189     initted = true;
00190     rbfs_hash = hash_value(rbfs);
00191   }

Here is the call graph for this function:

template<typename RBF>
void bvp::interpolator< RBF >::not_initted ( int  line,
string  file 
) const [inline, private]

00368                                                                 {
00369     badArgument exc;
00370     exc.reason = 
00371       "Interpolator can't interpolate without initialisation data.";
00372     exc.line = line;
00373     exc.file = file;
00374     throw exc;
00375   }

template<typename RBF>
size_t bvp::interpolator< RBF >::hash_value ( const std::vector< RBF > &  rbfs_in  )  [inline, private]

00195   {
00196     size_t seed = 0;  
00197 
00198     for(size_t i = 0; i < rbfs_in.size(); i++)
00199       boost::hash_combine(seed,rbfs[i]);
00200 
00201     return seed;
00202   }


Member Data Documentation

template<typename RBF>
shared_ptr<linear_BVP2> bvp::interpolator< RBF >::thebvp [private]

template<typename RBF>
size_t bvp::interpolator< RBF >::n [private]

template<typename RBF>
size_t bvp::interpolator< RBF >::m [private]

template<typename RBF>
matrix bvp::interpolator< RBF >::M [private]

template<typename RBF>
bool bvp::interpolator< RBF >::initted [private]

template<typename RBF>
linalg::vector bvp::interpolator< RBF >::coeffs [private]

template<typename RBF>
std::vector<RBF> bvp::interpolator< RBF >::rbfs [private]

template<typename RBF>
size_t bvp::interpolator< RBF >::rbfs_hash [private]

template<typename RBF>
map<diff_data, double> bvp::interpolator< RBF >::remtable [mutable, private]


The documentation for this class was generated from the following files:

Generated on Sat Jun 28 00:33:26 2008 by  doxygen 1.5.6