AFEPack
|
00001 00011 #ifndef _Geometry_h_ 00012 #define _Geometry_h_ 00013 00014 #ifndef NULL 00015 #define NULL 0 00016 #endif // NULL 00017 00018 #include <stdarg.h> 00019 #include <dlfcn.h> 00020 00021 #include <cmath> 00022 #include <iostream> 00023 #include <fstream> 00024 #include <string> 00025 #include <vector> 00026 #include <list> 00027 #include <iterator> 00028 #include <algorithm> 00029 00030 #include <base/exceptions.h> 00031 00032 #include "Miscellaneous.h" 00033 00034 namespace afepack { 00035 template <int DIM> class Point; 00036 } 00037 class Geometry; 00038 class GeometryBM; 00039 template <int DIM, int DOW> class Mesh; 00040 template <int DIM> class QuadratureInfo; 00041 template <int DIM> class QuadratureInfoAdmin; 00042 template <int DIM> class TemplateGeometry; 00043 template <int DIM, int DOW> class SimplestMesh; 00044 00045 template <int DIM> afepack::Point<DIM> midpoint(const afepack::Point<DIM>&, const afepack::Point<DIM>&); 00046 template <int DIM> double distance(const afepack::Point<DIM>&, const afepack::Point<DIM>&); 00047 template <int DIM> afepack::Point<DIM> barycenter(const std::vector<afepack::Point<DIM> >&, const double * = NULL); 00048 template <int DIM> afepack::Point<DIM> operator+(const afepack::Point<DIM>&, const afepack::Point<DIM>&); 00049 template <int DIM> afepack::Point<DIM> operator-(const afepack::Point<DIM>&, const afepack::Point<DIM>&); 00050 template <int DIM> std::istream& operator>>(std::istream&, afepack::Point<DIM>&); 00051 template <int DIM> std::ostream& operator<<(std::ostream&, const afepack::Point<DIM>&); 00052 00053 std::istream& operator>>(std::istream&, Geometry&); 00054 std::ostream& operator<<(std::ostream&, const Geometry&); 00055 bool is_same(const Geometry&, const Geometry&); 00056 00057 template <int DIM> filtering_istream& operator>>(filtering_istream&, QuadratureInfo<DIM>&); 00058 template <int DIM> std::ostream& operator<<(std::ostream&, const QuadratureInfo<DIM>&); 00059 00060 template <int DIM> filtering_istream& operator>>(filtering_istream&, QuadratureInfoAdmin<DIM>&); 00061 template <int DIM> std::ostream& operator<<(std::ostream&, const QuadratureInfoAdmin<DIM>&); 00062 00063 typedef std::vector<Geometry>::iterator GeometryIterator; 00064 00065 template <int DIM, int DOW> std::istream& operator>>(std::istream&, Mesh<DIM,DOW>&); 00066 template <int DIM, int DOW> std::ostream& operator<<(std::ostream&, const Mesh<DIM,DOW>&); 00067 00068 template <int DIM> filtering_istream& operator>>(filtering_istream&, TemplateGeometry<DIM>&); 00069 template <int DIM> std::ostream& operator<<(std::ostream&, const TemplateGeometry<DIM>&); 00070 00072 00077 namespace afepack { 00078 template <int DIM> 00079 class Point 00080 { 00081 public: 00082 enum { dim = DIM }; 00083 private: 00084 double x[DIM]; 00085 public: 00086 Point(); 00087 Point(const double *); 00088 Point(const Point&); 00089 Point(double, ...); 00090 ~Point(); 00091 public: 00092 Point<DIM>& operator=(const Point<DIM>&); 00093 operator const double *() const; 00094 operator double *(); 00095 const double& operator[](int) const; 00096 double& operator[](int); 00097 double length() const; 00098 Point<DIM>& operator+=(const Point<DIM>&); 00099 Point<DIM>& operator-=(const Point<DIM>&); 00100 Point<DIM>& operator*=(const double&); 00101 Point<DIM>& operator/=(const double&); 00102 public: 00103 friend Point<DIM> midpoint<>(const Point<DIM>&, const Point<DIM>&); 00104 friend double distance<>(const Point<DIM>&, const Point<DIM>&); 00105 friend Point<DIM> barycenter<>(const std::vector<Point<DIM> >&, const double *); 00106 friend Point<DIM> operator+ <>(const Point<DIM>&, const Point<DIM>&); 00107 friend Point<DIM> operator- <>(const Point<DIM>&, const Point<DIM>&); 00108 friend std::istream& operator>> <>(std::istream&, Point<DIM>&); 00109 friend std::ostream& operator<< <>(std::ostream&, const Point<DIM>&); 00110 }; 00111 } 00112 00118 class Geometry 00119 { 00120 private: 00121 int ind; 00122 std::vector<int> vtx; 00123 std::vector<int> bnd; 00124 public: 00125 Geometry(); 00126 Geometry(const Geometry&); 00127 ~Geometry(); 00128 public: 00129 Geometry& operator=(const Geometry&); 00130 int index() const; 00131 int& index(); 00132 int n_vertex() const; 00133 int n_boundary() const; 00134 const std::vector<int>& vertex() const; 00135 std::vector<int>& vertex(); 00136 const std::vector<int>& boundary() const; 00137 std::vector<int>& boundary(); 00138 int vertex(int) const; 00139 int& vertex(int); 00140 int boundary(int) const; 00141 int& boundary(int); 00142 public: 00143 friend bool isSame(const Geometry&, const Geometry&); 00144 template <int, int> friend class Mesh; 00145 friend std::istream& operator>>(std::istream&, Geometry&); 00146 friend std::ostream& operator<<(std::ostream&, const Geometry&); 00147 }; 00148 00154 class GeometryBM : public Geometry 00155 { 00156 public: 00157 typedef int bmark_t; 00158 private: 00159 bmark_t bm; 00160 public: 00161 GeometryBM(); 00162 GeometryBM(const GeometryBM&); 00163 ~GeometryBM(); 00164 public: 00165 GeometryBM& operator=(const GeometryBM&); 00166 int boundaryMark() const; 00167 int& boundaryMark(); 00168 public: 00169 template <int, int> friend class Mesh; 00170 friend std::istream& operator>>(std::istream&, GeometryBM&); 00171 friend std::ostream& operator<<(std::ostream&, const GeometryBM&); 00172 }; 00173 00179 template <int DIM, int DOW=DIM> 00180 class Mesh 00181 { 00182 public: 00183 enum { dim = DIM, dow = DOW }; 00184 typedef afepack::Point<DOW> point_t; 00185 public: 00186 typedef GeometryBM::bmark_t bmark_t; 00187 private: 00188 typedef Mesh<DIM,DOW> mesh_t; 00189 00190 std::vector<point_t > pnt; 00191 std::vector<std::vector<GeometryBM> > geo; 00194 public: 00195 Mesh(); 00196 Mesh(const mesh_t&); 00197 virtual ~Mesh(); 00198 public: 00199 mesh_t& operator=(const mesh_t&); 00200 unsigned int n_point() const; 00201 unsigned int n_geometry(int) const; 00202 const std::vector<point_t >& point() const; 00203 std::vector<point_t >& point(); 00204 const point_t& point(int) const; 00205 point_t& point(int); 00206 const std::vector<std::vector<GeometryBM> >& geometry() const; 00207 std::vector<std::vector<GeometryBM> >& geometry(); 00208 const std::vector<GeometryBM>& geometry(int) const; 00209 std::vector<GeometryBM>& geometry(int); 00210 const GeometryBM& geometry(int, int) const; 00211 GeometryBM& geometry(int, int); 00212 bmark_t boundaryMark(int, int) const; 00213 bmark_t& boundaryMark(int, int); 00214 void renumerateElement(); 00216 void renumerateElementHSFC(void (*f)(double,double,double,double&,double&,double&)=NULL); 00217 void readData(const std::string&); 00218 void writeData(const std::string&) const; 00219 void readData1d(const std::string&); 00220 virtual void writeEasyMesh(const std::string&) const {}; 00221 virtual void writeTecplotData(const std::string&) const {}; 00222 00223 DeclException1(ExcMeshData, 00224 char *, 00225 << "Mesh data error: " << arg1); 00226 public: 00227 friend std::istream& operator>> <>(std::istream&, mesh_t&); 00228 friend std::ostream& operator<< <>(std::ostream&, const mesh_t&); 00229 friend filtering_istream& operator>> <>(filtering_istream&, TemplateGeometry<DIM>&); 00230 friend std::ostream& operator<< <>(std::ostream&, const TemplateGeometry<DIM>&); 00231 }; 00232 00237 template <int DIM> 00238 class QuadratureInfo 00239 { 00240 public: 00241 enum { dim = DIM }; 00242 private: 00243 int alg_acc; 00244 std::vector<afepack::Point<DIM> > pnt; 00245 std::vector<double> wei; 00246 public: 00247 QuadratureInfo(); 00248 QuadratureInfo(const QuadratureInfo<DIM>&); 00249 ~QuadratureInfo(); 00250 public: 00251 QuadratureInfo<DIM>& operator=(const QuadratureInfo<DIM>&); 00252 int n_quadraturePoint() const; 00253 int algebricAccuracy() const; 00254 int& algebricAccuracy(); 00255 const std::vector<afepack::Point<DIM> >& quadraturePoint() const; 00256 std::vector<afepack::Point<DIM> >& quadraturePoint(); 00257 const afepack::Point<DIM>& quadraturePoint(int) const; 00258 afepack::Point<DIM>& quadraturePoint(int); 00259 const std::vector<double>& weight() const; 00260 std::vector<double>& weight(); 00261 const double& weight(int) const; 00262 double& weight(int); 00263 public: 00264 friend filtering_istream& operator>> <>(filtering_istream&, QuadratureInfo<DIM>&); 00265 friend std::ostream& operator<< <>(std::ostream&, const QuadratureInfo<DIM>&); 00266 }; 00267 00272 template <int DIM> 00273 class QuadratureInfoAdmin : public std::vector<QuadratureInfo<DIM> > 00274 { 00275 public: 00276 enum { dim = DIM }; 00277 private: 00278 std::vector<int> acc_tbl; 00279 public: 00280 QuadratureInfoAdmin(); 00281 QuadratureInfoAdmin(const QuadratureInfoAdmin<DIM>&); 00282 ~QuadratureInfoAdmin(); 00283 public: 00284 QuadratureInfoAdmin<DIM>& operator=(const QuadratureInfoAdmin<DIM>&); 00285 const QuadratureInfo<DIM>& find(int) const; 00286 QuadratureInfo<DIM>& find(int); 00287 public: 00288 friend filtering_istream& operator>> <>(filtering_istream&, QuadratureInfoAdmin<DIM>&); 00289 friend std::ostream& operator<< <>(std::ostream&, const QuadratureInfoAdmin<DIM>&); 00290 }; 00291 00299 template <int DIM> 00300 class TemplateGeometry : public Mesh<DIM,DIM> 00301 { 00302 public: 00303 enum { dim = DIM }; 00304 private: 00305 std::string library_path; 00306 void * handle; 00307 std::string library_name; 00308 std::string volume_function_name; 00309 double (*volume_function)(const double **); 00310 QuadratureInfoAdmin<DIM> quad_info; 00311 public: 00312 TemplateGeometry(); 00313 TemplateGeometry(const TemplateGeometry<DIM>&); 00314 ~TemplateGeometry(); 00315 public: 00316 TemplateGeometry& operator=(const TemplateGeometry<DIM>&); 00317 void loadFunction(); 00318 void unloadFunction(); 00319 const std::vector<afepack::Point<DIM> >& vertexArray() const; 00320 const QuadratureInfoAdmin<DIM>& quadratureInfo() const; 00321 QuadratureInfoAdmin<DIM>& quadratureInfo(); 00322 const QuadratureInfo<DIM>& findQuadratureInfo(int) const; 00323 double volume() const; 00324 void readData(const std::string&); 00325 void writeData(const std::string&) const; 00327 DeclException1(ExcTemplateGeometryData, 00328 char *, 00329 << "Template geometry data error: " << arg1); 00330 DeclException1(ExcFileOpen, char *, 00331 << "Can't open library " 00332 << arg1); 00333 DeclException2(ExcLoadFunction, char *, char *, 00334 << "Can't load function " 00335 << arg1 00336 << " from library " 00337 << arg2); 00338 00339 public: 00340 friend filtering_istream& operator>> <>(filtering_istream&, TemplateGeometry<DIM>&); 00341 friend std::ostream& operator<< <>(std::ostream&, const TemplateGeometry<DIM>&); 00342 }; 00343 00353 template <int DIM, int DOW=DIM> 00354 class SimplestMesh 00355 { 00356 public: 00357 enum { dim = DIM, dow = DOW }; 00358 protected: 00359 struct SimplestGeometry { 00360 public: 00361 int template_element; 00362 std::vector<int> vertex; 00363 }; 00364 std::vector<afepack::Point<DOW> > pnt; 00365 std::vector<SimplestGeometry> ele; 00366 std::vector<TemplateGeometry<DIM> > * tg; 00367 public: 00368 SimplestMesh() {}; 00369 virtual ~SimplestMesh() {}; 00370 public: 00371 int n_point() const {return pnt.size();} 00372 int n_element() const {return ele.size();} 00373 const std::vector<afepack::Point<DOW> >& point() const {return pnt;} 00374 std::vector<afepack::Point<DOW> >& point() {return pnt;} 00375 const afepack::Point<DOW>& point(int i) const {return pnt[i];} 00376 afepack::Point<DOW>& point(int i) {return pnt[i];} 00377 const std::vector<SimplestGeometry>& element() const {return ele;} 00378 std::vector<SimplestGeometry>& element() {return ele;} 00379 const SimplestGeometry& element(int i) const {return ele[i];} 00380 SimplestGeometry& element(int i) {return ele[i];} 00381 const std::vector<int>& elementVertex(int i) const {return ele[i].vertex;} 00382 std::vector<int>& elementVertex(int i) {return ele[i].vertex;} 00383 int elementVertex(int i, int j) const {return ele[i].vertex[j];} 00384 int& elementVertex(int i, int j) {return ele[i].vertex[j];} 00385 const std::vector<TemplateGeometry<DIM> >& templateGeometry() const {return *tg;} 00386 std::vector<TemplateGeometry<DIM> >& templateGeometry() {return *tg;} 00387 const TemplateGeometry<DIM>& templateGeometry(int i) const {return (*tg)[i];} 00388 TemplateGeometry<DIM>& templateGeometry(int i) {return (*tg)[i];} 00389 void setTemplateGeometry(std::vector<TemplateGeometry<DIM> >& t) {tg = &t;} 00390 public: 00391 virtual void generateMesh(Mesh<DIM,DOW>&); 00392 }; 00393 00394 #endif //_Geometry_h_ 00395