Go to the documentation of this file.00001 #ifndef VIENNAGRID_ALGORITHM_CENTROID_HPP
00002 #define VIENNAGRID_ALGORITHM_CENTROID_HPP
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include <iostream>
00022 #include <sstream>
00023 #include <string>
00024 #include <stdexcept>
00025
00026 #include "viennagrid/forwards.h"
00027 #include "viennagrid/topology/all.hpp"
00028 #include "viennagrid/algorithm/volume.hpp"
00029
00035 namespace viennagrid
00036 {
00037 namespace detail
00038 {
00039
00040
00041
00043 template <typename ElementType>
00044 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00045 centroid(ElementType const & cell, viennagrid::triangle_tag)
00046 {
00047 typedef typename ElementType::config_type Config;
00048 typedef typename ElementType::tag ElementTag;
00049 typedef typename viennagrid::result_of::point<Config>::type PointType;
00050 typedef typename viennagrid::result_of::ncell<Config, 0>::type VertexType;
00051 typedef typename viennagrid::result_of::ncell<Config, 1>::type EdgeType;
00052
00053 typedef typename viennagrid::result_of::const_ncell_range<ElementType, 0>::type VertexOnCellRange;
00054 typedef typename viennagrid::result_of::iterator<VertexOnCellRange>::type VertexOnCellIterator;
00055
00056 PointType p0(0.0, 0.0);
00057
00058 VertexOnCellRange vertices = viennagrid::ncells<0>(cell);
00059 for (VertexOnCellIterator vocit = vertices.begin();
00060 vocit != vertices.end();
00061 ++vocit)
00062 {
00063 p0 += vocit->point();
00064 }
00065
00066 p0 /= viennagrid::topology::bndcells<ElementTag, 0>::num;
00067
00068 return p0;
00069 }
00070
00071
00073 template <typename ElementType>
00074 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00075 centroid(ElementType const & cell, viennagrid::tetrahedron_tag)
00076 {
00077 return centroid(cell, viennagrid::triangle_tag());
00078 }
00079
00080
00081
00082
00084 template <typename ElementType>
00085 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00086 centroid(ElementType const & cell, viennagrid::quadrilateral_tag)
00087 {
00088 return centroid(cell, viennagrid::triangle_tag());
00089 }
00090
00091
00092
00093
00095 template <typename ElementType>
00096 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00097 centroid(ElementType const & cell, viennagrid::hexahedron_tag)
00098 {
00099 return centroid(cell, viennagrid::triangle_tag());
00100 }
00101
00102
00103
00105 template <typename ElementType>
00106 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00107 centroid(ElementType const & cell, viennagrid::simplex_tag<1>)
00108 {
00109 return centroid(cell, viennagrid::triangle_tag());
00110 }
00111
00113 template <typename ElementType>
00114 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00115 centroid(ElementType const & cell, viennagrid::hypercube_tag<1>)
00116 {
00117 return centroid(cell, viennagrid::triangle_tag());
00118 }
00119
00120
00122 template <typename ElementType>
00123 typename viennagrid::result_of::point<typename ElementType::config_type>::type
00124 centroid(ElementType const & cell, viennagrid::point_tag)
00125 {
00126 return cell.point();
00127 }
00128
00129
00130
00131
00132 template <typename DomainSegmentType>
00133 typename viennagrid::result_of::point<typename DomainSegmentType::config_type>::type
00134 centroid_domseg(DomainSegmentType const & domseg)
00135 {
00136 typedef typename DomainSegmentType::config_type ConfigType;
00137 typedef typename ConfigType::cell_tag CellTag;
00138
00139 typedef typename viennagrid::result_of::point<ConfigType>::type PointType;
00140 typedef typename viennagrid::result_of::ncell<ConfigType, 0>::type VertexType;
00141 typedef typename viennagrid::result_of::ncell<ConfigType, 1>::type EdgeType;
00142
00143 typedef typename viennagrid::result_of::const_ncell_range<DomainSegmentType,
00144 CellTag::dim>::type CellRange;
00145 typedef typename viennagrid::result_of::iterator<CellRange>::type CellIterator;
00146
00147 PointType result = 0;
00148 double volume = 0;
00149
00150 CellRange cells = viennagrid::ncells(domseg);
00151 for (CellIterator cit = cells.begin(); cit != cells.end(); ++cit)
00152 {
00153 double vol_cell = viennagrid::volume(*cit);
00154 result += vol_cell * centroid(*cit);
00155 volume += vol_cell;
00156 }
00157
00158 return result / volume;
00159 }
00160
00161 }
00162
00163
00168 template <typename CellType>
00169 typename viennagrid::result_of::point<typename CellType::config_type>::type
00170 centroid(CellType const & cell)
00171 {
00172 return detail::centroid(cell, typename CellType::tag());
00173 }
00174
00179 template <typename ConfigType>
00180 typename viennagrid::result_of::point<ConfigType>::type
00181 centroid(domain_t<ConfigType> const & domain)
00182 {
00183 return detail::centroid_domseg(domain);
00184 }
00185
00190 template <typename ConfigType>
00191 typename viennagrid::result_of::point<ConfigType>::type
00192 centroid(segment_t<ConfigType> const & segment)
00193 {
00194 return detail::centroid_domseg(segment);
00195 }
00196
00197 }
00198 #endif