vtkwriter.hh

Go to the documentation of this file.
00001 // $Id: vtkwriter.hh 3963 2008-01-28 12:52:56Z christi $
00002 
00003 #ifndef DUNE_VTKWRITER_HH
00004 #define DUNE_VTKWRITER_HH
00005 
00006 #include <iostream>
00007 #include <fstream>
00008 #include <vector>
00009 #include <list>
00010 #include <string.h>
00011 #include <dune/common/exceptions.hh>
00012 #include <dune/common/iteratorfacades.hh>
00013 #include <dune/grid/common/mcmgmapper.hh>
00014 #include <dune/grid/common/referenceelements.hh>
00015 
00016 
00017 // namespace base64
00018 // {
00019 // #include"cencode.c"
00020 // }
00021 
00032 namespace Dune
00033 {
00036   struct VTKOptions
00037   {
00038     enum OutputType {
00040       ascii,
00042       binary, 
00044       binaryappended
00045     };
00046     enum DataMode {
00047       conforming, nonconforming
00048     };
00049   };
00050 
00051 
00052   // map type to name in data array
00053   template<class T>
00054   struct VTKTypeNameTraits {
00055     std::string operator () (){
00056       return "";
00057     }
00058   };
00059 
00060   template<>
00061   struct VTKTypeNameTraits<char> {
00062     std::string operator () () {
00063       return "Int8";
00064     }
00065     typedef int PrintType;
00066   };
00067 
00068   template<>
00069   struct VTKTypeNameTraits<unsigned char> {
00070     std::string operator () () {
00071       return "UInt8";
00072     }   
00073     typedef int PrintType;
00074   };
00075 
00076   template<>
00077   struct VTKTypeNameTraits<short> {
00078     std::string operator () () {
00079       return "Int16";
00080     }   
00081     typedef short PrintType;
00082   };
00083 
00084   template<>
00085   struct VTKTypeNameTraits<unsigned short> {
00086     std::string operator () () {
00087       return "UInt16";
00088     }   
00089     typedef unsigned short PrintType;
00090   };
00091 
00092   template<>
00093   struct VTKTypeNameTraits<int> {
00094     std::string operator () () {
00095       return "Int32";
00096     }   
00097     typedef int PrintType;
00098   };
00099 
00100   template<>
00101   struct VTKTypeNameTraits<unsigned int> {
00102     std::string operator () () {
00103       return "UInt32";
00104     }   
00105     typedef unsigned int PrintType;
00106   };
00107 
00108   template<>
00109   struct VTKTypeNameTraits<float> {
00110     std::string operator () () {
00111       return "Float32";
00112     }   
00113     typedef float PrintType;
00114   };
00115 
00116   template<>
00117   struct VTKTypeNameTraits<double> {
00118     std::string operator () () {
00119       return "Float64";
00120     }   
00121     typedef double PrintType;
00122   };
00123 
00124 
00133   template<class GridImp, class IS=typename GridImp::template Codim<0>::LeafIndexSet>
00134   class VTKWriter {
00135     template<int dim>
00136     struct P0Layout
00137     {
00138       bool contains (Dune::GeometryType gt)
00139         {
00140             return gt.dim()==dim;
00141         }
00142     };
00143     
00144     template<int dim>
00145     struct P1Layout
00146     {
00147       bool contains (Dune::GeometryType gt)
00148         {
00149             return gt.dim()==0;
00150         }
00151     };    
00152 
00153     // extract types
00154     enum {n=GridImp::dimension};
00155     enum {w=GridImp::dimensionworld};
00156     typedef typename GridImp::ctype DT;
00157     typedef typename GridImp::Traits::template Codim<0>::Entity Entity;
00158     typedef typename GridImp::Traits::template Codim<0>::Entity Cell;
00159     typedef typename GridImp::Traits::template Codim<n>::Entity Vertex;
00160     typedef IS IndexSet;
00161     static const PartitionIteratorType vtkPartition = InteriorBorder_Partition;
00162     typedef typename IS::template Codim<0>::template Partition<vtkPartition>::Iterator GridCellIterator;
00163     typedef typename IS::template Codim<n>::template Partition<vtkPartition>::Iterator GridVertexIterator;
00164     typedef MultipleCodimMultipleGeomTypeMapper<GridImp,IS,P1Layout> VertexMapper;
00165   public:
00166 
00172     class VTKFunction
00173     {
00174     public:
00176       virtual int ncomps () const = 0;
00177 
00179 
00185       virtual double evaluate (int comp, const Entity& e, const Dune::FieldVector<DT,n>& xi) const = 0;
00186 
00188       virtual std::string name () const = 0;
00189 
00191       virtual ~VTKFunction () {}
00192     };
00193 
00194   private:
00195     typedef typename std::list<VTKFunction*>::iterator FunctionIterator;
00196     
00197     class CellIterator :
00198       public ForwardIteratorFacade<CellIterator, Entity, Entity&, int>
00199     {
00200       GridCellIterator git;
00201       GridCellIterator gend;
00202     public:
00203       CellIterator(const GridCellIterator & x, const GridCellIterator & end) : git(x), gend(end) {};
00204       void increment ()
00205         {
00206           ++git;
00207           while (git!=gend && git->partitionType()!=InteriorEntity) ++git;
00208         }
00209       bool equals (const CellIterator & cit) const
00210         {
00211           return git == cit.git;
00212         }
00213       Entity& dereference() const
00214         {
00215           return *git;
00216         }
00217       const FieldVector<DT,n> position() const
00218         {
00219           return ReferenceElements<DT,n>::general(git->type()).position(0,0);
00220         }
00221     };
00222     CellIterator cellBegin() const
00223       {
00224         return CellIterator(is.template begin<0,vtkPartition>(), is.template end<0,vtkPartition>());
00225       }
00226     CellIterator cellEnd() const
00227       {
00228         return CellIterator(is.template end<0,vtkPartition>(), is.template end<0,vtkPartition>());
00229       }
00230     
00231     class VertexIterator :
00232       public ForwardIteratorFacade<VertexIterator, Entity, Entity&, int>
00233     {
00234       GridCellIterator git;
00235       GridCellIterator gend;
00236       VTKOptions::DataMode datamode;
00237       int index;
00238       const VertexMapper & vertexmapper;
00239       std::vector<bool> visited;
00240       const std::vector<int> & number;
00241       int offset;
00242     protected:
00243       void basicIncrement ()
00244         {
00245           if (git == gend) return;
00246           index++;
00247           if (index == git->template count<n>()) {
00248             offset += git->template count<n>();
00249             index = 0;
00250             ++git;
00251             if(git == gend) return;
00252             while (git->partitionType()!=InteriorEntity) ++git;
00253           }
00254         }
00255     public:
00256       VertexIterator(const GridCellIterator & x,
00257                      const GridCellIterator & end,
00258                      const VTKOptions::DataMode & dm,
00259                      const VertexMapper & vm,
00260                      const std::vector<int> & num) :
00261         git(x), gend(end), datamode(dm), index(0),
00262         vertexmapper(vm), visited(vm.size(), false),
00263         number(num), offset(0)
00264         {
00265           if (datamode == VTKOptions::conforming && git != gend)
00266             visited[vertexmapper.template map<n>(*git,index)] = true;
00267         };
00268       void increment ()
00269         {
00270           switch (datamode)
00271           {
00272           case VTKOptions::conforming:
00273             while(visited[vertexmapper.template map<n>(*git,index)])
00274             {
00275               basicIncrement();
00276               if (git == gend) return;
00277             }
00278             visited[vertexmapper.template map<n>(*git,index)] = true;
00279             break;
00280           case VTKOptions::nonconforming:
00281             basicIncrement();
00282             break;
00283           }
00284        }
00285       bool equals (const VertexIterator & cit) const
00286         {
00287           return git == cit.git
00288             && index == cit.index && datamode == cit.datamode;
00289         }
00290       Entity& dereference() const
00291         {
00292           return *git;
00293         }
00294       int id () const
00295         {
00296           switch (datamode)
00297           {
00298           case VTKOptions::conforming:
00299             return
00300               number[vertexmapper.template map<n>(*git,renumber(*git,index))];
00301           case VTKOptions::nonconforming:
00302             return offset + renumber(*git,index);
00303           default:
00304             DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
00305           }
00306         }
00307       int localindex () const
00308         {
00309           return index;
00310         }
00311       const FieldVector<DT,n> & position () const
00312         {
00313           return ReferenceElements<DT,n>::general(git->type()).position(index,n);
00314         }
00315     };    
00316     VertexIterator vertexBegin() const
00317       {
00318         return VertexIterator(is.template begin<0,vtkPartition>(),
00319                               is.template end<0,vtkPartition>(),
00320                               datamode, *vertexmapper, number);
00321       }
00322     VertexIterator vertexEnd() const
00323       {
00324         return VertexIterator(is.template end<0,vtkPartition>(),
00325                               is.template end<0,vtkPartition>(),
00326                               datamode, *vertexmapper, number);
00327       }    
00328     
00329     class CornerIterator :
00330       public ForwardIteratorFacade<CornerIterator, Entity, Entity&, int>
00331     {
00332       GridCellIterator git;
00333       GridCellIterator gend;
00334       VTKOptions::DataMode datamode;
00335       int index;
00336       const VertexMapper & vertexmapper;
00337       std::vector<bool> visited;
00338       const std::vector<int> & number;
00339       int offset;
00340     protected:
00341       void basicIncrement ()
00342         {
00343           if (git == gend) return;
00344           index++;
00345           if (index == git->template count<n>()) {
00346             offset += git->template count<n>();
00347             index = 0;
00348             ++git;
00349             if (git == gend) return;
00350             while (git->partitionType()!=InteriorEntity) ++git;
00351           }
00352         }
00353     public:
00354       CornerIterator(const GridCellIterator & x,
00355                      const GridCellIterator & end,
00356                      const VTKOptions::DataMode & dm,
00357                      const VertexMapper & vm,
00358                      const std::vector<int> & num) :
00359         git(x), gend(end), datamode(dm), index(0),
00360         vertexmapper(vm),
00361         number(num), offset(0) {};
00362       void increment ()
00363         {
00364           basicIncrement();
00365         }
00366       bool equals (const CornerIterator & cit) const
00367         {
00368           return git == cit.git
00369             && index == cit.index && datamode == cit.datamode;
00370         }
00371       Entity& dereference() const
00372         {
00373           return *git;
00374         }
00375       int id () const
00376         {
00377           switch (datamode)
00378           {
00379           case VTKOptions::conforming:
00380             return
00381               number[vertexmapper.template map<n>(*git,renumber(*git,index))];
00382           case VTKOptions::nonconforming:
00383             return offset + renumber(*git,index);
00384           default:
00385             DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
00386           }
00387         }
00388       int localindex () const
00389         {
00390           return index;
00391         }
00392     };    
00393     CornerIterator cornerBegin() const
00394       {
00395         return CornerIterator(is.template begin<0,vtkPartition>(),
00396                               is.template end<0,vtkPartition>(),
00397                               datamode, *vertexmapper, number);
00398       }
00399     CornerIterator cornerEnd() const
00400       {
00401         return CornerIterator(is.template end<0,vtkPartition>(),
00402                               is.template end<0,vtkPartition>(),
00403                               datamode, *vertexmapper, number);
00404       }    
00405         
00409     template<class V>
00410     class P0VectorWrapper : public VTKFunction  
00411     {
00412       typedef MultipleCodimMultipleGeomTypeMapper<GridImp,IS,P0Layout> VM0;
00413     public:
00415       virtual int ncomps () const
00416         {
00417           return 1;
00418         }
00419 
00421       virtual double evaluate (int comp, const Entity& e, const Dune::FieldVector<DT,n>& xi) const
00422         {
00423           return v[mapper.map(e)];
00424         }
00425       
00427       virtual std::string name () const
00428         {
00429           return s;
00430         }
00431 
00433       P0VectorWrapper (const GridImp& g_, const IS& is_, const V& v_, std::string s_) 
00434         : g(g_), is(is_), v(v_), s(s_), mapper(g_,is_)
00435         {
00436           if (v.size()!=(unsigned int)mapper.size())
00437             DUNE_THROW(IOError,"VTKWriter::P0VectorWrapper: size mismatch");
00438         }
00439 
00440       virtual ~P0VectorWrapper() {}
00441       
00442     private:
00443       const GridImp& g;
00444       const IS& is;
00445       const V& v;
00446       std::string s;
00447       VM0 mapper;
00448     };
00449 
00453     template<class V>
00454     class P1VectorWrapper : public VTKFunction  
00455     {
00456       typedef MultipleCodimMultipleGeomTypeMapper<GridImp,IS,P1Layout> VM1;
00457     public:
00459       virtual int ncomps () const
00460         {
00461           return 1;
00462         }
00463 
00465       virtual double evaluate (int comp, const Entity& e, const Dune::FieldVector<DT,n>& xi) const
00466         {
00467           double min=1E100;
00468           int imin=-1;
00469           Dune::GeometryType gt = e.type();
00470           for (int i=0; i<e.template count<n>(); ++i)
00471           {
00472             Dune::FieldVector<DT,n> 
00473               local = Dune::ReferenceElements<DT,n>::general(gt).position(i,n);
00474             local -= xi;
00475             if (local.infinity_norm()<min)
00476             {
00477               min = local.infinity_norm();
00478               imin = i;
00479             }
00480           }
00481           return v[mapper.template map<n>(e,imin)];
00482         }
00483       
00485       virtual std::string name () const
00486         {
00487           return s;
00488         }
00489 
00491       P1VectorWrapper (const GridImp& g_, const IS& is_, const V& v_, std::string s_) 
00492         : g(g_), is(is_), v(v_), s(s_), mapper(g_,is_)
00493         {
00494           if (v.size()!=(unsigned int)mapper.size())
00495             DUNE_THROW(IOError,"VTKWriter::P1VectorWrapper: size mismatch");
00496         }
00497 
00498       virtual ~P1VectorWrapper() {}
00499       
00500     private:
00501       const GridImp& g;
00502       const IS& is;
00503       const V& v;
00504       std::string s;
00505       VM1 mapper;
00506     };
00507 
00508   public:
00518     VTKWriter (const GridImp& g, VTKOptions::DataMode dm = VTKOptions::conforming) :
00519       grid(g), is(grid.leafIndexSet()), datamode(dm)
00520       {
00521         indentCount = 0;
00522         numPerLine = 4*3; //should be a multiple of 3 !
00523       }
00524  
00533     VTKWriter (const GridImp& g, const IndexSet& i, VTKOptions::DataMode dm = VTKOptions::conforming) :
00534       grid(g), is(i), datamode(dm)
00535       {
00536         indentCount = 0;
00537         numPerLine = 4*3; //should be a multiple of 3 !
00538       }
00539 
00544     void addCellData (VTKFunction* p)
00545       {
00546         celldata.push_back(p);
00547       }
00548 
00560     template<class V>
00561     void addCellData (const V& v, std::string name)
00562       {
00563         VTKFunction* p = new P0VectorWrapper<V>(grid,is,v,name);
00564         celldata.push_back(p);
00565       }
00566 
00571     void addVertexData (VTKFunction* p)
00572       {
00573         vertexdata.push_back(p);
00574       }
00575 
00587     template<class V>
00588     void addVertexData (const V& v, std::string name)
00589       {
00590         VTKFunction* p = new P1VectorWrapper<V>(grid,is,v,name);
00591         vertexdata.push_back(p);
00592       }
00593 
00595     void clear ()
00596       {
00597         for (FunctionIterator it=celldata.begin();
00598              it!=celldata.end(); ++it)
00599           delete *it;
00600         celldata.clear();
00601         for (FunctionIterator it=vertexdata.begin();
00602              it!=vertexdata.end(); ++it)
00603           delete *it;
00604         vertexdata.clear();
00605       }
00606 
00608     ~VTKWriter ()
00609       {
00610         this->clear();
00611       }
00612 
00618     void write (const char* name, VTKOptions::OutputType ot = VTKOptions::ascii)
00619       {
00620         // make data mode visible to private functions
00621         outputtype=ot;
00622 
00623         // reset byte counter for binary appended output
00624         bytecount = 0;
00625 
00626         if (grid.comm().size()==1)
00627         {
00628           std::ofstream file;
00629           char fullname[128];
00630           if (n>1)
00631             sprintf(fullname,"%s.vtu",name);
00632           else
00633             sprintf(fullname,"%s.vtp",name);
00634           if (outputtype==VTKOptions::binaryappended)
00635             file.open(fullname,std::ios::binary);
00636           else
00637             file.open(fullname);
00638           writeDataFile(file);
00639           file.close();
00640         }
00641         else
00642         {
00643           std::ofstream file;
00644           char fullname[128];
00645           if (n>1)
00646             sprintf(fullname,"s%04d:p%04d:%s.vtu",grid.comm().size(),grid.comm().rank(),name);
00647           else
00648             sprintf(fullname,"s%04d:p%04d:%s.vtp",grid.comm().size(),grid.comm().rank(),name);
00649           if (outputtype==VTKOptions::binaryappended)
00650             file.open(fullname,std::ios::binary);
00651           else
00652             file.open(fullname);
00653           writeDataFile(file);
00654           file.close();
00655           grid.comm().barrier();
00656           if (grid.comm().rank()==0)
00657           {
00658             if (n>1)
00659               sprintf(fullname,"s%04d:%s.pvtu",grid.comm().size(),name);
00660             else
00661               sprintf(fullname,"s%04d:%s.pvtp",grid.comm().size(),name);
00662             file.open(fullname);
00663             writeParallelHeader(file,name,"");
00664             file.close();
00665           }
00666           grid.comm().barrier();
00667         }
00668       }
00669 
00671     void pwrite (const char* name,  const char* path, const char* extendpath, 
00672                  VTKOptions::OutputType ot = VTKOptions::ascii)
00673       {
00674         // make data mode visible to private functions
00675         outputtype=ot;
00676 
00677         // reset byte counter for binary appended output
00678         bytecount = 0;
00679 
00680         // do some magic because paraview can only cope with relative pathes to piece files
00681         std::ofstream file;
00682         char piecepath[256];
00683         char relpiecepath[256];
00684         int n=strlen(path);
00685         int m=strlen(extendpath);
00686         if (n>0 && path[0]=='/' && path[n-1]=='/')
00687         {
00688           // 1) path is an absolute path to the directory where the pvtu file will be placed
00689           // 2) extendpath is an absolute path from "/" where the pieces are placed
00690           // 3) pieces are addressed relative in the pvtu files
00691           if (m==0)
00692           {
00693             // write pieces to root :-)
00694             piecepath[0] = '/';
00695             piecepath[1] = '\0';
00696           }
00697           else
00698           {
00699             // make piecepath absolute with trailing "/"
00700             char *p=piecepath;
00701             if (extendpath[0]!='/')
00702             {
00703               *p = '/';
00704               p++;
00705             }
00706             for (int i=0; i<m; i++)
00707             {
00708               *p = extendpath[i];
00709               p++;
00710             }
00711             if (*(p-1)!='/')
00712             {
00713               *p = '/';
00714               p++;
00715             }
00716             *p = '\0';
00717           }
00718           // path and piecepath are either "/" or have leading and trailing /
00719           // count slashes in path 
00720           int k=0;
00721           const char *p=path;
00722           while (*p!='\0')
00723           {
00724             if (*p=='/') k++;
00725             p++;
00726           }
00727           char *pp = relpiecepath;
00728           if (k>1)
00729           {
00730             for (int i=0; i<k; i++)
00731             {
00732               *pp='.'; pp++; *pp='.'; pp++; *pp='/'; pp++; 
00733             }
00734           }
00735           // now copy the extendpath
00736           for (int i=0; i<m; i++)
00737           {
00738             if (i==0 && extendpath[i]=='/') continue;
00739             *pp = extendpath[i];
00740             pp++;
00741           }
00742           if ( pp!=relpiecepath && (*(pp-1)!='/') )
00743           {
00744             *pp = '/';
00745             pp++;
00746           }
00747           *pp = '\0';
00748         }
00749         else
00750         {
00751           // 1) path is a relative path to the directory where pvtu files are placed
00752           // 2) extendpath is relative to where the pvtu files are and there the pieces are placed
00753           if (n==0 || m==0)
00754             sprintf(piecepath,"%s%s",path,extendpath);
00755           else
00756           {
00757             // both are non-zero
00758             if (path[n-1]!='/' && extendpath[0]!='/')
00759               sprintf(piecepath,"%s/%s",path,extendpath);
00760             else
00761               sprintf(piecepath,"%s%s",path,extendpath);
00762           }
00763           // the pieces are relative to the pvtu files
00764           sprintf(relpiecepath,"%s",extendpath);
00765         }
00766         char fullname[256];
00767         if (n>1)
00768           sprintf(fullname,"%s/s%04d:p%04d:%s.vtu",piecepath,grid.comm().size(),grid.comm().rank(),name);
00769         else
00770           sprintf(fullname,"%s/s%04d:p%04d:%s.vtp",piecepath,grid.comm().size(),grid.comm().rank(),name);
00771         if (outputtype==VTKOptions::binaryappended)
00772           file.open(fullname,std::ios::binary);
00773         else
00774           file.open(fullname);
00775         writeDataFile(file);
00776         file.close();
00777         grid.comm().barrier();
00778         if (grid.comm().rank()==0)
00779         {
00780           if (n>1)
00781             sprintf(fullname,"%s/s%04d:%s.pvtu",path,grid.comm().size(),name);
00782           else
00783             sprintf(fullname,"%s/s%04d:%s.pvtp",path,grid.comm().size(),name);
00784           file.open(fullname);
00785           writeParallelHeader(file,name,relpiecepath);
00786           file.close();
00787         }
00788         grid.comm().barrier();
00789       }
00790 
00791   private:
00792 
00793     enum VTKGeometryType
00794     {
00795       vtkLine = 3,
00796       vtkTriangle = 5,
00797       vtkQuadrilateral = 9,
00798       vtkTetrahedron = 10,
00799       vtkHexahedron = 12,
00800       vtkPrism = 13,
00801       vtkPyramid = 14
00802     };
00803     
00805     static VTKGeometryType vtkType(const Dune::GeometryType & t)
00806       {
00807         if (t.isLine())
00808           return vtkLine;
00809         if (t.isTriangle())
00810           return vtkTriangle;
00811         if (t.isQuadrilateral())
00812           return vtkQuadrilateral;
00813         if (t.isTetrahedron())
00814           return vtkTetrahedron;
00815         if (t.isPyramid())
00816           return vtkPyramid;
00817         if (t.isPrism())
00818           return vtkPrism;
00819         if (t.isHexahedron())
00820           return vtkHexahedron;
00821         DUNE_THROW(IOError,"VTKWriter: unsupported GeometryType " << t);
00822       }
00823 
00825     void writeParallelHeader (std::ostream& s, const char* piecename, const char* piecepath)
00826       {
00827         // xml header
00828         s << "<?xml version=\"1.0\"?>" << std::endl;
00829 
00830         // VTKFile
00831         if (n>1)
00832           s << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
00833         else
00834           s << "<VTKFile type=\"PPolyData\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
00835         indentUp();
00836 
00837         // PUnstructuredGrid
00838         indent(s); 
00839         if (n>1)
00840           s << "<PUnstructuredGrid GhostLevel=\"0\">" << std::endl;
00841         else
00842           s << "<PPolyData GhostLevel=\"0\">" << std::endl;
00843         indentUp();
00844 
00845         // PPointData
00846         indent(s); s << "<PPointData";
00847         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
00848           if ((*it)->ncomps()==1)
00849           {
00850             s << " Scalars=\"" << (*it)->name() << "\"" ;
00851             break;
00852           }
00853         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
00854           if ((*it)->ncomps()>1)
00855           {
00856             s << " Vectors=\"" << (*it)->name() << "\"" ;
00857             break;
00858           }
00859         s << ">" << std::endl;
00860         indentUp();
00861         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
00862         {
00863           indent(s); s << "<PDataArray type=\"Float32\" Name=\"" << (*it)->name() << "\" ";
00864           s << "NumberOfComponents=\"" << ((*it)->ncomps()>1?3:1) << "\" ";
00865           if (outputtype==VTKOptions::ascii)
00866             s << "format=\"ascii\"/>" << std::endl;
00867           if (outputtype==VTKOptions::binary)
00868             s << "format=\"binary\"/>" << std::endl;
00869           if (outputtype==VTKOptions::binaryappended)
00870             s << "format=\"appended\"/>" << std::endl;
00871         }
00872         indentDown();
00873         indent(s); s << "</PPointData>" << std::endl;
00874 
00875         // PCellData
00876         indent(s); s << "<PCellData";
00877         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
00878           if ((*it)->ncomps()==1)
00879           {
00880             s << " Scalars=\"" << (*it)->name() << "\"" ;
00881             break;
00882           }
00883         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
00884           if ((*it)->ncomps()>1)
00885           {
00886             s << " Vectors=\"" << (*it)->name() << "\"" ;
00887             break;
00888           }
00889         s << ">" << std::endl;
00890         indentUp();
00891         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
00892         {
00893           indent(s); s << "<PDataArray type=\"Float32\" Name=\"" << (*it)->name() << "\" ";
00894           s << "NumberOfComponents=\"" << ((*it)->ncomps()>1?3:1) << "\" ";
00895           if (outputtype==VTKOptions::ascii)
00896             s << "format=\"ascii\"/>" << std::endl;
00897           if (outputtype==VTKOptions::binary)
00898             s << "format=\"binary\"/>" << std::endl;
00899           if (outputtype==VTKOptions::binaryappended)
00900             s << "format=\"appended\"/>" << std::endl;
00901         }
00902         indentDown();
00903         indent(s); s << "</PCellData>" << std::endl;
00904 
00905         // PPoints
00906         indent(s); s << "<PPoints>" << std::endl;
00907         indentUp();
00908         indent(s); s << "<PDataArray type=\"Float32\" Name=\"Coordinates\" NumberOfComponents=\"" << "3" << "\" ";
00909         if (outputtype==VTKOptions::ascii)
00910           s << "format=\"ascii\"/>" << std::endl;
00911         if (outputtype==VTKOptions::binary)
00912           s << "format=\"binary\"/>" << std::endl;
00913         if (outputtype==VTKOptions::binaryappended)
00914           s << "format=\"appended\"/>" << std::endl;
00915         indentDown();
00916         indent(s); s << "</PPoints>" << std::endl;
00917 
00918         // Pieces
00919         for (int i=0; i<grid.comm().size(); i++)
00920         {
00921           char fullname[128];
00922           if (n>1)
00923             sprintf(fullname,"%s/s%04d:p%0d:%s.vtu",piecepath,grid.comm().size(),i,piecename);
00924           else
00925             sprintf(fullname,"%s/s%04d:p%0d:%s.vtp",piecepath,grid.comm().size(),i,piecename);
00926           indent(s); s << "<Piece Source=\"" << fullname << "\"/>" << std::endl;
00927         }
00928 
00929         // /PUnstructuredGrid
00930         indentDown();
00931         indent(s); 
00932         if (n>1)
00933           s << "</PUnstructuredGrid>" << std::endl;
00934         else
00935           s << "</PPolyData>" << std::endl;
00936 
00937         // /VTKFile
00938         indentDown();
00939         s << "</VTKFile>" << std::endl;   
00940       }
00941 
00943     void writeDataFile (std::ostream& s)
00944       {
00945         // xml header
00946         s << "<?xml version=\"1.0\"?>" << std::endl;
00947 
00948         // VTKFile
00949         if (n>1)
00950           s << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
00951         else
00952           s << "<VTKFile type=\"PolyData\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
00953         indentUp();
00954 
00955         // Grid characteristics
00956         vertexmapper = new VertexMapper(grid,is);
00957         if (datamode == VTKOptions::conforming)
00958         {
00959           number.resize(vertexmapper->size());
00960           for (std::vector<int>::size_type i=0; i<number.size(); i++) number[i] = -1;
00961         }
00962         nvertices = 0;
00963         ncells = 0;
00964         ncorners = 0;
00965         for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
00966         {
00967           ncells++;
00968           for (int i=0; i<it->template count<n>(); ++i)
00969           {
00970             ncorners++;
00971             if (datamode == VTKOptions::conforming)
00972             {
00973               int alpha = vertexmapper->template map<n>(*it,i);
00974               if (number[alpha]<0)
00975                 number[alpha] = nvertices++;
00976             }
00977             else
00978             {
00979               nvertices++;
00980             }
00981           }
00982         }
00983 
00984         // UnstructuredGrid
00985         indent(s);
00986         if (n>1)
00987           s << "<UnstructuredGrid>" << std::endl;
00988         else
00989           s << "<PolyData>" << std::endl;
00990         indentUp();
00991 
00992         // Piece
00993         indent(s);
00994         if (n>1)
00995           s << "<Piece NumberOfPoints=\"" << nvertices << "\" NumberOfCells=\"" << ncells << "\">" << std::endl;
00996         else
00997           s << "<Piece NumberOfPoints=\"" << nvertices << "\""
00998             << " NumberOfVerts=\"0\""
00999             << " NumberOfLines=\"" << ncells << "\">" 
01000             << " NumberOfPolys=\"0\"" << std::endl;
01001         indentUp();
01002 
01003         // PointData
01004         writeVertexData(s);
01005 
01006         // CellData
01007         writeCellData(s);
01008 
01009         // Points
01010         writeGridPoints(s);
01011 
01012         // Cells
01013         writeGridCells(s);
01014 
01015         // /Piece
01016         indentDown();
01017         indent(s); s << "</Piece>" << std::endl;
01018 
01019         // /UnstructuredGrid
01020         indentDown();
01021         indent(s); 
01022         if (n>1)
01023           s << "</UnstructuredGrid>" << std::endl;
01024         else
01025           s << "</PolyData>" << std::endl;
01026 
01027         // write appended binary dat section
01028         if (outputtype==VTKOptions::binaryappended)
01029           writeAppendedData(s);
01030 
01031         // /VTKFile
01032         indentDown();
01033         s << "</VTKFile>" << std::endl;
01034 
01035         delete vertexmapper; number.clear();
01036       }
01037 
01038     void writeCellData (std::ostream& s)
01039       {
01040         indent(s); s << "<CellData";
01041         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
01042           if ((*it)->ncomps()==1)
01043           {
01044             s << " Scalars=\"" << (*it)->name() << "\"" ;
01045             break;
01046           }
01047         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
01048           if ((*it)->ncomps()>1)
01049           {
01050             s << " Vectors=\"" << (*it)->name() << "\"" ;
01051             break;
01052           }
01053         s << ">" << std::endl;
01054         indentUp();
01055         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
01056         {
01057           VTKDataArrayWriter<float> *p=0;
01058           if (outputtype==VTKOptions::ascii)
01059             p = new VTKAsciiDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps());
01060           if (outputtype==VTKOptions::binary)       
01061             p = new VTKBinaryDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps(),(*it)->ncomps()*ncells); 
01062           if (outputtype==VTKOptions::binaryappended)       
01063             p = new VTKBinaryAppendedDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps(),bytecount); 
01064           for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
01065             for (int j=0; j<(*it)->ncomps(); j++)
01066               p->write((*it)->evaluate(j,*i,i.position()));
01067           delete p;
01068         }
01069         indentDown();
01070         indent(s); s << "</CellData>" << std::endl;
01071       }
01072 
01073     void writeVertexData (std::ostream& s)
01074       {
01075         indent(s); s << "<PointData";
01076         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
01077           if ((*it)->ncomps()==1)
01078           {
01079             s << " Scalars=\"" << (*it)->name() << "\"" ;
01080             break;
01081           }
01082         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
01083           if ((*it)->ncomps()>1)
01084           {
01085             s << " Vectors=\"" << (*it)->name() << "\"" ;
01086             break;
01087           }
01088         s << ">" << std::endl;
01089         indentUp();
01090         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
01091         {
01092           VTKDataArrayWriter<float> *p=0;
01093           if (outputtype==VTKOptions::ascii)
01094             p = new VTKAsciiDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps()); 
01095           if (outputtype==VTKOptions::binary)
01096             p = new VTKBinaryDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps(),(*it)->ncomps()*nvertices); 
01097           if (outputtype==VTKOptions::binaryappended)
01098             p = new VTKBinaryAppendedDataArrayWriter<float>(s,(*it)->name(),(*it)->ncomps(),bytecount);
01099           for (VertexIterator vit=vertexBegin(); vit!=vertexEnd(); ++vit)
01100           {
01101             for (int j=0; j<(*it)->ncomps(); j++)
01102               p->write((*it)->evaluate(j,*vit,vit.position()));
01103                         //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01104                         if((*it)->ncomps()==2)
01105                           p->write(0.0);
01106           }
01107           delete p;
01108         }
01109         indentDown();
01110         indent(s); s << "</PointData>" << std::endl;
01111       }
01112 
01113     void writeGridPoints (std::ostream& s)
01114       {
01115         indent(s); s << "<Points>" << std::endl;
01116         indentUp();
01117 
01118         VTKDataArrayWriter<float> *p=0;
01119         if (outputtype==VTKOptions::ascii)
01120           p = new VTKAsciiDataArrayWriter<float>(s,"Coordinates",3);
01121         if (outputtype==VTKOptions::binary)
01122           p = new VTKBinaryDataArrayWriter<float>(s,"Coordinates",3,3*nvertices);
01123         if (outputtype==VTKOptions::binaryappended)
01124           p = new VTKBinaryAppendedDataArrayWriter<float>(s,"Coordinates",3,bytecount);
01125         VertexIterator vEnd = vertexEnd();
01126         for (VertexIterator vit=vertexBegin(); vit!=vEnd; ++vit)
01127         {
01128           int dimw=w;
01129           for (int j=0; j<std::min(dimw,3); j++)
01130             p->write(vit->geometry()[vit.localindex()][j]);
01131           for (int j=std::min(dimw,3); j<3; j++)
01132             p->write(0.0);
01133         }
01134         delete p;
01135       
01136         indentDown();
01137         indent(s); s << "</Points>" << std::endl;
01138       }
01139 
01140     void writeGridCells (std::ostream& s)
01141       {
01142         indent(s); 
01143         if (n>1)
01144           s << "<Cells>" << std::endl;
01145         else
01146           s << "<Lines>" << std::endl;
01147         indentUp();
01148 
01149         // connectivity
01150         VTKDataArrayWriter<int> *p1=0;
01151         if (outputtype==VTKOptions::ascii)
01152           p1 = new VTKAsciiDataArrayWriter<int>(s,"connectivity",1); 
01153         if (outputtype==VTKOptions::binary)
01154           p1 = new VTKBinaryDataArrayWriter<int>(s,"connectivity",1,ncorners); 
01155         if (outputtype==VTKOptions::binaryappended)
01156           p1 = new VTKBinaryAppendedDataArrayWriter<int>(s,"connectivity",1,bytecount);
01157         for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
01158           p1->write(it.id());
01159         delete p1;
01160 
01161         // offsets
01162         VTKDataArrayWriter<int> *p2=0;
01163         if (outputtype==VTKOptions::ascii)
01164           p2 = new VTKAsciiDataArrayWriter<int>(s,"offsets",1);
01165         if (outputtype==VTKOptions::binary)
01166           p2 = new VTKBinaryDataArrayWriter<int>(s,"offsets",1,ncells); 
01167         if (outputtype==VTKOptions::binaryappended)
01168           p2 = new VTKBinaryAppendedDataArrayWriter<int>(s,"offsets",1,bytecount);
01169         {
01170           int offset = 0;
01171           for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
01172           {
01173             offset += it->template count<n>();
01174             p2->write(offset);
01175           }
01176         }
01177         delete p2;
01178 
01179         // types
01180         if (n>1)
01181         {
01182           VTKDataArrayWriter<unsigned char> *p3=0; 
01183           if (outputtype==VTKOptions::ascii)
01184             p3 = new VTKAsciiDataArrayWriter<unsigned char>(s,"types",1);
01185           if (outputtype==VTKOptions::binary)
01186             p3 = new VTKBinaryDataArrayWriter<unsigned char>(s,"types",1,ncells); 
01187           if (outputtype==VTKOptions::binaryappended)
01188             p3 = new VTKBinaryAppendedDataArrayWriter<unsigned char>(s,"types",1,bytecount); 
01189           for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
01190           {
01191             int vtktype = vtkType(it->type());
01192             p3->write(vtktype);
01193           }
01194           delete p3;
01195         }
01196 
01197         indentDown();
01198         indent(s); 
01199         if (n>1)
01200           s << "</Cells>" << std::endl;
01201         else
01202           s << "</Lines>" << std::endl;
01203       }
01204 
01205 
01206     void writeAppendedData (std::ostream& s)
01207       {
01208         indent(s); s << "<AppendedData encoding=\"raw\">" << std::endl;
01209         indentUp();
01210         indent(s); s << "_"; // indicates start of binary data
01211 
01212         SimpleStream stream(s);
01213 
01214         // write length before each data block
01215         unsigned int blocklength;
01216 
01217         // point data     
01218         for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
01219         {
01220                   
01221                   blocklength = nvertices * (*it)->ncomps() * sizeof(float);
01222                   //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01223                   if((*it)->ncomps()==2)
01224                         blocklength = nvertices * (3) * sizeof(float);
01225           stream.write(blocklength);
01226           std::vector<bool> visited(vertexmapper->size(), false);
01227           for (VertexIterator vit=vertexBegin(); vit!=vertexEnd(); ++vit)
01228           {
01229             for (int j=0; j<(*it)->ncomps(); j++)
01230             {
01231               float data = (*it)->evaluate(j,*vit,vit.position());
01232               stream.write(data);
01233             }
01234                         //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01235                         if((*it)->ncomps()==2){
01236                           float data=0.0;
01237                           stream.write(data);}
01238           }
01239         }
01240 
01241         // cell data
01242         for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
01243         {
01244           blocklength = ncells * (*it)->ncomps() * sizeof(float);
01245           stream.write(blocklength);
01246           for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
01247             for (int j=0; j<(*it)->ncomps(); j++)
01248             {
01249               float data = (*it)->evaluate(j,*i,i.position());
01250               stream.write(data);
01251             }
01252         }
01253       
01254         // point coordinates
01255         blocklength = nvertices * 3 * sizeof(float);
01256         stream.write(blocklength);
01257         std::vector<bool> visited(vertexmapper->size(), false);
01258         for (VertexIterator vit=vertexBegin(); vit!=vertexEnd(); ++vit)
01259         {
01260           int dimw=w;
01261           float data;
01262           for (int j=0; j<std::min(dimw,3); j++)
01263           {
01264             data = vit->geometry()[vit.localindex()][j];
01265             stream.write(data);
01266           }
01267           data = 0;
01268           for (int j=std::min(dimw,3); j<3; j++)
01269             stream.write(data);
01270         }
01271       
01272         // connectivity
01273         blocklength = ncorners * sizeof(unsigned int);
01274         stream.write(blocklength);
01275         for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
01276         {
01277           stream.write(it.id());
01278         }
01279 
01280         // offsets
01281         blocklength = ncells * sizeof(unsigned int);
01282         stream.write(blocklength);
01283         {
01284           int offset = 0;
01285           for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
01286           {
01287             offset += it->template count<n>();
01288             stream.write(offset);
01289           }
01290         }
01291 
01292         // cell types
01293         if (n>1)
01294         {
01295           blocklength = ncells * sizeof(unsigned char);
01296           stream.write(blocklength);
01297           for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
01298           {
01299             unsigned char vtktype = vtkType(it->type());
01300             stream.write(vtktype);
01301           }
01302         }
01303 
01304         s << std::endl;
01305         indentDown();
01306         indent(s); s << "</AppendedData>" << std::endl;
01307       }
01308     
01309     // base class for data array writers
01310     template<class T>
01311     class VTKDataArrayWriter
01312     {
01313     public:
01314       virtual void write (T data) = 0;
01315       virtual ~VTKDataArrayWriter () {}
01316     };
01317 
01318     // a streaming writer for data array tags
01319     template<class T>
01320     class VTKAsciiDataArrayWriter : public VTKDataArrayWriter<T>
01321     {
01322     public:
01324       VTKAsciiDataArrayWriter (std::ostream& theStream, std::string name, int ncomps) 
01325         : s(theStream), counter(0), numPerLine(12)
01326         {
01327           VTKTypeNameTraits<T> tn;
01328           s << "<DataArray type=\"" << tn() << "\" Name=\"" << name << "\" ";
01329                   //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01330           if (ncomps>3)
01331             DUNE_THROW(IOError, "VTKWriter does not support more than 3 components");
01332           s << "NumberOfComponents=\"" << (ncomps>1?3:1) << "\" ";
01333           s << "format=\"ascii\">" << std::endl;
01334         }
01335 
01337       void write (T data)
01338         {
01339           typedef typename VTKTypeNameTraits<T>::PrintType PT;
01340           s << (PT) data << " ";
01341           counter++;
01342           if (counter%numPerLine==0) s << std::endl;
01343         }
01344 
01346       ~VTKAsciiDataArrayWriter ()
01347         {
01348           if (counter%numPerLine!=0) s << std::endl;
01349           s << "</DataArray>" << std::endl;   
01350         }
01351 
01352     private:
01353       std::ostream& s;
01354       int counter;
01355       int numPerLine;
01356     };
01357 
01358     // a streaming writer for data array tags
01359     template<class T>
01360     class VTKBinaryDataArrayWriter : public VTKDataArrayWriter<T>
01361     {
01362     public:
01364       VTKBinaryDataArrayWriter (std::ostream& theStream, std::string name, int ncomps, int nitems) 
01365         : s(theStream),bufsize(4096),n(0)
01366         {
01367           DUNE_THROW(IOError, "binary does not work yet, use binaryappended!");
01368           VTKTypeNameTraits<T> tn;
01369           s << "<DataArray type=\"" << tn() << "\" Name=\"" << name << "\" ";
01370                   //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01371           if (ncomps>3)
01372             DUNE_THROW(IOError, "VTKWriter does not support more than 3 components");
01373           s << "NumberOfComponents=\"" << (ncomps>1?3:1) << "\" ";
01374           s << "format=\"binary\">" << std::endl;
01375           buffer = new char[bufsize*sizeof(T)];
01376           code = new char[2*bufsize*sizeof(T)];
01377           unsigned int size = nitems*sizeof(T);
01378           char* p = reinterpret_cast<char*>(&size);
01379           memcpy(buffer+n,p,sizeof(int));
01380           n += sizeof(int);
01381           //        base64::base64_init_encodestate(&_state);
01382         }
01383 
01385       void write (T data)
01386         {
01387           if (n+sizeof(T)>bufsize)
01388           {
01389             // flush buffer
01390             //          int codelength = base64::base64_encode_block(buffer,n,code,&_state);
01391             //          s.write(code,codelength);
01392             n=0;
01393           }
01394           char* p = reinterpret_cast<char*>(&data);
01395           memcpy(buffer+n,p,sizeof(T));
01396           n += sizeof(T);
01397         }
01398 
01400       ~VTKBinaryDataArrayWriter ()
01401         {
01402 //      int codelength;
01403           if (n>0)
01404           {
01405             //          codelength = base64::base64_encode_block(buffer,n,code,&_state);
01406             //          s.write(code,codelength);
01407           }
01408           //        codelength = base64::base64_encode_blockend(code,&_state);
01409 //      s.write(code,codelength);
01410           //        base64::base64_init_encodestate(&_state);
01411           s << std::endl;
01412           s << "</DataArray>" << std::endl;
01413           delete [] code;
01414           delete [] buffer;
01415         }
01416 
01417     private:
01418       std::ostream& s;
01419       //      base64::base64_encodestate _state;
01420       size_t bufsize;
01421       char* buffer;
01422       char* code;
01423       int n;
01424     };
01425 
01426     // a streaming writer for data array tags
01427     template<class T>
01428     class VTKBinaryAppendedDataArrayWriter : public VTKDataArrayWriter<T>
01429     {
01430     public:
01432       VTKBinaryAppendedDataArrayWriter (std::ostream& theStream, std::string name, int ncomps, unsigned int& bc) 
01433         : s(theStream),bytecount(bc)
01434         {
01435           VTKTypeNameTraits<T> tn;
01436           s << "<DataArray type=\"" << tn() << "\" Name=\"" << name << "\" ";
01437                   //vtk file format: a vector data always should have 3 comps(with 3rd comp = 0 in 2D case)
01438           if (ncomps>3)
01439             DUNE_THROW(IOError, "VTKWriter does not support more than 3 components");
01440           s << "NumberOfComponents=\"" << (ncomps>1?3:1) << "\" ";
01441           s << "format=\"appended\" offset=\""<< bytecount << "\" />" << std::endl;
01442           bytecount += 4; // header
01443         }
01444 
01446       void write (T data)
01447         {
01448           bytecount += sizeof(T);
01449         }
01450 
01451     private:
01452       std::ostream& s;
01453       unsigned int& bytecount;
01454     };
01455 
01456 
01457     // write out data in binary
01458     class SimpleStream
01459     {
01460     public:
01461       SimpleStream (std::ostream& theStream)
01462         : s(theStream)
01463         {}
01464       template<class T>
01465       void write (T data)
01466         {
01467           char* p = reinterpret_cast<char*>(&data);
01468           s.write(p,sizeof(T));
01469         }
01470     private:
01471       std::ostream& s;
01472     };
01473 
01474     void indentUp ()
01475       {
01476         indentCount++;
01477       }
01478 
01479     void indentDown ()
01480       {
01481         indentCount--;
01482       }
01483 
01484     void indent (std::ostream& s)
01485       {
01486         for (int i=0; i<indentCount; i++) 
01487           s << "  ";
01488       }
01489 
01490     // renumber VTK -> Dune
01491     static int renumber (const Entity& e, int i)
01492       {
01493         static const int quadRenumbering[4] = {0,1,3,2};
01494         static const int cubeRenumbering[8] = {0,1,3,2,4,5,7,6};
01495         static const int prismRenumbering[6] = {0,2,1,3,5,4};
01496         switch (vtkType(e.type()))
01497         {
01498         case vtkQuadrilateral:
01499           return quadRenumbering[i];
01500         case vtkHexahedron:
01501           return cubeRenumbering[i];
01502         case vtkPrism:
01503           return prismRenumbering[i];
01504         default:
01505           return i;
01506         }
01507       }
01508 
01509     // the list of registered functions
01510     std::list<VTKFunction*> celldata;
01511     std::list<VTKFunction*> vertexdata;
01512 
01513     // the grid
01514     const GridImp& grid;
01515 
01516     // the indexset
01517     const IndexSet& is;
01518 
01519     // intend counter
01520     int indentCount;
01521     int numPerLine;
01522 
01523     // temporary grid information
01524     int ncells;
01525     int nvertices;
01526     int ncorners;
01527     VertexMapper* vertexmapper;
01528     std::vector<int> number;
01529     VTKOptions::DataMode datamode;
01530     VTKOptions::OutputType outputtype;
01531     unsigned int bytecount;
01532   };
01533 
01537   template<class G>
01538   class LeafVTKWriter : public VTKWriter<G,typename G::template Codim<0>::LeafIndexSet>
01539   {
01540   public:
01542     LeafVTKWriter (const G& grid, VTKOptions::DataMode dm = VTKOptions::conforming)
01543       : VTKWriter<G,typename G::template Codim<0>::LeafIndexSet>(grid,grid.leafIndexSet(),dm)
01544       {}
01545   };
01546 
01550   template<class G>
01551   class LevelVTKWriter : public VTKWriter<G, typename G::template Codim<0>::LevelIndexSet>
01552   {
01553   public:
01555     LevelVTKWriter (const G& grid, int level, VTKOptions::DataMode dm = VTKOptions::conforming)
01556       : VTKWriter<G,typename G::template Codim<0>::LevelIndexSet>(grid,grid.levelIndexSet(level),dm)
01557       {}
01558   };
01559 }
01560 #endif

Generated on 9 Apr 2008 with Doxygen (ver 1.5.2) [logfile].