• Main Page
  • Namespaces
  • Data Structures
  • Files
  • File List
  • Globals

/export/development/ViennaGrid/release/ViennaGrid-1.0.0/viennagrid/algorithm/voronoi.hpp

Go to the documentation of this file.
00001 #ifndef VIENNAGRID_ALGORITHM_VORONOI_HPP
00002 #define VIENNAGRID_ALGORITHM_VORONOI_HPP
00003 
00004 /* =======================================================================
00005    Copyright (c) 2011, Institute for Microelectronics,
00006                        Institute for Analysis and Scientific Computing,
00007                        TU Wien.
00008 
00009                             -----------------
00010                      ViennaGrid - The Vienna Grid Library
00011                             -----------------
00012 
00013    Authors:      Karl Rupp                           rupp@iue.tuwien.ac.at
00014                  Josef Weinbub                    weinbub@iue.tuwien.ac.at
00015                
00016    (A list of additional contributors can be found in the PDF manual)
00017 
00018    License:      MIT (X11), see file LICENSE in the base directory
00019 ======================================================================= */
00020 
00021 #include <iostream>
00022 #include <sstream>
00023 #include <string>
00024 #include <stdexcept>
00025 
00026 #include "viennagrid/forwards.h"
00027 #include "viennagrid/algorithm/circumcenter.hpp"
00028 #include "viennagrid/algorithm/spanned_volume.hpp"
00029 #include "viennagrid/algorithm/volume.hpp"
00030 
00031 #include "viennadata/api.hpp"
00032 
00038 namespace viennagrid
00039 {
00040   namespace detail
00041   {
00043     template <typename DomainType,
00044               typename InterfaceAreaKey,
00045               typename BoxVolumeKey>
00046     void write_voronoi_info(DomainType const & domain,
00047                             InterfaceAreaKey const & interface_key,
00048                             BoxVolumeKey const & box_volume_key,
00049                             viennagrid::simplex_tag<1>)
00050     {
00051       typedef typename DomainType::config_type           Config;
00052       typedef typename Config::cell_tag                  CellTag;
00053       
00054       typedef typename viennagrid::result_of::point<Config>::type                            PointType;
00055       typedef typename viennagrid::result_of::ncell<Config, 0>::type                         VertexType;
00056       typedef typename viennagrid::result_of::ncell<Config, CellTag::dim>::type              CellType;
00057       
00058       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim>::type    CellRange;
00059       typedef typename viennagrid::result_of::iterator<CellRange>::type                            CellIterator;
00060       
00061       
00062       typedef typename viennagrid::result_of::const_ncell_range<CellType, 0>::type           VertexOnCellRange;
00063       typedef typename viennagrid::result_of::iterator<VertexOnCellRange>::type              VertexOnCellIterator;
00064       
00065       CellRange cells = viennagrid::ncells<CellTag::dim>(domain);
00066       for (CellIterator cit  = cells.begin();
00067                         cit != cells.end();
00068                       ++cit)
00069       {
00070         viennadata::access<InterfaceAreaKey, double>(interface_key)(*cit) = 1;
00071         viennadata::access<BoxVolumeKey, double>(box_volume_key)(*cit) += volume(*cit) / 2.0;
00072         
00073         VertexOnCellRange vertices_on_cell = viennagrid::ncells<0>(*cit);
00074         for (VertexOnCellIterator vocit  = vertices_on_cell.begin();
00075                                   vocit != vertices_on_cell.end();
00076                                 ++vocit)
00077         {
00078           viennadata::access<BoxVolumeKey, double>(box_volume_key)(*vocit) += volume(*cit) / 2.0;
00079         }
00080       }
00081     }
00082     
00084     template <typename DomainType,
00085               typename InterfaceAreaKey,
00086               typename BoxVolumeKey>
00087     void write_voronoi_info(DomainType const & domain,
00088                             InterfaceAreaKey const & interface_key,
00089                             BoxVolumeKey const & box_volume_key,
00090                             viennagrid::hypercube_tag<1>)
00091     {
00092       write_voronoi_info(domain, interface_key, box_volume_key, viennagrid::simplex_tag<1>());
00093     }
00094     
00095     //
00096     // Voronoi information in two (topological) dimensions
00097     //
00099     template <typename DomainType,
00100               typename InterfaceAreaKey,
00101               typename BoxVolumeKey>
00102     void write_voronoi_info(DomainType const & domain,
00103                             InterfaceAreaKey const & interface_key,
00104                             BoxVolumeKey const & box_volume_key,
00105                             viennagrid::quadrilateral_tag)
00106     {
00107       //std::cout << "Warning: Voronoi info for quadrilaterals is only correct when having rectangles only." << std::endl;
00108       typedef typename DomainType::config_type           Config;
00109       typedef typename Config::cell_tag                  CellTag;
00110       typedef typename viennagrid::result_of::point<Config>::type                            PointType;
00111       typedef typename viennagrid::result_of::ncell<Config, 0>::type                         VertexType;
00112       typedef typename viennagrid::result_of::ncell<Config, 1>::type                         EdgeType;
00113       typedef typename viennagrid::result_of::ncell<Config, CellTag::dim>::type   CellType;
00114       
00115       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim>::type    CellRange;
00116       typedef typename viennagrid::result_of::iterator<CellRange>::type                                       CellIterator;
00117 
00118       typedef typename viennagrid::result_of::const_ncell_range<CellType, 1>::type                            EdgeOnCellRange;
00119       typedef typename viennagrid::result_of::iterator<EdgeOnCellRange>::type                                 EdgeOnCellIterator;
00120       
00121       typedef typename viennagrid::result_of::const_ncell_range<EdgeType, 0>::type                            VertexOnEdgeRange;
00122       typedef typename viennagrid::result_of::iterator<VertexOnEdgeRange>::type                               VertexOnEdgeIterator;
00123 
00124       //
00125       // Algorithm: Iterate over all cells, compute circumcenter and add interface area to edge, box volume to vertex.
00126       //
00127       
00128       CellRange cells = viennagrid::ncells<CellTag::dim>(domain);
00129       for (CellIterator cit  = cells.begin();
00130                         cit != cells.end();
00131                       ++cit)
00132       {
00133         PointType circ_center = circumcenter(*cit);
00134         
00135         //iterate over edges:
00136         EdgeOnCellRange edges_on_cell = viennagrid::ncells<1>(*cit);
00137         for (EdgeOnCellIterator eocit  = edges_on_cell.begin();
00138                                 eocit != edges_on_cell.end();
00139                               ++eocit)
00140         {
00141           PointType edge_midpoint = circumcenter(*eocit);
00142 
00143           // interface contribution:
00144           viennadata::access<InterfaceAreaKey, double>(interface_key)(*eocit) += spanned_volume(circ_center, edge_midpoint);
00145           
00146           //box volume contribution:
00147           double edge_contribution = 0;
00148           VertexOnEdgeRange vertices_on_edge = viennagrid::ncells<0>(*eocit);
00149           for (VertexOnEdgeIterator voeit  = vertices_on_edge.begin();
00150                                     voeit != vertices_on_edge.end();
00151                                   ++voeit)
00152           {
00153             double contribution = spanned_volume(circ_center, edge_midpoint, voeit->point());
00154             edge_contribution += contribution;
00155             viennadata::access<BoxVolumeKey, double>(box_volume_key)(*voeit) += contribution;
00156           }
00157           viennadata::access<BoxVolumeKey, double>(box_volume_key)(*eocit) += edge_contribution;
00158         } //for edges on cells
00159         
00160       } //for cells
00161       
00162     } //write_voronoi_info(triangle_tag)
00163 
00164 
00166     template <typename DomainType,
00167               typename InterfaceAreaKey,
00168               typename BoxVolumeKey>
00169     void write_voronoi_info(DomainType const & domain,
00170                             InterfaceAreaKey const & interface_key,
00171                             BoxVolumeKey const & box_volume_key,
00172                             viennagrid::triangle_tag)
00173     {
00174       typedef typename DomainType::config_type           Config;
00175       typedef typename Config::cell_tag                  CellTag;
00176       typedef typename viennagrid::result_of::point<Config>::type                            PointType;
00177       typedef typename viennagrid::result_of::ncell<Config, 0>::type                         VertexType;
00178       typedef typename viennagrid::result_of::ncell<Config, 1>::type                         EdgeType;
00179       typedef typename viennagrid::result_of::ncell<Config, CellTag::dim>::type   CellType;
00180       
00181       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim>::type    CellRange;
00182       typedef typename viennagrid::result_of::iterator<CellRange>::type                                       CellIterator;
00183 
00184       typedef typename viennagrid::result_of::const_ncell_range<DomainType, 1>::type                          EdgeRange;
00185       typedef typename viennagrid::result_of::iterator<EdgeRange>::type                                       EdgeIterator;
00186       
00187       typedef typename viennagrid::result_of::const_ncell_range<CellType, 1>::type                            EdgeOnCellRange;
00188       typedef typename viennagrid::result_of::iterator<EdgeOnCellRange>::type                                 EdgeOnCellIterator;
00189       
00190       typedef typename viennagrid::result_of::const_ncell_range<EdgeType, 0>::type                            VertexOnEdgeRange;
00191       typedef typename viennagrid::result_of::iterator<VertexOnEdgeRange>::type                               VertexOnEdgeIterator;
00192 
00193       //
00194       // Phase 1: Compute circumcenter of cells and store them on each of the edges (avoids coboundary operations)
00195       //
00196       
00197       CellRange cells = viennagrid::ncells<CellTag::dim>(domain);
00198       for (CellIterator cit  = cells.begin();
00199                         cit != cells.end();
00200                       ++cit)
00201       {
00202         PointType circ_center = circumcenter(*cit);
00203         
00204         // store circumcenter on edges
00205         EdgeOnCellRange edges_on_cell = viennagrid::ncells<1>(*cit);
00206         for (EdgeOnCellIterator eocit  = edges_on_cell.begin();
00207                                 eocit != edges_on_cell.end();
00208                               ++eocit)
00209         {
00210           viennadata::access<InterfaceAreaKey, std::vector<PointType> >(interface_key)(*eocit).push_back(circ_center);
00211         } //for edges on cells
00212       } //for cells
00213 
00214 
00215       //
00216       // Phase 2: Compute edge lengths, interface areas and box volumes from the information stored on edges:
00217       //
00218       EdgeRange edges = viennagrid::ncells<1>(domain);
00219       for (EdgeIterator eit  = edges.begin();
00220                         eit != edges.end();
00221                       ++eit)
00222       {
00223         double edge_length = volume(*eit);
00224         
00225         {
00226           std::vector<PointType> & circ_centers = viennadata::access<InterfaceAreaKey, std::vector<PointType> >(interface_key)(*eit);
00227           
00228           double interface_area = -42.0;
00229           if (circ_centers.size() == 1)  //edge is located on the domain boundary
00230             interface_area = spanned_volume(circ_centers[0], circumcenter(*eit));
00231           else if (circ_centers.size() == 2)
00232             interface_area = spanned_volume(circ_centers[0], circ_centers[1]);
00233           else
00234             throw "More than two circum centers for an edge in two dimensions!"; 
00235           
00236           viennadata::access<InterfaceAreaKey, double>(interface_key)(*eit) = interface_area;
00237           double volume_contribution = interface_area * edge_length / 4.0;
00238           viennadata::access<BoxVolumeKey, double>(box_volume_key)(*eit) = 2.0 * volume_contribution; //volume contribution of both box volumes associated with the edge
00239           viennadata::access<BoxVolumeKey, double>(box_volume_key)(viennagrid::ncells<0>(*eit)[0]) += volume_contribution;
00240           viennadata::access<BoxVolumeKey, double>(box_volume_key)(viennagrid::ncells<0>(*eit)[1]) += volume_contribution;
00241         }
00242         
00243         //delete circ_centers:
00244         viennadata::erase<InterfaceAreaKey, std::vector<PointType> >()(*eit);
00245       }
00246 
00247     }
00248 
00249 
00250     //
00251     // Voronoi information in three dimensions
00252     //
00253 
00255     template <typename DomainType,
00256               typename InterfaceAreaKey,
00257               typename BoxVolumeKey>
00258     void write_voronoi_info(DomainType const & domain,
00259                             InterfaceAreaKey const & interface_key,
00260                             BoxVolumeKey const & box_volume_key,
00261                             viennagrid::tetrahedron_tag)
00262     {
00263       typedef typename DomainType::config_type           Config;
00264       typedef typename Config::cell_tag                  CellTag;
00265       typedef typename viennagrid::result_of::point<Config>::type                            PointType;
00266       typedef typename viennagrid::result_of::ncell<Config, 0>::type                         VertexType;
00267       typedef typename viennagrid::result_of::ncell<Config, 1>::type                         EdgeType;
00268       typedef typename viennagrid::result_of::ncell<Config, 2>::type                         FacetType;
00269       typedef typename viennagrid::result_of::ncell<Config, CellTag::dim>::type   CellType;
00270       
00271       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim>::type    CellRange;
00272       typedef typename viennagrid::result_of::iterator<CellRange>::type                                       CellIterator;
00273 
00274       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim-1>::type  FacetRange;
00275       typedef typename viennagrid::result_of::iterator<FacetRange>::type                                      FacetIterator;
00276       
00277       typedef typename viennagrid::result_of::const_ncell_range<DomainType, 1>::type                          EdgeRange;
00278       typedef typename viennagrid::result_of::iterator<EdgeRange>::type                                       EdgeIterator;
00279       
00280       typedef typename viennagrid::result_of::const_ncell_range<CellType, CellTag::dim-1>::type    FacetOnCellRange;
00281       typedef typename viennagrid::result_of::iterator<FacetOnCellRange>::type                                FacetOnCellIterator;
00282 
00283       typedef typename viennagrid::result_of::const_ncell_range<FacetType, 1>::type                           EdgeOnFacetRange;
00284       typedef typename viennagrid::result_of::iterator<EdgeOnFacetRange>::type                                EdgeOnFacetIterator;
00285       
00286       typedef typename viennagrid::result_of::const_ncell_range<EdgeType, 0>::type                            VertexOnEdgeRange;
00287       typedef typename viennagrid::result_of::iterator<VertexOnEdgeRange>::type                               VertexOnEdgeIterator;
00288 
00289       //
00290       // Step one: Write circumcenters to facets
00291       //
00292       
00293       CellRange cells = viennagrid::ncells<CellTag::dim>(domain);
00294       for (CellIterator cit  = cells.begin();
00295                         cit != cells.end();
00296                       ++cit)
00297       {
00298         PointType circ_center = circumcenter(*cit);
00299         
00300         FacetOnCellRange facets_on_cell = viennagrid::ncells<CellTag::dim-1>(*cit);
00301         for (FacetOnCellIterator focit  = facets_on_cell.begin();
00302                                 focit != facets_on_cell.end();
00303                               ++focit)
00304         {
00305           viennadata::access<InterfaceAreaKey, std::vector<PointType> >(interface_key)(*focit).push_back(circ_center);
00306         } //for edges on cells
00307       } //for cells
00308 
00309 
00310       //
00311       // Step two: Write lines connecting circumcenters to edges
00312       //
00313       FacetRange facets = viennagrid::ncells<CellTag::dim-1>(domain);
00314       for (FacetIterator fit  = facets.begin();
00315                          fit != facets.end();
00316                        ++fit)
00317       {
00318         std::vector<PointType> & circ_centers = viennadata::access<InterfaceAreaKey, std::vector<PointType> >(interface_key)(*fit);
00319         
00320         EdgeOnFacetRange edges_on_facet = viennagrid::ncells<1>(*fit);
00321         for (EdgeOnFacetIterator eofit  = edges_on_facet.begin();
00322                                 eofit != edges_on_facet.end();
00323                               ++eofit)
00324         {
00325           viennadata::access<InterfaceAreaKey, std::vector<PointType> >(interface_key)(*eofit).push_back(circ_centers[0]);
00326           
00327           if (circ_centers.size() == 1)
00328           {
00329             viennadata::access<InterfaceAreaKey,
00330                                std::vector<std::pair<PointType, PointType> > 
00331                               >(interface_key)(*eofit).push_back(std::make_pair(circ_centers[0], circumcenter(*fit)));
00332             viennadata::access<InterfaceAreaKey,
00333                                std::vector<std::pair<PointType, PointType> >
00334                               >(interface_key)(*eofit).push_back(std::make_pair(circumcenter(*eofit), circumcenter(*fit)));
00335           }
00336           else if (circ_centers.size() == 2)
00337             viennadata::access<InterfaceAreaKey,
00338                                std::vector<std::pair<PointType, PointType> >
00339                               >(interface_key)(*eofit).push_back(std::make_pair(circ_centers[0], circ_centers[1]));
00340           else
00341           {
00342             std::cerr << "circ_centers.size() = " << circ_centers.size() << std::endl;
00343             std::cerr << "*fit: " << *fit << std::endl;
00344             throw "More than two circumcenters for a facet in three dimensions!"; 
00345           }
00346             
00347         } //for edges on cells
00348       }
00349 
00350 
00351       //
00352       // Step three: Compute Voronoi information:
00353       //
00354       EdgeRange edges = viennagrid::ncells<1>(domain);
00355       for (EdgeIterator eit  = edges.begin();
00356                         eit != edges.end();
00357                       ++eit)
00358       {
00359         //get vertices of edge:
00360         VertexOnEdgeRange vertices_on_edge = viennagrid::ncells<0>(*eit);
00361         VertexOnEdgeIterator voeit = vertices_on_edge.begin();
00362         
00363         VertexType const & v0 = *voeit;
00364         ++voeit;
00365         VertexType const & v1 = *voeit;
00366         
00367         PointType const & p0 = v0.point();
00368         PointType const & p1 = v1.point();
00369 
00370         //write edge length
00371         double edge_length = spanned_volume(p0, p1);
00372         
00373         std::vector<std::pair<PointType, PointType> > & interface_segments =
00374           viennadata::access<InterfaceAreaKey, std::vector<std::pair<PointType, PointType> > >(interface_key)(*eit);
00375 
00376         //determine inner point of convex interface polygon:
00377         PointType inner_point = (interface_segments[0].first + interface_segments[0].second) / 2.0;
00378         for (size_t i=1; i<interface_segments.size(); ++i)
00379         {
00380           inner_point += (interface_segments[i].first + interface_segments[i].second) / 2.0;
00381         }
00382         inner_point /= interface_segments.size();
00383         //std::cout << "Inner point: " << inner_point << std::endl;
00384         
00385         //compute interface area
00386         double interface_area = 0.0;
00387         for (size_t i=0; i<interface_segments.size(); ++i)
00388         {
00389           //std::cout << "Interface line segment: " << interface_segments[i].first << " to " << interface_segments[i].second << std::endl;
00390           interface_area += spanned_volume(interface_segments[i].first, interface_segments[i].second, inner_point);
00391         }
00392         //std::cout << "Interface area: " << interface_area << std::endl;
00393         viennadata::access<InterfaceAreaKey, double>(interface_key)(*eit) = interface_area;
00394         double volume_contribution = interface_area * edge_length / 6.0;
00395         //std::cout << "Volume contribution: " << volume_contribution << std::endl;
00396         viennadata::access<BoxVolumeKey, double>(box_volume_key)(*eit) = 2.0 * volume_contribution; //volume contribution of both box volumes associated with the edge
00397         viennadata::access<BoxVolumeKey, double>(box_volume_key)(v0) += volume_contribution;
00398         viennadata::access<BoxVolumeKey, double>(box_volume_key)(v1) += volume_contribution;
00399       }
00400 
00401     } //write_voronoi_info(tetrahedron_tag)
00402 
00403 
00404 
00406     template <typename DomainType,
00407               typename InterfaceAreaKey,
00408               typename BoxVolumeKey>
00409     void write_voronoi_info(DomainType const & domain,
00410                             InterfaceAreaKey const & interface_key,
00411                             BoxVolumeKey const & box_volume_key,
00412                             viennagrid::hexahedron_tag)
00413     {
00414       //std::cout << "Warning: Voronoi info for hexahedron is only correct when having regular cuboids only." << std::endl;
00415       typedef typename DomainType::config_type           Config;
00416       typedef typename Config::cell_tag                  CellTag;
00417       typedef typename viennagrid::result_of::point<Config>::type                            PointType;
00418       typedef typename viennagrid::result_of::ncell<Config, 0>::type                         VertexType;
00419       typedef typename viennagrid::result_of::ncell<Config, 1>::type                         EdgeType;
00420       typedef typename viennagrid::result_of::ncell<Config, 2>::type                         FacetType;
00421       typedef typename viennagrid::result_of::ncell<Config, CellTag::dim>::type   CellType;
00422       
00423       typedef typename viennagrid::result_of::const_ncell_range<DomainType, CellTag::dim>::type    CellRange;
00424       typedef typename viennagrid::result_of::iterator<CellRange>::type                                       CellIterator;
00425 
00426       typedef typename viennagrid::result_of::const_ncell_range<CellType, 2>::type                            FacetOnCellRange;
00427       typedef typename viennagrid::result_of::iterator<FacetOnCellRange>::type                                FacetOnCellIterator;
00428 
00429       typedef typename viennagrid::result_of::const_ncell_range<FacetType, 1>::type                           EdgeOnFacetRange;
00430       typedef typename viennagrid::result_of::iterator<EdgeOnFacetRange>::type                                EdgeOnFacetIterator;
00431       
00432       typedef typename viennagrid::result_of::const_ncell_range<EdgeType, 0>::type                            VertexOnEdgeRange;
00433       typedef typename viennagrid::result_of::iterator<VertexOnEdgeRange>::type                               VertexOnEdgeIterator;
00434 
00435       //
00436       // Algorithm: Iterate over all cells, compute circumcenter and add interface area to edge, box volume to vertex.
00437       //
00438       
00439       CellRange cells = viennagrid::ncells<CellTag::dim>(domain);
00440       for (CellIterator cit  = cells.begin();
00441                         cit != cells.end();
00442                       ++cit)
00443       {
00444         PointType cell_center = circumcenter(*cit);
00445 
00446         FacetOnCellRange facets_on_cell = viennagrid::ncells<2>(*cit);
00447         for (FacetOnCellIterator focit  = facets_on_cell.begin();
00448                                 focit != facets_on_cell.end();
00449                               ++focit)
00450         {
00451           PointType facet_center = circumcenter(*focit);
00452           
00453           //iterate over edges:
00454           EdgeOnFacetRange edges_on_facet = viennagrid::ncells<1>(*focit);
00455           for (EdgeOnFacetIterator eocit  = edges_on_facet.begin();
00456                                   eocit != edges_on_facet.end();
00457                                 ++eocit)
00458           {
00459             PointType edge_midpoint = viennagrid::circumcenter(*eocit);
00460             
00461             // interface contribution:
00462             viennadata::access<InterfaceAreaKey, double>(interface_key)(*eocit) += spanned_volume(cell_center, facet_center, edge_midpoint);
00463             
00464             //box volume contribution:
00465             double edge_contribution = 0;
00466             VertexOnEdgeRange vertices_on_edge = viennagrid::ncells<0>(*eocit);
00467             for (VertexOnEdgeIterator voeit  = vertices_on_edge.begin();
00468                                       voeit != vertices_on_edge.end();
00469                                     ++voeit)
00470             {
00471               double contribution = spanned_volume(cell_center, facet_center, edge_midpoint, voeit->point());
00472               edge_contribution += contribution;
00473               viennadata::access<BoxVolumeKey, double>(box_volume_key)(*voeit) += contribution;
00474             }
00475             viennadata::access<BoxVolumeKey, double>(box_volume_key)(*eocit) += edge_contribution;
00476           } //for edges on facet
00477 
00478         } //for facets on cell
00479 
00480       } //for cells
00481       
00482     }
00483 
00484   } //namespace detail
00485 
00486   //
00487   // The public interface
00488   //
00489   
00496   template <typename DomainType,
00497             typename InterfaceAreaKey,
00498             typename BoxVolumeKey>
00499   void apply_voronoi(DomainType const & domain,
00500                      InterfaceAreaKey const & interface_area_key = viennagrid::voronoi_interface_area_key(),
00501                      BoxVolumeKey const & box_volume_key = viennagrid::voronoi_box_volume_key())
00502   {
00503     detail::write_voronoi_info(domain,
00504                                interface_area_key,
00505                                box_volume_key,
00506                                typename DomainType::config_type::cell_tag());
00507   }
00508     
00510   template <typename DomainType>
00511   void apply_voronoi(DomainType const & domain)
00512   {
00513     apply_voronoi(domain, viennagrid::voronoi_interface_area_key(), viennagrid::voronoi_box_volume_key());
00514   }
00515     
00516 } //namespace viennagrid
00517 #endif

Generated on Wed Sep 14 2011 19:21:30 for ViennaGrid by  doxygen 1.7.1