NGSolve  4.9
comp/gridfunction.hpp
00001 #ifndef FILE_GRIDFUNCTION
00002 #define FILE_GRIDFUNCTION
00003 
00004 /*********************************************************************/
00005 /* File:   gridfunction.hpp                                           */
00006 /* Author: Joachim Schoeberl                                         */
00007 /* Date:   10. Jul. 2000                                             */
00008 /*********************************************************************/
00009 
00010 namespace ngcomp
00011 {
00012 
00016   class NGS_DLL_HEADER GridFunction : public NGS_Object
00017   {
00018   protected:
00020     const FESpace & fespace;
00022     bool nested;
00024     bool visual;
00026     int multidim;
00028     int level_updated;
00029 
00031     int cacheblocksize;
00033     Array<BaseVector*> vec;
00035     Array<GridFunction*> compgfs;
00036 
00037   public:
00039     GridFunction (const FESpace & afespace, 
00040                   const string & name = "gfu", 
00041                   const Flags & flags = Flags());
00043     virtual ~GridFunction ();
00045     virtual void Update ();
00047     virtual BaseVector & GetVector (int comp = 0) { return *(vec[comp]); }
00049     virtual const BaseVector & GetVector (int comp = 0) const  { return *(vec[comp]); }
00050 
00051 
00052     /*
00053     operator BaseVector& () { return GetVector(); }
00054     template <typename T> 
00055     GridFunction & operator= (const VVecExpr<T> & v) { GetVector() = v; return *this; }
00056     GridFunction & operator= (const BaseVector & v) { GetVector() = v; return *this; }
00057     template <typename T> 
00058     GridFunction & operator+= (const VVecExpr<T> & v) { GetVector() += v; return *this; }
00059     GridFunction & operator+= (const BaseVector & v) { GetVector() += v; return *this; }
00060     template <typename T> 
00061     GridFunction & operator-= (const VVecExpr<T> & v) { GetVector() -= v; return *this; }
00062     GridFunction & operator-= (const BaseVector & v) { GetVector() -= v; return *this; }
00063     */
00064 
00066     void SetNested (int anested = 1) { nested = anested; }
00068     void SetVisual (bool avisual = 1) { visual = avisual; }
00069 
00070     int GetMultiDim () const { return multidim; }
00071   
00072     int GetLevelUpdated() const { return level_updated; }
00074     const FESpace & GetFESpace() const
00075     { return fespace; }
00076 
00078     virtual string GetClassName () const
00079     {
00080       return "GridFunction";
00081     }
00082 
00084     virtual void PrintReport (ostream & ost);
00086     virtual void MemoryUsage (Array<MemoryUsageStruct*> & mu) const;
00087 
00089     void Visualize(const string & name);
00090 
00092     virtual void SetCacheBlockSize (const int size) {cacheblocksize = size;}
00094     virtual int GetCacheBlockSize (void) const { return cacheblocksize;}
00096     virtual bool IsUpdated () const; 
00097 
00098     
00099     int GetNComponents () const { return compgfs.Size(); }
00100     GridFunction * GetComponent (int compound_comp) const;
00101 
00102 
00104     virtual void GetElementVector (const FlatArray<int> & dnums,
00105                                    FlatVector<double> & elvec) const
00106     { vec[0] -> GetIndirect (dnums, elvec); }
00107       
00109     virtual void SetElementVector (const FlatArray<int> & dnums,
00110                                    const FlatVector<double> & elvec) 
00111     { vec[0] -> SetIndirect (dnums, elvec); }
00112 
00113 
00115     virtual void GetElementVector (int comp,
00116                                    const FlatArray<int> & dnums,
00117                                    FlatVector<double> & elvec) const 
00118     { vec[comp] -> GetIndirect (dnums, elvec); }
00119 
00120 
00122     virtual void SetElementVector (int comp,
00123                                    const FlatArray<int> & dnums,
00124                                    const FlatVector<double> & elvec) 
00125     { vec[comp] -> SetIndirect (dnums, elvec); }
00126 
00128     virtual void GetElementVector (const FlatArray<int> & dnums,
00129                                    FlatVector<Complex> & elvec) const
00130     { vec[0] -> GetIndirect (dnums, elvec); }
00131       
00133     virtual void SetElementVector (const FlatArray<int> & dnums,
00134                                    const FlatVector<Complex> & elvec) 
00135     { vec[0] -> SetIndirect (dnums, elvec); }
00136 
00137 
00139     virtual void GetElementVector (int comp,
00140                                    const FlatArray<int> & dnums,
00141                                    FlatVector<Complex> & elvec) const 
00142     { vec[comp] -> GetIndirect (dnums, elvec); }
00143 
00144 
00146     virtual void SetElementVector (int comp,
00147                                    const FlatArray<int> & dnums,
00148                                    const FlatVector<Complex> & elvec) 
00149     { vec[comp] -> SetIndirect (dnums, elvec); }
00150 
00151 
00152 
00153     virtual void Load (istream & ist) = 0;
00154     virtual void Save (ostream & ost) const = 0;
00155   };
00156 
00157 
00158 
00159   template <class SCAL>
00160   class NGS_DLL_HEADER S_GridFunction : public GridFunction
00161   {
00162   public:
00163     S_GridFunction (const FESpace & afespace, const string & aname, const Flags & flags)
00164       : GridFunction (afespace, aname, flags) { ; }
00165   
00166 
00167     // parallel Load/Save by Martin Huber and Lothar Nannen 
00168     virtual void Load (istream & ist);
00169     virtual void Save (ostream & ost) const;
00170 
00171   private:
00172     template <int N, NODE_TYPE NT> void LoadNodeType (istream & ist);
00173 
00174     template <int N, NODE_TYPE NT> void SaveNodeType (ostream & ost) const;
00175   };
00176 
00177 
00178 
00179   template <class TV>
00180   class NGS_DLL_HEADER T_GridFunction : public S_GridFunction<typename mat_traits<TV>::TSCAL>
00181   {
00182     using S_GridFunction<typename mat_traits<TV>::TSCAL>::vec;
00183     using S_GridFunction<typename mat_traits<TV>::TSCAL>::compgfs;
00184 
00185   public:
00186     typedef typename mat_traits<TV>::TSCAL TSCAL;
00187     enum { VDIM = mat_traits<TV>::HEIGHT };
00188 
00189     T_GridFunction (const FESpace & afespace, 
00190                     const string & aname = "gfu", 
00191                     const Flags & flags = Flags());
00192     virtual ~T_GridFunction ();
00193 
00194     virtual void Update ();
00195   };
00196 
00197 
00198   extern NGS_DLL_HEADER GridFunction * CreateGridFunction (const FESpace * space,
00199                                             const string & name, const Flags & flags);
00200 
00201 
00202 
00203 
00204   template <class SCAL>
00205   class NGS_DLL_HEADER S_ComponentGridFunction : public S_GridFunction<SCAL>
00206   {
00207     const S_GridFunction<SCAL> & gf_parent;
00208     int comp;
00209   public:
00210     S_ComponentGridFunction (const S_GridFunction<SCAL> & agf_parent, int acomp);
00211     virtual ~S_ComponentGridFunction ();
00212     virtual void Update ();
00213   };
00214   
00215 
00216 
00217 
00218   class NGS_DLL_HEADER GridFunctionCoefficientFunction : public CoefficientFunction
00219   {
00220   protected:
00221     GridFunction & gf;
00222     const DifferentialOperator * diffop;
00223     int comp;
00224   public:
00225     GridFunctionCoefficientFunction (GridFunction & agf, int acomp = 0);
00226     GridFunctionCoefficientFunction (GridFunction & agf, DifferentialOperator * adiffop, int acomp = 0);
00227     
00228     virtual ~GridFunctionCoefficientFunction ();
00230     virtual bool IsComplex() const;
00231     virtual int Dimension() const;
00232     virtual double Evaluate (const BaseMappedIntegrationPoint & ip) const;
00233 
00234     virtual void Evaluate(const BaseMappedIntegrationPoint & ip,
00235                           FlatVector<> result) const;
00236     virtual void Evaluate(const BaseMappedIntegrationPoint & ip,
00237                           FlatVector<Complex> result) const;
00238     
00239     virtual void Evaluate (const BaseMappedIntegrationRule & ir, 
00240                            FlatMatrix<double> values) const;
00241   };
00242 
00243 
00244 
00245   template <class SCAL>
00246   class NGS_DLL_HEADER VisualizeGridFunction : public netgen::SolutionData
00247   {
00248     const MeshAccess & ma;
00249     const S_GridFunction<SCAL> * gf;
00250     Array<const BilinearFormIntegrator *> bfi2d;
00251     Array<const BilinearFormIntegrator *> bfi3d;
00252     bool applyd;
00253     //
00254     int cache_elnr;
00255     bool cache_bound;
00256     LocalHeap lh;
00257     // ElementTransformation eltrans;
00258     const FiniteElement * fel;
00259     Array<int> dnums;
00260     FlatVector<SCAL> elu;
00261 
00262   public:
00263     VisualizeGridFunction (const MeshAccess & ama,
00264                            const GridFunction * agf,
00265                            const BilinearFormIntegrator * abfi2d,
00266                            const BilinearFormIntegrator * abfi3d,
00267                            bool aapplyd);
00268     VisualizeGridFunction (const MeshAccess & ama,
00269                            const GridFunction * agf,
00270                            const Array<BilinearFormIntegrator *> & abfi2d,
00271                            const Array<BilinearFormIntegrator *> & abfi3d,
00272                            bool aapplyd);
00273 
00274     virtual ~VisualizeGridFunction ();
00275   
00276     virtual bool GetValue (int elnr, 
00277                            double lam1, double lam2, double lam3,
00278                            double * values) ;
00279 
00280     virtual bool GetValue (int elnr, 
00281                            const double xref[], const double x[], const double dxdxref[],
00282                            double * values) ;
00283 
00284     virtual bool GetMultiValue (int elnr, int facetnr, int npts,
00285                                 const double * xref, int sxref,
00286                                 const double * x, int sx,
00287                                 const double * dxdxref, int sdxdxref,
00288                                 double * values, int svalues);
00289 
00290     virtual bool GetSurfValue (int elnr, int facetnr,
00291                                double lam1, double lam2, 
00292                                double * values) ;
00293 
00294     virtual bool GetSurfValue (int selnr, int facetnr, 
00295                                const double xref[], const double x[], const double dxdxref[],
00296                                double * values);
00297 
00298     virtual bool GetMultiSurfValue (int selnr, int facetnr, int npts,
00299                                     const double * xref, int sxref,
00300                                     const double * x, int sx,
00301                                     const double * dxdxref, int sdxdxref,
00302                                     double * values, int svalues);
00303 
00304 
00305     virtual int GetNumMultiDimComponents ()
00306     {
00307       return gf->GetMultiDim();
00308     }
00309 
00310 
00311     void Analyze(Array<double> & minima, Array<double> & maxima, Array<double> & averages, int component = -1);
00312     void Analyze(Array<double> & minima, Array<double> & maxima, Array<double> & averages_times_volumes, Array<double> & volumes, int component = -1);
00313     
00314 
00315   };
00316 
00317 
00318 
00319 
00320 
00321 
00322 
00323 
00324 
00325 
00326 
00327 
00328 
00329   class NGS_DLL_HEADER VisualizeCoefficientFunction : public netgen::SolutionData
00330   {
00331     const MeshAccess & ma;
00332     const CoefficientFunction * cf;
00333     LocalHeap lh;
00334   public:
00335     VisualizeCoefficientFunction (const MeshAccess & ama,
00336                                   const CoefficientFunction * acf);
00337     
00338     virtual ~VisualizeCoefficientFunction ();
00339   
00340     virtual bool GetValue (int elnr, 
00341                            double lam1, double lam2, double lam3,
00342                            double * values) ;
00343 
00344     virtual bool GetValue (int elnr, 
00345                            const double xref[], const double x[], const double dxdxref[],
00346                            double * values) ;
00347 
00348     virtual bool GetMultiValue (int elnr, int npts,
00349                                 const double * xref, int sxref,
00350                                 const double * x, int sx,
00351                                 const double * dxdxref, int sdxdxref,
00352                                 double * values, int svalues);
00353 
00354     virtual bool GetSurfValue (int elnr, int facetnr,
00355                                double lam1, double lam2, 
00356                                double * values) ;
00357 
00358     virtual bool GetSurfValue (int selnr, int facetnr, 
00359                                const double xref[], const double x[], const double dxdxref[],
00360                                double * values);
00361 
00362     virtual bool GetMultiSurfValue (int selnr, int facetnr, int npts,
00363                                     const double * xref, int sxref,
00364                                     const double * x, int sx,
00365                                     const double * dxdxref, int sdxdxref,
00366                                     double * values, int svalues);
00367 
00368 
00369     void Analyze(Array<double> & minima, Array<double> & maxima, Array<double> & averages, int component = -1);
00370     void Analyze(Array<double> & minima, Array<double> & maxima, Array<double> & averages_times_volumes, Array<double> & volumes, int component = -1);
00371   };
00372 
00373 
00374 
00375 }
00376 
00377 
00378 
00379 
00380 #endif