dune-grid  2.2.0
dgfoned.hh
Go to the documentation of this file.
00001 #ifndef DUNE_GRID_FILE_IO_DGFPARSER_DGFONED_HH
00002 #define DUNE_GRID_FILE_IO_DGFPARSER_DGFONED_HH
00003 
00004 //- C++ includes
00005 #include <algorithm>
00006 #include <fstream>
00007 #include <iostream>
00008 #include <istream>
00009 #include <vector>
00010 
00011 //- dune-common includes
00012 #include <dune/common/exceptions.hh>
00013 
00014 //- dune-grid includes
00015 #include <dune/grid/common/intersection.hh>
00016 #include <dune/grid/onedgrid.hh>
00017 
00018 //- local includes
00019 #include "dgfparser.hh"
00020 
00021 
00022 namespace
00023 {
00024   // helper method used below
00025   double getfirst ( std::vector< double > v )
00026   {
00027     return v[ 0 ];
00028   }
00029 }  // end anonymous namespace
00030 
00031 
00032 
00033 namespace Dune 
00034 {
00035 
00036   // DGFGridInfo
00037   // -----------
00038 
00039   template< >
00040   struct DGFGridInfo< OneDGrid > 
00041   {
00042     static int refineStepsForHalf ()
00043     {
00044       return 1;
00045     }
00046 
00047     static double refineWeight ()
00048     {
00049       return 0.5;
00050     }
00051   };
00052 
00053 
00054 
00055   // DGFGridFactory< OneDGrid >
00056   // --------------------------
00057 
00058   template< >
00059   struct DGFGridFactory< OneDGrid >
00060   {
00062     typedef OneDGrid Grid;
00064     const static int dimension = Grid::dimension;
00066     typedef MPIHelper::MPICommunicator MPICommunicatorType;
00067 
00069     explicit DGFGridFactory ( std::istream &input,
00070                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00071     : grid_( 0 ),
00072       emptyParameters_( 0 )
00073     {
00074       generate( input, comm );
00075     }
00076 
00078     explicit DGFGridFactory ( const std::string &filename,
00079                               MPICommunicatorType comm = MPIHelper::getCommunicator() )
00080     : grid_( 0 ),
00081       emptyParameters_( 0 )
00082     {
00083       std::ifstream input( filename.c_str() );
00084       generate( input, comm );
00085     }
00086 
00088     Grid *grid () const
00089     {
00090       return grid_;
00091     }
00092 
00094     template < class GG, template< class > class II >
00095     bool wasInserted ( const Dune::Intersection< GG, II > &intersection ) const
00096     {
00097       return false;
00098     }
00099 
00100     template < class GG, template< class > class II >
00101     int boundaryId ( const Dune::Intersection< GG, II > &intersection ) const
00102     {
00103       // OneDGrid returns boundary segment index; 
00104       // we return the index as the method boundaryId is deprecated
00105       return intersection.boundarySegmentIndex();
00106     }
00107 
00109     template< class Entity >
00110     int numParameters ( const Entity & ) const
00111     {
00112       return 0;
00113     }
00114 
00116     template< int codim >
00117     int numParameters () const
00118     {
00119       return 0;
00120     }
00121 
00122     template< class Entity >
00123     std::vector< double >& parameter ( const Entity &entity )
00124     {
00125       return parameter< Entity::codimension >( entity );
00126     }
00127 
00129     template< int codim >
00130     std::vector< double > &parameter ( const typename Grid::Codim< codim >::Entity &element )
00131     {
00132       return emptyParameters_;
00133     }
00134 
00136     bool haveBoundaryParameters () const
00137     {
00138       return false;
00139     }
00140 
00142     template < class GG, template< class > class II >
00143     const DGFBoundaryParameter::type &boundaryParameter ( const Dune::Intersection< GG, II > &intersection ) const
00144     {
00145       return DGFBoundaryParameter::defaultValue();
00146     }
00147 
00148   private:
00149     // generate grid
00150     void generate ( std::istream &input, MPICommunicatorType comm );
00151 
00152     Grid *grid_;
00153     std::vector< double > emptyParameters_;
00154   };
00155 
00156 
00157 
00158   // Implementation of DGFGridFactory< OneDGrid >
00159   // --------------------------------------------
00160 
00161   inline void DGFGridFactory< OneDGrid >::generate ( std::istream &input, MPICommunicatorType comm )
00162   {
00163     // try to open interval block 
00164     dgf::IntervalBlock intervalBlock( input );
00165 
00166     // try to open vertex block 
00167     int dimensionworld = Grid::dimensionworld;
00168     dgf::VertexBlock vertexBlock( input, dimensionworld );
00169 
00170     // check at least one block is active
00171     if( !( vertexBlock.isactive() || intervalBlock.isactive() ))
00172       DUNE_THROW( DGFException, "No readable block found" );
00173     
00174     std::vector< std::vector< double > > vertices;
00175 
00176     // read vertices first
00177     if( vertexBlock.isactive() )
00178     {
00179       int nparameter = 0;
00180       std::vector< std::vector< double > > parameter;
00181       vertexBlock.get( vertices, parameter, nparameter );
00182 
00183       if( nparameter > 0 )
00184         std::cerr << "Warning: vertex parameters will be ignored" << std::endl;
00185     }
00186 
00187     // get vertices from interval block
00188     if ( intervalBlock.isactive() )
00189     {
00190       if( intervalBlock.dimw() != dimensionworld )
00191       {
00192         DUNE_THROW( DGFException, "Error: wrong coordinate dimension in interval block \
00193                                    (got " << intervalBlock.dimw() << ", expected " << dimensionworld << ")" );
00194       }
00195 
00196       int nintervals = intervalBlock.numIntervals();
00197       for( int i = 0; i < nintervals; ++i )
00198         intervalBlock.getVtx( i, vertices );
00199     }
00200 
00201     // copy to vector of doubles
00202     std::vector< double > vtx( vertices.size() );
00203     transform( vertices.begin(), vertices.end(), vtx.begin(), getfirst );
00204 
00205     // remove duplicates
00206     std::sort( vtx.begin(), vtx.end() );
00207     std::vector< double >::iterator it = std::unique( vtx.begin(), vtx.end() );
00208     vtx.erase( it, vtx.end() );
00209     if( vertices.size() != vtx.size() )
00210       std::cerr << "Warning: removed duplicate vertices" << std::endl;
00211 
00212     // create grid
00213     grid_ = new OneDGrid( vtx );
00214   }
00215 
00216 } // end namespace Dune
00217 
00218 #endif // #ifndef DUNE_GRID_FILE_IO_DGFPARSER_DGFONED_HH