Home | Doxygen Documentation | Tutorials | Developer Tools (restricted)

bem/bform.hh
Go to the documentation of this file.
00001 /* Basic bilinear forms for the boundary element method
00002  */
00003 
00004 #ifndef bform_hh
00005 #define bform_hh
00006 
00007 #include "operator/bilinearForm.hh"
00008 #include "bemInt.hh"
00009 #include "bem/element.hh"
00010 
00011 namespace bem {
00012 
00013   // ************************************************************* LaplaceSL **
00014 
00018   template <class F = concepts::Real>
00019   class LaplaceSL : public concepts::BilinearForm<F> {
00020 
00022     LplCol006<F> qrA_;
00023     LplGal006<F> qrB_;
00024     LplGal018<F> qrC_;
00025 
00027     uint stroud_, gauss_;
00031     concepts::Real dist_;
00032 
00033   public:
00039     inline LaplaceSL(uint stroud = 0, uint gauss = 0,
00040          concepts::Real dist = 0.0);
00041 
00048     void operator()(const concepts::Element<F>& elmX,
00049         const concepts::Element<F>& elmY,
00050         concepts::ElementMatrix<F>& em);
00051 
00052     inline void operator()(const Dirac3d000<F>& elmX,
00053          const Linear3d000<F>& elmY,
00054          concepts::ElementMatrix<F>& em);
00055     inline void operator()(const Constant3d000<F>& elmX,
00056          const Constant3d000<F>& elmY,
00057          concepts::ElementMatrix<F>& em);
00058     inline void operator()(const Constant3d001<F>& elmX,
00059          const Constant3d001<F>& elmY,
00060          concepts::ElementMatrix<F>& em);
00061     inline void operator()(const Constant3d002<F>& elmX,
00062          const Constant3d002<F>& elmY,
00063          concepts::ElementMatrix<F>& em);
00064     inline void operator()(const Linear3d000<F>& elmX,
00065          const Linear3d000<F>& elmY,
00066          concepts::ElementMatrix<F>& em);
00067 
00068     virtual LaplaceSL* clone() const {
00069       return new LaplaceSL(stroud_, gauss_, dist_); }
00070   };
00071 
00072   template <class F>
00073   inline LaplaceSL<F>::LaplaceSL(uint stroud, uint gauss, concepts::Real dist)
00074   : stroud_(stroud), gauss_(gauss), dist_(dist) {
00075   }
00076 
00077   template <class F>
00078   inline void LaplaceSL<F>::operator()(const Dirac3d000<F>& elmX,
00079                const Linear3d000<F>& elmY,
00080                concepts::ElementMatrix<F>& em) {
00081     F m[9];
00082     qrA_(elmX, elmY, gauss_, dist_, m);
00083 
00084     uint n = elmX.T().n();
00085     F*   mm = m;
00086 
00087     em.resize(n, 3);
00088     for(uint i = 0; i < n; ++i) {
00089       em(i, 0) = *mm++; em(i, 1) = *mm++; em(i, 2) = *mm++;
00090     }
00091   }
00092 
00093   template <class F>
00094   inline void LaplaceSL<F>::operator()(const Constant3d000<F>& elmX,
00095                const Constant3d000<F>& elmY,
00096                concepts::ElementMatrix<F>& em) {
00097     F m;
00098     qrC_(elmX, elmY, stroud_, gauss_, dist_, &m);
00099 
00100     em.resize(1, 1);
00101     em(0, 0) = m;
00102   }
00103 
00104   template <class F>
00105   inline void LaplaceSL<F>::operator()(const Constant3d001<F>& elmX,
00106                const Constant3d001<F>& elmY,
00107                concepts::ElementMatrix<F>& em) {
00108     F m;
00109     qrC_(elmX, elmY, stroud_, gauss_, dist_, &m);
00110 
00111     em.resize(1, 1);
00112     em(0, 0) = m;
00113   }
00114 
00115   template <class F>
00116   inline void LaplaceSL<F>::operator()(const Constant3d002<F>& elmX,
00117                const Constant3d002<F>& elmY,
00118                concepts::ElementMatrix<F>& em) {
00119     F m;
00120     qrC_(elmX, elmY, stroud_, gauss_, dist_, &m);
00121 
00122     em.resize(1, 1);
00123     em(0, 0) = m;
00124   }
00125 
00126   template <class F>
00127   inline void LaplaceSL<F>::operator()(const Linear3d000<F>& elmX,
00128                const Linear3d000<F>& elmY,
00129                concepts::ElementMatrix<F>& em) {
00130     F m[9];
00131     qrB_(elmX, elmY, stroud_, gauss_, dist_, m);
00132 
00133     em.resize(3, 3);
00134     em(0, 0) = m[0]; em(0, 1) = m[1]; em(0, 2) = m[2];
00135     em(1, 0) = m[3]; em(1, 1) = m[4]; em(1, 2) = m[5];
00136     em(2, 0) = m[6]; em(2, 1) = m[7]; em(2, 2) = m[8];
00137   }
00138 
00139   // ************************************************************* LaplaceDL **
00140 
00144   template <class F = concepts::Real>
00145   class LaplaceDL : public concepts::BilinearForm<F> {
00146 
00148     LplCol007<F> qrA_;
00149     LplGal007<F> qrB_;
00150     LplGal014<F> qrC_;
00151 
00153     uint stroud_, gauss_;
00157     concepts::Real dist_;
00158 
00159   public:
00165     inline LaplaceDL(uint stroud = 0, uint gauss = 0,
00166          concepts::Real dist = 0.0);
00167 
00174     void operator()(const concepts::Element<F>& elmX,
00175         const concepts::Element<F>& elmY,
00176         concepts::ElementMatrix<F>& em);
00177 
00178     inline void operator()(const Dirac3d000<F>& elmX,
00179          const Linear3d000<F>& elmY,
00180          concepts::ElementMatrix<F>& em);
00181     inline void operator()(const Constant3d000<F>& elmX,
00182          const Constant3d000<F>& elmY,
00183          concepts::ElementMatrix<F>& em);
00184     inline void operator()(const Constant3d002<F>& elmX,
00185          const Constant3d002<F>& elmY,
00186          concepts::ElementMatrix<F>& em);
00187     inline void operator()(const Linear3d000<F>& elmX,
00188          const Linear3d000<F>& elmY,
00189          concepts::ElementMatrix<F>& em);
00190 
00191     virtual LaplaceDL* clone() const {
00192       return new LaplaceDL(stroud_, gauss_, dist_); }
00193   };
00194 
00195   template <class F>
00196   inline LaplaceDL<F>::LaplaceDL(uint stroud, uint gauss, concepts::Real dist)
00197     : stroud_(stroud), gauss_(gauss), dist_(dist) {
00198   }
00199 
00200   template <class F>
00201   inline void LaplaceDL<F>::operator()(const Dirac3d000<F>& elmX,
00202                const Linear3d000<F>& elmY,
00203                concepts::ElementMatrix<F>& em) {
00204     F m[9];
00205     qrA_(elmX, elmY, gauss_, dist_, m);
00206 
00207     uint n = elmX.T().n();
00208     F* mm = m;
00209 
00210     em.resize(n, 3);
00211     for(uint i = 0; i < n; ++i) {
00212       em(i, 0) = *mm++; em(i, 1) = *mm++; em(i, 2) = *mm++;
00213     }
00214   }
00215 
00216   template <class F>
00217   inline void LaplaceDL<F>::operator()(const Constant3d000<F>& elmX,
00218                const Constant3d000<F>& elmY,
00219                concepts::ElementMatrix<F>& em) {
00220     F m;
00221     qrC_(elmX, elmY, stroud_, gauss_, dist_, &m);
00222 
00223     em.resize(1, 1);
00224     em(0, 0) = m;
00225   }
00226 
00227   template <class F>
00228   inline void LaplaceDL<F>::operator()(const Constant3d002<F>& elmX,
00229                const Constant3d002<F>& elmY,
00230                concepts::ElementMatrix<F>& em) {
00231     F m;
00232     qrC_(elmX, elmY, stroud_, gauss_, dist_, &m);
00233 
00234     em.resize(1, 1);
00235     em(0, 0) = m;
00236   }
00237 
00238   template <class F>
00239   inline void LaplaceDL<F>::operator()(const Linear3d000<F>& elmX,
00240                const Linear3d000<F>& elmY,
00241                concepts::ElementMatrix<F>& em) {
00242     F m[9];
00243     qrB_(elmX, elmY, stroud_, gauss_, dist_, m);
00244 
00245     em.resize(3, 3);
00246     em(0, 0) = m[0]; em(0, 1) = m[1]; em(0, 2) = m[2];
00247     em(1, 0) = m[3]; em(1, 1) = m[4]; em(1, 2) = m[5];
00248     em(2, 0) = m[6]; em(2, 1) = m[7]; em(2, 2) = m[8];
00249   }
00250 
00251   // ************************************************************ LaplacePLD **
00252 
00257   template <class F = concepts::Real>
00258   class LaplacePLD : public concepts::BilinearForm<F> {
00259 
00261     LplGal007<F> qrA_;
00262     LplGal014<F> qrB_;
00263 
00265     uint stroud_, gauss_;
00269     concepts::Real dist_;
00270 
00271   public:
00277     inline LaplacePLD(uint stroud = 0, uint gauss = 0,
00278           concepts::Real dist = 0.0);
00279 
00286     void operator()(const concepts::Element<F>& elmX,
00287         const concepts::Element<F>& elmY,
00288         concepts::ElementMatrix<F>& em);
00289 
00290     inline void operator()(const Constant3d000<F>& elmX,
00291          const Constant3d000<F>& elmY,
00292          concepts::ElementMatrix<F>& em);
00293     inline void operator()(const Linear3d000<F>& elmX,
00294          const Linear3d000<F>& elmY,
00295          concepts::ElementMatrix<F>& em);
00296 
00297     virtual LaplacePLD* clone() const {
00298       return new LaplacePLD(stroud_, gauss_, dist_); }
00299   };
00300 
00301 
00302   template <class F>
00303   inline LaplacePLD<F>::LaplacePLD(uint stroud, uint gauss,
00304            concepts::Real dist)
00305     : stroud_(stroud), gauss_(gauss), dist_(dist) {
00306   }
00307 
00308   template <class F>
00309   inline void LaplacePLD<F>::operator()(const Constant3d000<F>& elmX,
00310           const Constant3d000<F>& elmY,
00311           concepts::ElementMatrix<F>& em) {
00312     F m;
00313     qrB_(elmY, elmX, stroud_, gauss_, dist_, &m);
00314 
00315     em.resize(1, 1);
00316     em(0, 0) = m;
00317   }
00318 
00319   template <class F>
00320   inline void LaplacePLD<F>::operator()(const Linear3d000<F>& elmX,
00321           const Linear3d000<F>& elmY,
00322           concepts::ElementMatrix<F>& em) {
00323     F m[9];
00324     qrA_(elmY, elmX, stroud_, gauss_, dist_, m);
00325 
00326     em.resize(3, 3);
00327     em(0, 0) = m[0]; em(0, 1) = m[1]; em(0, 2) = m[2];
00328     em(1, 0) = m[3]; em(1, 1) = m[4]; em(1, 2) = m[5];
00329     em(2, 0) = m[6]; em(2, 1) = m[7]; em(2, 2) = m[8];
00330   }
00331 
00332   // ************************************************************ LaplaceHyp **
00333 
00337   template <class F = concepts::Real>
00338   class LaplaceHyp : public concepts::BilinearForm<F> {
00339 
00341     LplGal008<F> qrA_;
00342 
00344     uint stroud_, gauss_;
00348     concepts::Real dist_;
00349 
00350   public:
00356     inline LaplaceHyp(uint stroud = 0, uint gauss = 0,
00357           concepts::Real dist = 0.0);
00358 
00365     void operator()(const concepts::Element<F>& elmX,
00366         const concepts::Element<F>& elmY,
00367         concepts::ElementMatrix<F>& em);
00368 
00369     inline void operator()(const Linear3d000<F>& elmX,
00370          const Linear3d000<F>& elmY,
00371          concepts::ElementMatrix<F>& em);
00372 
00373     virtual LaplaceHyp* clone() const {
00374       return new LaplaceHyp(stroud_, gauss_, dist_); }
00375   };
00376 
00377   template <class F>
00378   inline LaplaceHyp<F>::LaplaceHyp(uint stroud, uint gauss,
00379            concepts::Real dist)
00380     : stroud_(stroud), gauss_(gauss), dist_(dist) {
00381   }
00382 
00383   template <class F>
00384   inline void LaplaceHyp<F>::operator()(const Linear3d000<F>& elmX,
00385           const Linear3d000<F>& elmY,
00386           concepts::ElementMatrix<F>& em) {
00387     F m[9];
00388     qrA_(elmX, elmY, stroud_, gauss_, dist_, m);
00389 
00390     em.resize(3, 3);
00391     em(0, 0) = m[0]; em(0, 1) = m[1]; em(0, 2) = m[2];
00392     em(1, 0) = m[3]; em(1, 1) = m[4]; em(1, 2) = m[5];
00393     em(2, 0) = m[6]; em(2, 1) = m[7]; em(2, 2) = m[8];
00394   }
00395 
00396   // ************************************************************** Identity **
00397 
00401   template <class F = concepts::Real>
00402   class Identity : public concepts::BilinearForm<F> {
00403 
00404   public:
00411     void operator()(const concepts::Element<F>& elmX,
00412         const concepts::Element<F>& elmY,
00413         concepts::ElementMatrix<F>& em);
00414 
00415     inline void operator()(const Constant3d000<F>& elmX,
00416          const Constant3d000<F>& elmY,
00417          concepts::ElementMatrix<F>& em);
00418     inline void operator()(const Constant3d001<F>& elmX,
00419          const Constant3d001<F>& elmY,
00420          concepts::ElementMatrix<F>& em);
00421     inline void operator()(const Constant3d002<F>& elmX,
00422          const Constant3d002<F>& elmY,
00423          concepts::ElementMatrix<F>& em);
00424     void operator()(const Linear3d000<F>& elmX, const Linear3d000<F>& elmY,
00425         concepts::ElementMatrix<F>& em);
00426 
00427     virtual Identity* clone() const { return new Identity(); }
00428   };
00429 
00430   template <class F>
00431   inline void Identity<F>::operator()(const Constant3d000<F>& elmX,
00432               const Constant3d000<F>& elmY,
00433               concepts::ElementMatrix<F>& em) {
00434     em.resize(1, 1);
00435     em(0, 0) = (elmX.support().key() != elmY.support().key())
00436       ? 0.0 : elmX.n().length() * 0.5;
00437   }
00438 
00439   template <class F>
00440   inline void Identity<F>::operator()(const Constant3d001<F>& elmX,
00441               const Constant3d001<F>& elmY,
00442               concepts::ElementMatrix<F>& em) {
00443     em.resize(1, 1);
00444     em(0, 0) = (elmX.support().key() != elmY.support().key())
00445       ? 0.0 : elmX.n().length() * 0.5;
00446   }
00447 
00448   template <class F>
00449   inline void Identity<F>::operator()(const Constant3d002<F>& elmX,
00450               const Constant3d002<F>& elmY,
00451               concepts::ElementMatrix<F>& em) {
00452     em.resize(1, 1);
00453     em(0, 0) = (elmX.support().key() != elmY.support().key())
00454       ? 0.0 : 1.0;
00455   }
00456 
00457   // ************************************************************ ArbKrnl000 **
00458 
00465   template <class F, class K>
00466   class ArbKrnl000 : public concepts::BilinearForm<F> {
00467 
00469     ArbKrnlGal000<F, K> qrA_;
00470     ArbKrnlGal001<F, K> qrB_;
00471     K& k_;
00472   public:
00476     inline ArbKrnl000(K& k) : qrA_(k), qrB_(k), k_(k) {}
00477 
00484     void operator()(const concepts::Element<F>& elmX,
00485         const concepts::Element<F>& elmY,
00486         concepts::ElementMatrix<F>& em);
00487 
00488     inline void operator()(const Constant3d000<F>& elmX,
00489          const Constant3d000<F>& elmY,
00490          concepts::ElementMatrix<F>& em);
00491     inline void operator()(const Constant3d001<F>& elmX,
00492          const Constant3d001<F>& elmY,
00493          concepts::ElementMatrix<F>& em);
00494     inline void operator()(const Constant3d002<F>& elmX,
00495          const Constant3d002<F>& elmY,
00496          concepts::ElementMatrix<F>& em);
00497     inline void operator()(const Linear3d000<F>& elmX,
00498          const Linear3d000<F>& elmY,
00499          concepts::ElementMatrix<F>& em);
00500 
00501     virtual ArbKrnl000* clone() const { return new ArbKrnl000(k_); }
00502   };
00503 
00504   template <class F, class K>
00505   inline void ArbKrnl000<F, K>::operator()(const Constant3d000<F>& elmX,
00506              const Constant3d000<F>& elmY,
00507              concepts::ElementMatrix<F>& em) {
00508     F m;
00509     qrA_(elmX, elmY, &m);
00510 
00511     em.resize(1, 1);
00512     em(0, 0) = m;
00513   }
00514 
00515   template <class F, class K>
00516   inline void ArbKrnl000<F, K>::operator()(const Constant3d001<F>& elmX,
00517              const Constant3d001<F>& elmY,
00518              concepts::ElementMatrix<F>& em) {
00519     F m;
00520     qrA_(elmX, elmY, &m);
00521 
00522     em.resize(1, 1);
00523     em(0, 0) = m;
00524   }
00525 
00526   template <class F, class K>
00527   inline void ArbKrnl000<F, K>::operator()(const Constant3d002<F>& elmX,
00528              const Constant3d002<F>& elmY,
00529              concepts::ElementMatrix<F>& em) {
00530     F m;
00531     qrA_(elmX, elmY, &m);
00532 
00533     em.resize(1, 1);
00534     em(0, 0) = m;
00535   }
00536 
00537   template <class F, class K>
00538   inline void ArbKrnl000<F, K>::operator()(const Linear3d000<F>& elmX,
00539              const Linear3d000<F>& elmY,
00540              concepts::ElementMatrix<F>& em) {
00541     F m[9];
00542     qrB_(elmX, elmY, m);
00543 
00544     em.resize(3, 3);
00545     em(0, 0) = m[0];  em(0, 1) = m[1];  em(0, 2) = m[2];
00546     em(1, 0) = m[3];  em(1, 1) = m[4];  em(1, 2) = m[5];
00547     em(2, 0) = m[6];  em(2, 1) = m[7];  em(2, 2) = m[8];
00548   }
00549 
00550 } // namespace bem
00551 
00552 #endif // bform_hh

Home | Doxygen Documentation | Tutorials | Developer Tools (restricted)