dune-localfunctions  2.2.0
monomlocalinterpolation.hh
Go to the documentation of this file.
00001 // -*- tab-width: 4; indent-tabs-mode: nil -*-
00002 #ifndef DUNE_MONOMLOCALINTERPOLATION_HH
00003 #define DUNE_MONOMLOCALINTERPOLATION_HH
00004 
00005 #include <vector>
00006 
00007 #include <dune/common/deprecated.hh>
00008 #include <dune/common/fvector.hh>
00009 #include <dune/common/fmatrix.hh>
00010 
00011 #include <dune/geometry/type.hh>
00012 #include <dune/geometry/quadraturerules.hh>
00013 
00014 namespace Dune 
00015 {
00016 
00017   template<class LB, unsigned int size>
00018   class MonomLocalInterpolation 
00019   {
00020     typedef typename LB::Traits::DomainType D;
00021     typedef typename LB::Traits::DomainFieldType DF;
00022     static const int dimD=LB::Traits::dimDomain;
00023     typedef typename LB::Traits::RangeType R;
00024     typedef typename LB::Traits::RangeFieldType RF;
00025 
00026     typedef QuadratureRule<DF,dimD> QR;
00027     typedef typename QR::iterator QRiterator;
00028 
00029     void init() {
00030       if(size != lb.size())
00031         DUNE_THROW(Exception, "size template parameter does not match size of "
00032                    "local basis");
00033 
00034       const QRiterator qrend = qr.end();
00035       for(QRiterator qrit = qr.begin(); qrit != qrend; ++qrit) {
00036         std::vector<R> base;
00037         lb.evaluateFunction(qrit->position(),base);
00038 
00039         for(unsigned int i = 0; i < size; ++i)
00040           for(unsigned int j = 0; j < size; ++j)
00041             Minv[i][j] += qrit->weight() * base[i] * base[j];
00042       }
00043       Minv.invert();
00044     }
00045 
00046   public:
00047         MonomLocalInterpolation (const GeometryType::BasicType &bt_,
00048                              const LB &lb_) DUNE_DEPRECATED
00049       : gt(bt_, dimD), lb(lb_), Minv(0)
00050       , qr(QuadratureRules<DF,dimD>::rule(gt, 2*lb.order()))
00051     { init(); }
00052 
00053     MonomLocalInterpolation (const GeometryType &gt_,
00054                              const LB &lb_)
00055       : gt(gt_), lb(lb_), Minv(0)
00056       , qr(QuadratureRules<DF,dimD>::rule(gt, 2*lb.order()))
00057     { init(); }
00058 
00060         template<typename F, typename C>
00061         void interpolate (const F& f, std::vector<C>& out) const
00062         {
00063           out.clear();
00064       out.resize(size, 0);
00065       
00066       const QRiterator qrend = qr.end();
00067       for(QRiterator qrit = qr.begin(); qrit != qrend; ++qrit) {
00068         //TODO: mass matrix
00069         R y;
00070         f.evaluate(qrit->position(),y);
00071 
00072         std::vector<R> base;
00073         lb.evaluateFunction(qrit->position(),base);
00074 
00075         for(unsigned int i = 0; i < size; ++i)
00076           for(unsigned int j = 0; j < size; ++j)
00077             out[i] += Minv[i][j] * qrit->weight() * y * base[j];
00078       }
00079         }
00080 
00081   private:
00082     GeometryType gt;
00083     const LB &lb;
00084     FieldMatrix<RF, size, size> Minv;
00085     const QR &qr;
00086   };
00087 
00088 }
00089 
00090 #endif //DUNE_MONOMLOCALINTERPOLATION_HH