AFEPack
DGFEMSpace.h
浏览该文件的文档。
00001 
00002 // DGFEMSpace.h : by R.Lie, July 7, 2003
00003 //
00004 
00005 #ifndef _DGFEMSpace_h_
00006 #define _DGFEMSpace_h_
00007 
00008 #include <lac/full_matrix.h>
00009 
00010 #include "Miscellaneous.h"
00011 #include "Geometry.h"
00012 #include "TemplateElement.h"
00013 #include "FEMSpace.h"
00014 
00015 template <int TDIM, int DIM> class TemplateDGElement;
00016 template <class value_type, int DIM, int DOW, int TDIM, int TDIM1> class DGElement;
00017 template <class value_type, int DIM, int DOW, int TDIM, int TDIM1> class DGFEMSpace;
00018 
00019 template <class value_type, int DIM, int DOW, int TDIM, int TDIM1> std::vector<double> unitOutNormal(const afepack::Point<DIM>&, const Element<value_type,DIM,DOW,TDIM>&, const DGElement<value_type,DIM,DOW,TDIM,TDIM1>&);
00020 template <class value_type, int DIM, int DOW, int TDIM, int TDIM1> std::vector<std::vector<double> > unitOutNormal(const std::vector<afepack::Point<DIM> >&, const Element<value_type,DIM,DOW,TDIM>&, const DGElement<value_type,DIM,DOW,TDIM,TDIM1>&);
00021 
00023 
00024 template <int TDIM, int DIM=TDIM+1>
00025   class TemplateDGElement
00026   {
00027   private:
00028   TemplateGeometry<TDIM> * geo;
00029   CoordTransform<TDIM,DIM> * ct;
00030   public:
00031   TemplateDGElement(TemplateGeometry<TDIM>& _geo = *((TemplateGeometry<TDIM> *)NULL),
00032                     CoordTransform<TDIM,DIM>& _ct = *((CoordTransform<TDIM,DIM> *)NULL)) :
00033   geo(&_geo), ct(&_ct) {};
00034   TemplateDGElement(const TemplateDGElement<TDIM,DIM>& other) :
00035   geo(other.geo), ct(other.ct) {};
00036   ~TemplateDGElement() {};
00037   public:
00038   void reinit(TemplateGeometry<TDIM>& _geo = *((TemplateGeometry<TDIM> *)NULL),
00039               CoordTransform<TDIM,DIM>& _ct = *((CoordTransform<TDIM,DIM> *)NULL)) {geo = &_geo; ct = &_ct;};
00040   const TemplateGeometry<TDIM>& geometry() const {return *geo;};
00041   TemplateGeometry<TDIM>& geometry() {return *geo;};
00042   const CoordTransform<TDIM,DIM>& coordTransform() const {return *ct;};
00043   CoordTransform<TDIM,DIM>& coordTransform() {return *ct;};
00044   const std::vector<afepack::Point<TDIM> >& vertexArray() const {return geo->vertexArray();};
00045   const QuadratureInfoAdmin<TDIM>& quadratureInfo() const {return geo->quadratureInfo();};
00046   QuadratureInfoAdmin<TDIM>& quadratureInfo() {return geo->quadratureInfo();};
00047   const QuadratureInfo<TDIM>& findQuadratureInfo(const int& i) const {return geo->findQuadratureInfo(i);};
00048   double volume() const {return geo->volume();};
00049   };
00050 
00051 template <class value_type, int DIM, int DOW=DIM, int TDIM=DIM, int TDIM1=DIM-1>
00052   class DGElement
00053   {
00054   public:
00055   enum { dim = DIM, dow = DOW, tdim = TDIM, tdim1 = TDIM1 };
00056   typedef value_type value_t;
00057   typedef DGFEMSpace<value_t,DIM,DOW,TDIM,TDIM1> fe_space_t;
00058   typedef Element<value_t,DIM,DOW,TDIM> element_t;
00059   typedef TemplateDGElement<TDIM1,DOW> dg_template_t;
00060   typedef DGElement<value_t,DIM,DOW,TDIM,TDIM1> dg_element_t;
00061   typedef typename Mesh<DIM,DOW>::bmark_t bmark_t;
00062   typedef CoordTransform<TDIM1,DOW> coord_trans_t;
00063 
00064   private:
00065   fe_space_t * sp; 
00066   int geometry_index; 
00067   int template_element_index; 
00068   std::vector<std::vector<int> > geo_img; 
00069   element_t * neigh[2]; 
00070   int bnd_idx[2]; 
00072   public:
00073   DGElement(fe_space_t& = *((fe_space_t *)NULL)); 
00074   DGElement(const dg_element_t&);  
00075   ~DGElement(); 
00076   public:
00077   void reinit(fe_space_t&, const int&, const int&); 
00078   void reinit(const int&, const int&); 
00079   dg_element_t& operator=(const dg_element_t&); 
00080   const fe_space_t& femSpace() const; 
00081   fe_space_t& femSpace(); 
00082   const int& index() const; 
00083   int& index(); 
00084   const int& boundaryIndex(const u_int&) const;
00085   int& boundaryIndex(const u_int&);
00086   GeometryBM& geometry() const; 
00087   bmark_t boundaryMark() const; 
00088   void geometry(const GeometryBM&); 
00089   void geometry(const int&); 
00090   TemplateDGElement<TDIM1,DOW>& templateElement() const; 
00091   void templateElement(const int&); 
00092   TemplateGeometry<TDIM1>& templateGeometry() const; 
00093   const std::vector<std::vector<int> >& geometryImage() const; 
00094   std::vector<std::vector<int> >& geometryImage(); 
00095   void buildGeometryImage(); 
00096   void buildVertexArray(std::vector<afepack::Point<DOW> >&) const; 
00097   const double ** buildVertexArray() const; 
00099   afepack::Point<DOW> local_to_global(const afepack::Point<TDIM1>&) const; 
00100   afepack::Point<TDIM1> global_to_local(const afepack::Point<DOW>&) const; 
00101   double local_to_global_jacobian(const afepack::Point<TDIM1>&) const; 
00102   double global_to_local_jacobian(const afepack::Point<DOW>&) const; 
00103   std::vector<afepack::Point<DOW> > local_to_global(const std::vector<afepack::Point<TDIM1> >&) const; 
00104   std::vector<afepack::Point<TDIM1> > global_to_local(const std::vector<afepack::Point<DOW> >&) const; 
00105   std::vector<double> local_to_global_jacobian(const std::vector<afepack::Point<TDIM1> >&) const; 
00106   std::vector<double> global_to_local_jacobian(const std::vector<afepack::Point<DOW> >&) const; 
00108   const element_t * p_neighbourElement(const int& i) const;
00109   element_t * p_neighbourElement(const int& i);
00110   const element_t& neighbourElement(const int& i) const;
00111   element_t& neighbourElement(const int& i);
00112 
00113   const QuadratureInfo<TDIM1>& findQuadratureInfo(const int&) const; 
00114   std::vector<double> unitNormal() const;
00115 
00116   public:
00117   DeclException1(ExcMeshData, std::string, << "Mesh data uncompatible: " << arg1);
00118   public:
00119   friend std::vector<double> unitOutNormal <>(const afepack::Point<DIM>&, 
00120                                               const element_t&, 
00121                                               const dg_element_t&);
00122   friend std::vector<std::vector<double> > unitOutNormal <>(const std::vector<afepack::Point<DIM> >&, 
00123                                                             const element_t&, 
00124                                                             const dg_element_t&);
00125   friend class DGFEMSpace<value_type, DIM, DOW, TDIM, TDIM1>;
00126 #ifdef __SERIALIZATION__
00127   template <class T> T *
00128   new_property(const property_id_t<T>& pid) const {
00129     return femSpace().new_property(*this, pid);
00130   }
00131   template <class T> T *
00132   get_property(const property_id_t<T>& pid) const {
00133     return femSpace().get_property(*this, pid);
00134   }
00135   template <class T> void
00136   free_property(const property_id_t<T>& pid) const {
00137     femSpace().free_property(*this, pid);
00138   }
00139 #endif // __SERIALIZATION__
00140   };
00141 
00142 template <class value_type, int DIM, int DOW=DIM, int TDIM=DIM, int TDIM1=DIM-1>
00143   class DGFEMSpace : public FEMSpace<value_type, DIM, DOW, TDIM>
00144   {
00145   public:
00146   enum { dim = DIM, dow = DOW, tdim = TDIM, tdim1 = TDIM1 };
00147   typedef value_type value_t;
00148   typedef DGFEMSpace<value_t,DIM,DOW,TDIM,TDIM1> fe_space_t;
00149   typedef DGElement<value_t,DIM,DOW,TDIM,TDIM1> dg_element_t;
00150   typedef FEMSpace<value_t, DIM, DOW, TDIM> base_t;
00151   typedef TemplateDGElement<TDIM1,DOW> dg_template_t;
00152   typedef typename base_t::template_t template_t;
00153   typedef typename base_t::element_t element_t;
00154 
00155   typedef typename std::vector<dg_element_t>::iterator DGElementIterator;
00156   typedef typename std::vector<dg_element_t>::const_iterator ConstDGElementIterator;
00157 
00158   private:
00159   std::vector<dg_template_t> * tmp_dgele;
00160   std::vector<dg_element_t> dgele;
00161   public:
00162   explicit DGFEMSpace(Mesh<DIM,DOW>& = *((Mesh<DIM,DOW> *)NULL),
00163                       std::vector<template_t>& = *((std::vector<template_t> *)NULL),
00164                       std::vector<dg_template_t>& = *((std::vector<dg_template_t> *)NULL));
00165   DGFEMSpace(const fe_space_t&);
00166   virtual ~DGFEMSpace();
00167   public:
00168   void reinit(Mesh<DIM,DOW>&, 
00169               std::vector<template_t>& = *((std::vector<template_t> *)NULL),
00170               std::vector<dg_template_t>& = *((std::vector<dg_template_t> *)NULL));
00171   fe_space_t& operator=(const fe_space_t&);
00172   const std::vector<dg_template_t>& templateDGElement() const {return *tmp_dgele;};
00173   std::vector<dg_template_t>& templateDGElement() {return *tmp_dgele;}
00174   const dg_template_t& templateDGElement(const int& i) const {return (*tmp_dgele)[i];};
00175   dg_template_t& templateDGElement(const int& i) {return (*tmp_dgele)[i];};
00176   int n_DGElement() const {return dgele.size();};
00177   const std::vector<dg_element_t>& dgElement() const {return dgele;};
00178   std::vector<dg_element_t>& dgElement() {return dgele;};
00179   const dg_element_t& dgElement(const int& i) const {return dgele[i];};
00180   dg_element_t& dgElement(const int& i) {return dgele[i];};
00181   virtual void buildDGElement();
00182   DGElementIterator beginDGElement() {return dgElement().begin();};
00183   DGElementIterator endDGElement() {return dgElement().end();};
00184   ConstDGElementIterator beginDGElement() const {return dgElement().begin();};
00185   ConstDGElementIterator endDGElement() const {return dgElement().end();};
00186 
00187 #ifdef __SERIALIZATION__
00188   template <class T> T *
00189   new_property(const dg_element_t& ele,
00190                const property_id_t<T>& pid) const {
00191     typedef RegularMesh<DIM,DOW> reg_mesh_t;
00192     const reg_mesh_t& mesh = dynamic_cast<const reg_mesh_t&>(this->mesh());
00193     return mesh.template new_property<T,DIM-1>(ele.geometry().index(), pid);
00194   }
00195   template <class T> T *
00196   get_property(const dg_element_t& ele,
00197                const property_id_t<T>& pid) const {
00198     typedef RegularMesh<DIM,DOW> reg_mesh_t;
00199     const reg_mesh_t& mesh = dynamic_cast<const reg_mesh_t&>(this->mesh());
00200     return mesh.template get_property<T,DIM-1>(ele.geometry().index(), pid);
00201   }
00202   template <class T> void
00203   free_property(const dg_element_t& ele,
00204                 const property_id_t<T>& pid) const {
00205     typedef RegularMesh<DIM,DOW> reg_mesh_t;
00206     const reg_mesh_t& mesh = dynamic_cast<const reg_mesh_t&>(this->mesh());
00207     mesh.template free_property<T,DIM-1>(ele.geometry().index(), pid);
00208   }
00209 
00210   using base_t::new_property;
00211   using base_t::get_property;
00212   using base_t::free_property;
00213 #endif // __SERIALIZATION__
00214   public:
00215   friend class DGElement<value_t,DIM,DOW,TDIM,TDIM1>;
00216   };
00217 
00218 
00220 
00229 template <typename value_type, int DIM>
00230   struct DGElementAdditionalData : public GeometryAdditionalData<DIM>
00231 {
00232  public:
00233   std::vector<std::vector<value_type> > basis_value0;
00234   std::vector<std::vector<value_type> > basis_value1;
00235   std::vector<std::vector<double> > unit_normal;
00236 };
00237 
00238 #endif // _DGFEMSpace_h_
00239 
00240 //
00241 // end of file