- Home
- About DUNE
- Download
- Documentation
- Community
- Development
00001 #ifndef DUNE_AMIRAMESH_WRITER_HH 00002 #define DUNE_AMIRAMESH_WRITER_HH 00003 00004 #include <string> 00005 00006 #include <dune/common/array.hh> 00007 00008 #if HAVE_AMIRAMESH 00009 #include <amiramesh/AmiraMesh.h> 00010 #endif 00011 00012 namespace Dune { 00013 00018 template<class GridView> 00019 class AmiraMeshWriter { 00020 00021 enum {dim = GridView::dimension}; 00022 00023 public: 00024 00033 void addGrid(const GridView& gridView, bool splitAll=false); 00034 00044 template <class GridType2> 00045 void addLevelGrid(const GridType2& grid, int level, bool splitAll=false); 00046 00055 template <class GridType2> 00056 void addLeafGrid(const GridType2& grid, bool splitAll=false); 00057 00064 template <class DataContainer> 00065 void addCellData(const DataContainer& data, const GridView& gridView, bool GridSplitUp=false); 00066 00073 template <class DataContainer> 00074 void addVertexData(const DataContainer& data, const GridView& gridView, bool GridSplitUp=false); 00075 00080 void write(const std::string& filename, bool ascii=false) const; 00081 00084 template <class DataContainer> 00085 void addUniformData(const GridView& gridView, 00086 const array<unsigned int, dim>& n, 00087 const DataContainer& data); 00088 00089 protected: 00090 00091 #if HAVE_AMIRAMESH // better: use a pointer here and forward-declare AmiraMesh 00092 AmiraMesh amiramesh_; 00093 #endif 00094 }; 00095 00100 template<class GridType> 00101 class LevelAmiraMeshWriter 00102 : public AmiraMeshWriter<typename GridType::LevelGridView> 00103 { 00104 00105 public: 00106 00108 LevelAmiraMeshWriter() {} 00109 00111 LevelAmiraMeshWriter(const GridType& grid, int level) { 00112 this->addGrid(grid.levelView(level)); 00113 } 00114 00121 static void writeGrid(const GridType& grid, 00122 const std::string& filename, 00123 int level) { 00124 LevelAmiraMeshWriter amiramesh(grid, level); 00125 amiramesh.write(filename); 00126 } 00127 00136 template <class VectorType> 00137 static void writeBlockVector(const GridType& grid, 00138 const VectorType& f, 00139 const std::string& filename, 00140 int level, 00141 bool GridSplitUp=false) { 00142 LevelAmiraMeshWriter amiramesh; 00143 if (f.size()==grid.size(level,GridType::dimension)) 00144 amiramesh.addVertexData(f, grid.levelView(level),GridSplitUp); 00145 else 00146 amiramesh.addCellData(f, grid.levelView(level),GridSplitUp); 00147 amiramesh.write(filename); 00148 } 00149 00150 }; 00151 00156 template<class GridType> 00157 class LeafAmiraMeshWriter 00158 : public AmiraMeshWriter<typename GridType::LeafGridView> 00159 { 00160 00161 public: 00162 00164 LeafAmiraMeshWriter() {} 00165 00167 LeafAmiraMeshWriter(const GridType& grid) { 00168 this->addLeafGrid(grid); 00169 } 00170 00176 static void writeGrid(const GridType& grid, 00177 const std::string& filename) { 00178 LeafAmiraMeshWriter amiramesh(grid); 00179 amiramesh.write(filename); 00180 } 00181 00188 template <class VectorType> 00189 static void writeBlockVector(const GridType& grid, 00190 const VectorType& f, 00191 const std::string& filename, 00192 bool GridSplitUp = false) { 00193 LeafAmiraMeshWriter amiramesh; 00194 if (f.size()==grid.size(GridType::dimension)) 00195 amiramesh.addVertexData(f, grid.leafView(),GridSplitUp); 00196 else 00197 amiramesh.addCellData(f, grid.leafView(),GridSplitUp); 00198 00199 amiramesh.write(filename); 00200 } 00201 00202 }; 00203 00204 } 00205 00206 // implementation 00207 #if HAVE_AMIRAMESH 00208 #include "amiramesh/amirameshwriter.cc" 00209 #endif 00210 00211 #endif
Generated on Fri Apr 29 2011 with Doxygen (ver 1.7.1) [doxygen-log,error-log].