AFEPack
BilinearOperator.h
浏览该文件的文档。
00001 
00002 //  BilinearOperator.h : by R.Lie
00003 //
00004 
00005 #ifndef _BilinearOperator_h_
00006 #define _BilinearOperator_h_
00007 
00008 #include <iostream>
00009 #include <vector>
00010 #include <iterator>
00011 #include <algorithm>
00012 
00013 #include <base/exceptions.h>
00014 #include <lac/full_matrix.h>
00015 #include <lac/sparsity_pattern.h>
00016 #include <lac/sparse_matrix.h>
00017 
00018 #include "Geometry.h"
00019 #include "TemplateElement.h"
00020 #include "FEMSpace.h"
00021 #include "HGeometry.h"
00022 
00023 
00072 template <int DIM, class value_type0, class value_type1=value_type0,
00073   int DOW=DIM, int TDIM0=DIM, int TDIM1 = DIM>
00074   class BilinearOperator : public SparseMatrix<double>
00075   {
00076   private:
00077   SparsityPattern                                       sparsity_pattern;
00078   FEMSpace<value_type0,DIM,DOW,TDIM0> *         fem_space0;
00079   FEMSpace<value_type1,DIM,DOW,TDIM1> *         fem_space1;
00080   int                                           n_dof0;
00081   int                                           n_dof1;
00082   int                                           n_max_coupling_dof;
00083   const std::vector<int> *                      element_dof0;
00084   const std::vector<int> *                      element_dof1;
00085   FullMatrix<double>                            element_matrix;
00086   int                                           algebric_accuracy;
00087   public:
00088   BilinearOperator();
00089   BilinearOperator(FEMSpace<value_type0,DIM,DOW,TDIM0>&, FEMSpace<value_type1,DIM,DOW,TDIM1>&);
00090   virtual ~BilinearOperator();
00091   public:
00092   void reinit(FEMSpace<value_type0,DIM,DOW,TDIM0>&, FEMSpace<value_type1,DIM,DOW,TDIM1>&);
00093   const FEMSpace<value_type0,DIM,DOW,TDIM0>& FEMSpace0() const {return *fem_space0;};
00094   const FEMSpace<value_type1,DIM,DOW,TDIM1>& FEMSpace1() const {return *fem_space1;};
00095   FEMSpace<value_type0,DIM,DOW,TDIM0> * FEMSpace0() {return fem_space0;};
00096   FEMSpace<value_type1,DIM,DOW,TDIM1> * FEMSpace1() {return fem_space1;};
00097   const SparsityPattern& getSparsityPattern() const {return sparsity_pattern;};
00098   SparsityPattern& getSparsityPattern() {return sparsity_pattern;};
00099   public:
00100   const int& nDof0() const {return n_dof0;};
00101   const int& nDof1() const {return n_dof1;};
00102   const int& nMaxCouplingDof() const {return n_max_coupling_dof;};
00103   int& nDof0() {return n_dof0;};
00104   int& nDof1() {return n_dof1;};
00105   int& nMaxCouplingDof() {return n_max_coupling_dof;};
00106 
00107   const std::vector<int>& elementDof0() const {return *element_dof0;};
00108   const std::vector<int>& elementDof1() const {return *element_dof1;};
00109   const int& elementDof0(const int& i) const {return (*element_dof0)[i];};
00110   const int& elementDof1(const int& i) const {return (*element_dof1)[i];};
00111 
00112   const FullMatrix<double>& elementMatrix() const {return element_matrix;};
00113   const double elementMatrix(const int& i, const int& j) const {return element_matrix(i,j);};
00114   FullMatrix<double>& elementMatrix() {return element_matrix;};
00115   double& elementMatrix(const int& i, const int& j) {return element_matrix(i,j);};
00116   void buildDofInfo();
00117   void addElementPattern();
00118   void addElementMatrix();
00119   void getElementPattern(const Element<value_type0,DIM,DOW,TDIM0>&, const Element<value_type1,DIM,DOW,TDIM1>&);
00120   public:
00121   const int& algebricAccuracy() const {return algebric_accuracy;};
00122   int& algebricAccuracy() {return algebric_accuracy;};
00123 
00124   virtual void getElementMatrix(const Element<value_type0,DIM,DOW,TDIM0>&, 
00125                                 const Element<value_type1,DIM,DOW,TDIM1>&,
00126                                 const typename ActiveElementPairIterator<DIM>::State state = ActiveElementPairIterator<DIM>::EQUAL) = 0;
00127   virtual void build();
00128   virtual void buildSparsityPattern();
00129   virtual void buildSparseMatrix();
00130   };
00131 
00137 template <int DIM, class value_type0, class value_type1=value_type0,
00138   int DOW=DIM, int TDIM0=DIM, int TDIM1 = DIM>
00139   class L2InnerProduct : public BilinearOperator<DIM, value_type0, value_type1, DOW, TDIM0, TDIM1>
00140   {
00141   public:
00142   L2InnerProduct() {};
00143   L2InnerProduct(FEMSpace<value_type0,DIM,DOW,TDIM0>& sp0, FEMSpace<value_type1,DIM,DOW,TDIM1>& sp1) :
00144   BilinearOperator<DIM, value_type0, value_type1, DOW, TDIM0, TDIM1>(sp0, sp1) {};
00145   virtual ~L2InnerProduct() {};
00146   public:
00147   void reinit(FEMSpace<value_type0,DIM,DOW,TDIM0>&, FEMSpace<value_type1,DIM,DOW,TDIM1>&);
00148   virtual void getElementMatrix(const Element<value_type0,DIM,DOW,TDIM0>&, 
00149                                 const Element<value_type1,DIM,DOW,TDIM1>&,
00150                                 const typename ActiveElementPairIterator<DIM>::State state = ActiveElementPairIterator<DIM>::EQUAL);
00151   };
00152 
00153 template <int DIM, class value_type, int DOW=DIM, int TDIM=DIM>
00154   class MassMatrix : public BilinearOperator<DIM, value_type, value_type, DOW, TDIM, TDIM>
00155   {
00156   public:
00157   MassMatrix() {};
00158   MassMatrix(FEMSpace<typename MassMatrix::value_type,DIM,DOW,TDIM>& sp) :
00159   BilinearOperator<DIM, typename MassMatrix::value_type, 
00160   typename MassMatrix::value_type, DOW, TDIM, TDIM>(sp, sp) {};
00161   virtual ~MassMatrix() {};
00162   public:
00163   void reinit(FEMSpace<typename MassMatrix::value_type,DIM,DOW,TDIM>& sp);
00164   virtual void getElementMatrix(const Element<typename MassMatrix::value_type,DIM,DOW,TDIM>& e0, 
00165                                 const Element<typename MassMatrix::value_type,DIM,DOW,TDIM>& e1,
00166                                 const typename ActiveElementPairIterator<DIM>::State state = ActiveElementPairIterator<DIM>::EQUAL);
00167   };
00168 
00169 template <int DIM, class value_type, int DOW=DIM, int TDIM=DIM>
00170   class StiffMatrix : public BilinearOperator<DIM, value_type, value_type, DOW, TDIM, TDIM>
00171   {
00172   public:
00173   StiffMatrix() {};
00174   StiffMatrix(FEMSpace<typename StiffMatrix::value_type,DIM,DOW,TDIM>& sp) :
00175   BilinearOperator<DIM, typename StiffMatrix::value_type, 
00176   typename StiffMatrix::value_type,DOW,TDIM,TDIM>(sp, sp) {};
00177   virtual ~StiffMatrix() {};
00178   public:
00179   void reinit(FEMSpace<typename StiffMatrix::value_type,DIM,DOW,TDIM>& sp);
00180   virtual void getElementMatrix(const Element<typename StiffMatrix::value_type,DIM,DOW,TDIM>&, 
00181                                 const Element<typename StiffMatrix::value_type,DIM,DOW,TDIM>&,
00182                                 const typename ActiveElementPairIterator<DIM>::State state = ActiveElementPairIterator<DIM>::EQUAL);
00183   };
00184 
00185 #endif
00186 
00187 //
00188 // end of file