00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #ifndef _FRAMETL_SIMPLE_ELLIPTIC_EQUATION_H
00011 #define _FRAMETL_SIMPLE_ELLIPTIC_EQUATION_H
00012
00013 #include <aggregated_frame.h>
00014 #include <numerics/bvp.h>
00015 #include <adaptive/compression.h>
00016 #include <interval/i_index.h>
00017 #include <galerkin/infinite_preconditioner.h>
00018 #include <frame_support.h>
00019 #include <index1D.h>
00020
00021 using FrameTL::AggregatedFrame;
00022 using MathTL::EllipticBVP;
00023 using WaveletTL::CompressionStrategy;
00024 using WaveletTL::IntervalIndex;
00025 using WaveletTL::FullyDiagonalEnergyNormPreconditioner;
00026
00027 namespace FrameTL
00028 {
00029
00069 template <class IBASIS, unsigned int DIM>
00070 class SimpleEllipticEquation
00071
00072 : public FullyDiagonalEnergyNormPreconditioner<typename AggregatedFrame<IBASIS,DIM>::Index>
00073 {
00074 public:
00075
00084 SimpleEllipticEquation(const EllipticBVP<DIM>* ell_bvp,
00085 const AggregatedFrame<IBASIS,DIM>* frame,
00086 const int jmax);
00087
00091 typedef AggregatedFrame<IBASIS,DIM> Frame;
00092
00097 typedef AggregatedFrame<IBASIS,DIM> WaveletBasis;
00098
00102 typedef typename Frame::Index Index;
00103
00107 const AggregatedFrame<IBASIS,DIM>& frame() const { return *frame_; }
00108
00112 const EllipticBVP<DIM>& get_bvp() const { return *ell_bvp_; }
00113
00118 const AggregatedFrame<IBASIS,DIM>& basis() const { return *frame_; }
00119
00123 static const int space_dimension = DIM;
00124
00128 static bool local_operator() { return true; }
00129
00133 double operator_order() const { return 1; }
00134
00138 double D(const typename AggregatedFrame<IBASIS,DIM>::Index& lambda) const;
00139
00143 void rescale(InfiniteVector<double, typename AggregatedFrame<IBASIS,DIM>::Index>& coeffs,
00144 const int n) const;
00145
00149 double a(const typename AggregatedFrame<IBASIS,DIM>::Index& lambda,
00150 const typename AggregatedFrame<IBASIS,DIM>::Index& nu) const;
00151
00155 double norm_A() const;
00156
00163 double norm_Ainv() const { return normAinv; };
00164
00168 void set_norm_A(const double _normA) { normA = _normA; }
00169
00173 void set_Ainv(const double nAinv) { normAinv = nAinv; };
00174
00178 double s_star() const;
00179
00184 double alphak(const unsigned int k) const {
00185
00186 return 2.*norm_A();
00187 }
00188
00192 double f(const typename AggregatedFrame<IBASIS,DIM>::Index& lambda) const;
00193
00198 void RHS(const double eta, InfiniteVector<double,
00199 typename AggregatedFrame<IBASIS,DIM>::Index>& coeffs) const;
00200
00205 void RHS(const double eta, const int p,
00206 InfiniteVector<double,
00207 typename AggregatedFrame<IBASIS,DIM>::Index>& coeffs) const;
00208
00212 double F_norm() const { return sqrt(fnorm_sqr); }
00213
00217 double F_norm_local(const int patch) const { return sqrt(fnorms_sqr_patch[patch]); }
00218
00222 void set_bvp(const EllipticBVP<DIM>*);
00223
00224
00229 void add_level (const Index& lambda,
00230 InfiniteVector<double, Index>& w, const int j,
00231 const double factor,
00232 const int J,
00233 const CompressionStrategy strategy) const;
00234
00235 protected:
00236
00240 const EllipticBVP<DIM>* ell_bvp_;
00241
00245 const AggregatedFrame<IBASIS,DIM>* frame_;
00246
00247
00248
00249
00250 typedef std::map<Index1D<IBASIS>,double > Column1D;
00251 typedef std::map<Index1D<IBASIS>,Column1D> One_D_IntegralCache;
00252
00253 mutable One_D_IntegralCache one_d_integrals;
00254
00255
00256 private:
00257
00275 double integrate(const Index1D<IBASIS>& lambda,
00276 const Index1D<IBASIS>& mu,
00277 const int N_Gauss,
00278 const int dir,
00279 const typename CubeBasis<IBASIS,DIM>::Support* supp_lambda,
00280 const typename CubeBasis<IBASIS,DIM>::Support* supp_mu) const;
00281
00283 void compute_rhs();
00284
00286 void compute_diagonal();
00287
00289 const int jmax_;
00290
00292 Array1D<std::pair<typename AggregatedFrame<IBASIS,DIM>::Index,double> > fcoeffs;
00293
00295 Array1D<Array1D<std::pair<typename AggregatedFrame<IBASIS,DIM>::Index,double> > > fcoeffs_patch;
00296
00298 Array1D<double> stiff_diagonal;
00299
00301 double fnorm_sqr;
00302
00307 Array1D<double> fnorms_sqr_patch;
00308
00309
00310
00311
00312
00313 mutable double normA, normAinv;
00314 };
00315 }
00316
00317 #include <simple_elliptic_equation.cpp>
00318
00319 #endif