dune-grid
2.2.0
|
00001 #ifndef DUNE_DGF_GRIDFACTORY_HH 00002 #define DUNE_DGF_GRIDFACTORY_HH 00003 00004 #include <iostream> 00005 #include <string> 00006 #include <vector> 00007 #include <map> 00008 #include <assert.h> 00009 00010 #include <dune/common/mpihelper.hh> 00011 #include <dune/grid/io/file/dgfparser/dgfexception.hh> 00012 #include <dune/grid/io/file/dgfparser/macrogrid.hh> 00013 00014 #include <dune/grid/io/file/dgfparser/parser.hh> 00015 #include <dune/grid/common/intersection.hh> 00016 00017 00018 namespace Dune 00019 { 00020 00021 // forward declaration 00022 template < class GridImp, template < class > class IntersectionImp > 00023 class Intersection; 00024 00025 template < class G > 00026 struct DGFGridFactory 00027 { 00028 typedef G Grid; 00029 const static int dimension = Grid::dimension; 00030 typedef MPIHelper::MPICommunicator MPICommunicatorType; 00031 00032 private: 00033 typedef typename Grid::template Codim< 0 >::Entity Element; 00034 00035 typedef typename Grid::template Codim< dimension >::Entity Vertex; 00036 00037 public: 00038 explicit DGFGridFactory ( std::istream &input, 00039 MPICommunicatorType comm = MPIHelper::getCommunicator() ) DUNE_DEPRECATED 00040 : macroGrid_( comm ) 00041 { 00042 DUNE_THROW( DGFException, "DGF factories using old MacroGrid implementation" 00043 "don't support creation from std::istream." ); 00044 } 00045 00046 explicit DGFGridFactory ( const std::string &filename, 00047 MPICommunicatorType comm = MPIHelper::getCommunicator() ) DUNE_DEPRECATED 00048 : macroGrid_( filename.c_str(), comm ) 00049 { 00050 grid_ = macroGrid_.template createGrid< Grid >(); 00051 00052 if( macroGrid_.nofelparams > 0 ) 00053 { 00054 const size_t nofElements = macroGrid_.elements.size(); 00055 for( size_t i = 0; i < nofElements; ++i ) 00056 { 00057 std::vector< double > coord; 00058 00059 DomainType p(0); 00060 const size_t nofCorners = macroGrid_.elements[i].size(); 00061 for (size_t k=0;k<nofCorners;++k) 00062 for (int j=0;j<DomainType::dimension;++j) 00063 p[j]+=macroGrid_.vtx[macroGrid_.elements[i][k]][j]; 00064 p/=double(nofCorners); 00065 00066 elInsertOrder_.insert( std::make_pair( p, i ) ); 00067 } 00068 } 00069 00070 if( macroGrid_.nofvtxparams > 0 ) 00071 { 00072 const size_t nofVertices = macroGrid_.vtx.size(); 00073 for( size_t i = 0; i < nofVertices; ++i ) 00074 { 00075 std::vector< double > coord; 00076 00077 DomainType p; 00078 for( int k = 0; k < DomainType::dimension; ++k ) 00079 p[ k ] = macroGrid_.vtx[i][k]; 00080 00081 vtxInsertOrder_.insert( std::make_pair( p, i ) ); 00082 } 00083 } 00084 } 00085 00086 Grid *grid() 00087 { 00088 return grid_; 00089 } 00090 00091 template <class Intersection> 00092 bool wasInserted(const Intersection &intersection) const 00093 { 00094 return intersection.boundary(); 00095 } 00096 00097 template <class Intersection> 00098 int boundaryId(const Intersection &intersection) const 00099 { 00100 return intersection.boundaryId(); 00101 } 00102 00103 template< int codim > 00104 int numParameters () const 00105 { 00106 if( codim == 0 ) 00107 return macroGrid_.nofelparams; 00108 else if( codim == dimension ) 00109 return macroGrid_.nofvtxparams; 00110 else 00111 return 0; 00112 } 00113 00114 template < class Entity > 00115 int numParameters ( const Entity & ) const 00116 { 00117 return numParameters< Entity::codimension >(); 00118 } 00119 00120 std::vector<double>& parameter(const Element &element) 00121 { 00122 const typename Element::Geometry &geo = element.geometry(); 00123 DomainType coord( geo.corner( 0 ) ); 00124 for( int i = 1; i < geo.corners(); ++i ) 00125 coord += geo.corner( i ); 00126 coord /= double( geo.corners() ); 00127 00128 InsertOrderIterator it = elInsertOrder_.find( coord ); 00129 if( it != elInsertOrder_.end() ) 00130 return macroGrid_.elParams[ it->second ]; 00131 assert(0); 00132 return emptyParam; 00133 } 00134 00135 std::vector<double>& parameter(const Vertex &vertex) 00136 { 00137 const typename Vertex::Geometry &geo = vertex.geometry(); 00138 DomainType coord( geo.corner( 0 ) ); 00139 00140 InsertOrderIterator it = vtxInsertOrder_.find( coord ); 00141 if( it != vtxInsertOrder_.end() ) 00142 return macroGrid_.vtxParams[ it->second ]; 00143 return emptyParam; 00144 } 00145 00146 // return true if boundary parameters found 00147 bool haveBoundaryParameters () const 00148 { 00149 return false; 00150 } 00151 00152 template < class GG, template < class > class II > 00153 const typename DGFBoundaryParameter::type & 00154 boundaryParameter ( const Intersection< GG, II > & intersection ) const 00155 { 00156 return DGFBoundaryParameter::defaultValue(); 00157 } 00158 00159 private: 00160 typedef FieldVector<typename Grid::ctype,Grid::dimensionworld> DomainType; 00161 struct Compare 00162 { 00163 bool operator() ( const DomainType &a, const DomainType &b ) const 00164 { 00165 // returns true, if a < b; c[i] < -eps; 00166 const DomainType c = a - b; 00167 const double eps = 1e-8; 00168 00169 for( int i = 0; i < DomainType::dimension; ++i ) 00170 { 00171 if( c[ i ] <= -eps ) 00172 return true; 00173 if( c[ i ] >= eps ) 00174 return false; 00175 } 00176 return false; 00177 } 00178 }; 00179 typedef std::map< DomainType, size_t, Compare > InsertOrderMap; 00180 typedef typename InsertOrderMap::const_iterator InsertOrderIterator; 00181 00182 MacroGrid macroGrid_; 00183 Grid *grid_; 00184 InsertOrderMap elInsertOrder_; 00185 InsertOrderMap vtxInsertOrder_; 00186 std::vector<double> emptyParam; 00187 }; 00188 00189 } // end namespace Dune 00190 00191 #endif