dune-geometry  2.2.0
generalvertexorder.hh
Go to the documentation of this file.
00001 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
00002 // vi: set et ts=8 sw=2 sts=2:
00003 
00004 #ifndef DUNE_GEOMETRY_GENERALVERTEXORDER_HH
00005 #define DUNE_GEOMETRY_GENERALVERTEXORDER_HH
00006 
00007 #include <algorithm>
00008 #include <cassert>
00009 #include <cstddef>
00010 #include <functional>
00011 #include <iterator>
00012 #include <vector>
00013 
00014 #include <dune/common/iteratorfacades.hh>
00015 
00016 #include "type.hh"
00017 #include "referenceelements.hh"
00018 
00019 namespace Dune {
00020 
00022 
00037   template<class InIterator, class OutIterator>
00038   void reduceOrder(const InIterator& inBegin, const InIterator& inEnd,
00039                    OutIterator outIt)
00040   {
00041     typedef std::less<
00042       typename std::iterator_traits<InIterator>::value_type
00043       > less_t;
00044     static const less_t less = less_t();
00045 
00046     for(InIterator inIt = inBegin; inIt != inEnd; ++inIt, ++outIt)
00047       *outIt = std::count_if(inBegin, inEnd, std::bind2nd(less, *inIt));
00048   }
00049 
00051 
00066   template<std::size_t dim, class Index_ = std::size_t>
00067   class GeneralVertexOrder {
00068     typedef GenericReferenceElement<double, dim> RefElem;
00069     typedef GenericReferenceElements<double, dim> RefElems;
00070 
00071     const RefElem& refelem;
00072     GeometryType gt;
00073     std::vector<Index_> vertexOrder;
00074 
00075   public:
00077     typedef Index_ Index;
00078 
00080     class iterator;
00081 
00083     static const std::size_t dimension = dim;
00085     const GeometryType &type() const { return gt; }
00086 
00088 
00096     template<class InIterator>
00097     GeneralVertexOrder(const GeometryType& gt_, const InIterator &inBegin,
00098                        const InIterator &inEnd) :
00099       refelem(RefElems::general(gt_)), gt(gt_),
00100       vertexOrder(refelem.size(dim))
00101     { reduceOrder(inBegin, inEnd, vertexOrder.begin()); }
00102 
00104 
00108     iterator begin(std::size_t codim, std::size_t subEntity) const
00109     { return iterator(*this, codim, subEntity); }
00111 
00115     iterator end(std::size_t codim, std::size_t subEntity) const {
00116       return iterator(*this, codim, subEntity,
00117                       refelem.size(subEntity, codim, dim));
00118     }
00119 
00121 
00128     void getReduced(std::size_t codim, std::size_t subEntity,
00129                     std::vector<Index>& order) const
00130     {
00131       order.resize(refelem.size(subEntity, codim, dim));
00132       reduceOrder(begin(codim, subEntity), end(codim, subEntity),
00133                   order.begin());
00134     }
00135   };
00136 
00138 
00141   template<std::size_t dim, class Index_>
00142   class GeneralVertexOrder<dim, Index_>::iterator :
00143     public Dune::RandomAccessIteratorFacade<iterator, const Index_>
00144   {
00145     const GeneralVertexOrder *order;
00146     std::size_t codim;
00147     std::size_t subEntity;
00148     std::size_t vertex;
00149 
00150     iterator(const GeneralVertexOrder &order_, std::size_t codim_,
00151              std::size_t subEntity_, std::size_t vertex_ = 0) :
00152       order(&order_), codim(codim_), subEntity(subEntity_), vertex(vertex_)
00153     { }
00154 
00155   public:
00156     const Index &dereference() const {
00157       return order->vertexOrder[order->refelem.subEntity(subEntity, codim,
00158                                                          vertex, dim)];
00159     }
00160     const Index &elementAt(std::ptrdiff_t n) const {
00161       return order->vertexOrder[order->refelem.subEntity(subEntity, codim,
00162                                                          vertex+n, dim)];
00163     }
00164     bool equals(const iterator &other) const {
00165       return order == other.order && codim == other.codim &&
00166         subEntity == other.subEntity && vertex == other.vertex;
00167     }
00168     void increment() { ++vertex; }
00169     void decrement() { --vertex; }
00170     void advance(std::ptrdiff_t n) { vertex += n; }
00171     std::ptrdiff_t distanceTo(const iterator &other) const {
00172       // make sure we reference the same container
00173       assert(order == other.order && codim == other.codim &&
00174              subEntity == other.subEntity);
00175       if(vertex < other.vertex)         return other.vertex - vertex;
00176       else return -static_cast<std::ptrdiff_t>(vertex - other.vertex);
00177     }
00178 
00179     friend class GeneralVertexOrder<dim, Index>;
00180 
00182 
00187     iterator() { }
00188   };
00189 } // namespace Dune
00190 
00191 #endif // DUNE_GEOMETRY_GENERALVERTEXORDER_HH