NGSolve
4.9
|
00001 #ifndef FILE_INTEGRATOR 00002 #define FILE_INTEGRATOR 00003 00004 /*********************************************************************/ 00005 /* File: integrator.hpp */ 00006 /* Author: Joachim Schoeberl */ 00007 /* Date: 25. Mar. 2000 */ 00008 /*********************************************************************/ 00009 00010 namespace ngfem 00011 { 00012 00013 /* 00014 bilinear-form and linear-form integrators 00015 */ 00016 00017 00022 class NGS_DLL_HEADER Integrator 00023 { 00024 protected: 00026 BitArray definedon; 00027 00029 int integration_order; 00030 00031 // lower bound for the integration order 00032 int higher_integration_order; 00033 00035 static int common_integration_order; 00036 00038 bool const_coef; 00039 00041 string name; 00042 00044 Array < FlatVector < double > * > curve_ips; 00045 Array < FlatVector < double > * > curve_ip_tangents; 00046 Array <int> continuous_curveparts; 00047 00048 int cachecomp; 00049 00050 00051 protected: 00052 void DeleteCurveIPs ( void ); 00053 00054 public: 00056 Integrator() throw (); 00057 00059 virtual ~Integrator(); 00060 00062 virtual bool BoundaryForm () const = 0; 00063 00065 virtual bool SkeletonForm () const {return 0;} 00066 00068 bool DefinedOn (int mat) const; 00069 00071 void SetDefinedOn (const BitArray & adefinedon); 00072 00073 bool DefinedOnSubdomainsOnly() const 00074 { return definedon.Size() != 0; } 00075 00076 00077 00079 static void SetCommonIntegrationOrder (int cio) 00080 { 00081 common_integration_order = cio; 00082 } 00083 00084 static int GetCommonIntegrationOrder () 00085 { 00086 return common_integration_order; 00087 } 00088 00090 void SetHigherIntegrationOrder (int io) 00091 { 00092 higher_integration_order = io; 00093 } 00094 00096 void SetIntegrationOrder (int io) 00097 { 00098 integration_order = io; 00099 } 00100 00102 int GetIntegrationOrder (void) const 00103 { 00104 return integration_order; 00105 } 00106 00108 void SetConstantCoefficient (bool acc = 1) 00109 { const_coef = acc; } 00110 00112 virtual int DimElement () const { return -1; } 00113 00115 virtual int DimSpace () const { return -1; } 00116 00118 void SetName (const string & aname); 00120 virtual string Name () const; 00121 00123 virtual void CheckElement (const FiniteElement & el) const { ; } 00124 00125 00126 00127 // special hacks by Markus 00128 bool IntegrationAlongCurve (void) const 00129 { return curve_ips.Size() > 0; } 00130 00131 void SetIntegrationAlongCurve ( const int npoints ); 00132 00133 void UnSetIntegrationAlongCurve ( void ); 00134 00135 int NumCurvePoints(void) const 00136 { return curve_ips.Size(); } 00137 00138 FlatVector<double> & CurvePoint(const int i) 00139 { return *(curve_ips[i]); } 00140 00141 const FlatVector<double> & CurvePoint(const int i) const 00142 { return *(curve_ips[i]); } 00143 00144 FlatVector<double> & CurvePointTangent(const int i) 00145 { return *(curve_ip_tangents[i]); } 00146 00147 const FlatVector<double> & CurvePointTangent(const int i) const 00148 { return *(curve_ip_tangents[i]); } 00149 00150 int GetNumCurveParts(void) const; 00151 int GetStartOfCurve(const int i) const; 00152 int GetEndOfCurve(const int i) const; 00153 00154 void AppendCurvePoint(const FlatVector<double> & point); 00155 void AppendCurvePoint(const FlatVector<double> & point, const FlatVector<double> & tangent); 00156 void SetCurveClearance(void); 00157 00158 00159 00160 virtual void SetCacheComp(const int comp) 00161 { cachecomp = comp; } 00162 00163 virtual int CacheComp(void) const 00164 { return cachecomp; } 00165 00166 virtual void SetFileName(const string & filename) 00167 { 00168 cerr << "SetFileName not defined for Integrator base class" << endl; 00169 } 00170 }; 00171 00172 00173 00174 00180 class NGS_DLL_HEADER BilinearFormIntegrator : public Integrator 00181 { 00182 public: 00183 // typedef double TSCAL; 00185 BilinearFormIntegrator () throw (); 00187 virtual ~BilinearFormIntegrator (); 00188 00190 virtual int DimFlux () const { return -1; } 00191 00195 virtual void 00196 CalcElementMatrix (const FiniteElement & fel, 00197 const ElementTransformation & eltrans, 00198 FlatMatrix<double> & elmat, 00199 LocalHeap & lh) const 00200 { 00201 AssembleElementMatrix (fel, eltrans, elmat, lh); 00202 } 00203 00208 virtual void 00209 CalcElementMatrix (const FiniteElement & fel, 00210 const ElementTransformation & eltrans, 00211 FlatMatrix<Complex> & elmat, 00212 LocalHeap & lh) const; 00213 00214 00215 /* 00216 Computes the element matrix. 00217 should be renamed to CalcElementMatrix 00218 */ 00219 virtual void 00220 AssembleElementMatrix (const FiniteElement & fel, 00221 const ElementTransformation & eltrans, 00222 FlatMatrix<double> & elmat, 00223 LocalHeap & lh) const; 00224 00225 /* 00226 virtual void 00227 AssembleElementMatrix (const FiniteElement & fel, 00228 const ElementTransformation & eltrans, 00229 FlatMatrix<Complex> & elmat, 00230 LocalHeap & lh) const; 00231 */ 00232 00233 00234 00235 virtual void 00236 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00237 const FiniteElement & bfel_master_element, 00238 const FiniteElement & bfel_slave, 00239 const ElementTransformation & eltrans_master, 00240 const ElementTransformation & eltrans_master_element, 00241 const ElementTransformation & eltrans_slave, 00242 const IntegrationPoint & ip_master, 00243 const IntegrationPoint & ip_master_element, 00244 const IntegrationPoint & ip_slave, 00245 FlatMatrix<double> & elmat, 00246 LocalHeap & lh) const 00247 {;} 00248 virtual void 00249 ApplyElementMatrixIndependent (const FiniteElement & bfel_master, 00250 const FiniteElement & bfel_master_element, 00251 const FiniteElement & bfel_slave, 00252 const ElementTransformation & eltrans_master, 00253 const ElementTransformation & eltrans_master_element, 00254 const ElementTransformation & eltrans_slave, 00255 const IntegrationPoint & ip_master, 00256 const IntegrationPoint & ip_master_element, 00257 const IntegrationPoint & ip_slave, 00258 const FlatVector<double> & elx, 00259 Vector<double> & result, 00260 LocalHeap & lh) const 00261 {;} 00262 virtual void 00263 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00264 const FiniteElement & bfel_master_element, 00265 const FiniteElement & bfel_slave, 00266 const ElementTransformation & eltrans_master, 00267 const ElementTransformation & eltrans_master_element, 00268 const ElementTransformation & eltrans_slave, 00269 const IntegrationPoint & ip_master, 00270 const IntegrationPoint & ip_master_element, 00271 const IntegrationPoint & ip_slave, 00272 FlatMatrix<Complex> & elmat, 00273 LocalHeap & lh) const 00274 { 00275 FlatMatrix<double> rmat; 00276 AssembleElementMatrixIndependent(bfel_master,bfel_master_element,bfel_slave, 00277 eltrans_master, eltrans_master_element, eltrans_slave, 00278 ip_master, ip_master_element, ip_slave, 00279 rmat, lh); 00280 elmat.AssignMemory(rmat.Height(), rmat.Width(), lh); 00281 elmat = rmat; 00282 } 00283 00284 virtual void 00285 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00286 const FiniteElement & bfel_slave, 00287 const ElementTransformation & eltrans_master, 00288 const ElementTransformation & eltrans_slave, 00289 const IntegrationPoint & ip_master, 00290 const IntegrationPoint & ip_slave, 00291 FlatMatrix<double> & elmat, 00292 LocalHeap & lh) const 00293 {;} 00294 virtual void 00295 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00296 const FiniteElement & bfel_slave, 00297 const ElementTransformation & eltrans_master, 00298 const ElementTransformation & eltrans_slave, 00299 const IntegrationPoint & ip_master, 00300 const IntegrationPoint & ip_slave, 00301 FlatMatrix<Complex> & elmat, 00302 LocalHeap & lh) const 00303 { 00304 FlatMatrix<double> rmat; 00305 AssembleElementMatrixIndependent(bfel_master,bfel_slave, 00306 eltrans_master, eltrans_slave, 00307 ip_master, ip_slave, 00308 rmat, lh); 00309 elmat.AssignMemory(rmat.Height(), rmat.Width(), lh); 00310 elmat = rmat; 00311 } 00312 00313 00314 virtual void 00315 CalcElementMatrixDiag (const FiniteElement & fel, 00316 const ElementTransformation & eltrans, 00317 FlatVector<double> & diag, 00318 LocalHeap & lh) const; 00319 00320 00321 00322 00323 virtual void 00324 CalcLinearizedElementMatrix (const FiniteElement & fel, 00325 const ElementTransformation & eltrans, 00326 FlatVector<double> & elveclin, 00327 FlatMatrix<double> & elmat, 00328 LocalHeap & lh) const; 00329 00330 virtual void 00331 CalcLinearizedElementMatrix (const FiniteElement & fel, 00332 const ElementTransformation & eltrans, 00333 FlatVector<Complex> & elveclin, 00334 FlatMatrix<Complex> & elmat, 00335 LocalHeap & lh) const; 00336 00337 00338 virtual void * 00339 PrecomputeData (const FiniteElement & fel, 00340 const ElementTransformation & eltrans, 00341 LocalHeap & lh) const { return 0; } 00342 00343 00344 virtual void 00345 ApplyElementMatrix (const FiniteElement & fel, 00346 const ElementTransformation & eltrans, 00347 const FlatVector<double> & elx, 00348 FlatVector<double> & ely, 00349 void * precomputed, 00350 LocalHeap & lh) const; 00351 00352 virtual void 00353 ApplyElementMatrix (const FiniteElement & fel, 00354 const ElementTransformation & eltrans, 00355 const FlatVector<Complex> & elx, 00356 FlatVector<Complex> & ely, 00357 void * precomputed, 00358 LocalHeap & lh) const; 00359 00360 /* 00361 template < int S, class T> 00362 void ApplyElementMatrix (const FiniteElement & fel, 00363 const ElementTransformation & eltrans, 00364 const FlatVector< Vec<S,T> > & elx, 00365 FlatVector< Vec<S,T> > & ely, 00366 void * precomputed, 00367 LocalHeap & lh) const 00368 { 00369 //cout << "call baseclass ApplyElementMatrix, type = " << typeid(*this).name() << endl; 00370 FlatMatrix<T> mat; 00371 CalcElementMatrix (fel, eltrans, mat, lh); 00372 ely = mat * elx; 00373 } 00374 */ 00375 00376 00377 virtual void 00378 ApplyLinearizedElementMatrix (const FiniteElement & fel, 00379 const ElementTransformation & eltrans, 00380 const FlatVector<double> & ellin, 00381 const FlatVector<double> & elx, 00382 FlatVector<double> & ely, 00383 LocalHeap & lh) const; 00384 00385 virtual void 00386 ApplyLinearizedElementMatrix (const FiniteElement & fel, 00387 const ElementTransformation & eltrans, 00388 const FlatVector<Complex> & ellin, 00389 const FlatVector<Complex> & elx, 00390 FlatVector<Complex> & ely, 00391 LocalHeap & lh) const; 00392 00393 00394 00395 virtual double Energy (const FiniteElement & fel, 00396 const ElementTransformation & eltrans, 00397 const FlatVector<double> & elx, 00398 LocalHeap & lh) const; 00399 00400 virtual double Energy (const FiniteElement & fel, 00401 const ElementTransformation & eltrans, 00402 const FlatVector<Complex> & elx, 00403 LocalHeap & lh) const; 00404 00405 00406 00407 00408 virtual void 00409 CalcFlux (const FiniteElement & fel, 00410 const BaseMappedIntegrationPoint & bmip, 00411 const FlatVector<double> & elx, 00412 FlatVector<double> & flux, 00413 bool applyd, 00414 LocalHeap & lh) const; 00415 00416 00417 virtual void 00418 CalcFlux (const FiniteElement & fel, 00419 const BaseMappedIntegrationPoint & bmip, 00420 const FlatVector<Complex> & elx, 00421 FlatVector<Complex> & flux, 00422 bool applyd, 00423 LocalHeap & lh) const; 00424 00425 00426 virtual void 00427 CalcFlux (const FiniteElement & fel, 00428 const BaseMappedIntegrationRule & mir, 00429 const FlatVector<double> & elx, 00430 FlatMatrix<double> & flux, 00431 bool applyd, 00432 LocalHeap & lh) const; 00433 00434 00435 virtual void 00436 CalcFlux (const FiniteElement & fel, 00437 const BaseMappedIntegrationRule & mir, 00438 const FlatVector<Complex> & elx, 00439 FlatMatrix<Complex> & flux, 00440 bool applyd, 00441 LocalHeap & lh) const; 00442 00443 00444 virtual void 00445 CalcFlux (const FiniteElement & fel, 00446 const FiniteElement & felflux, 00447 const ElementTransformation & eltrans, 00448 const FlatVector<> & elx, 00449 FlatVector<> & flux, 00450 bool applyd, 00451 LocalHeap & lh) const; 00452 00453 00454 virtual void 00455 CalcFluxMulti (const FiniteElement & fel, 00456 const BaseMappedIntegrationPoint & bmip, 00457 int m, 00458 const FlatVector<double> & elx, 00459 FlatVector<double> & flux, 00460 bool applyd, 00461 LocalHeap & lh) const; 00462 00463 00464 virtual void 00465 CalcFluxMulti (const FiniteElement & fel, 00466 const BaseMappedIntegrationPoint & bmip, 00467 int m, 00468 const FlatVector<Complex> & elx, 00469 FlatVector<Complex> & flux, 00470 bool applyd, 00471 LocalHeap & lh) const; 00472 00473 00474 virtual void 00475 ApplyBTrans (const FiniteElement & fel, 00476 const BaseMappedIntegrationPoint & bmip, 00477 const FlatVector<double> & elx, 00478 FlatVector<double> & ely, 00479 LocalHeap & lh) const; 00480 00481 virtual void 00482 ApplyBTrans (const FiniteElement & fel, 00483 const BaseMappedIntegrationPoint & bmip, 00484 const FlatVector<Complex> & elx, 00485 FlatVector<Complex> & ely, 00486 LocalHeap & lh) const; 00487 00488 00489 virtual void 00490 ApplyBTrans (const FiniteElement & fel, 00491 const BaseMappedIntegrationRule & mir, 00492 const FlatMatrix<double> & elx, 00493 FlatVector<double> & ely, 00494 LocalHeap & lh) const; 00495 00496 virtual void 00497 ApplyBTrans (const FiniteElement & fel, 00498 const BaseMappedIntegrationRule & mir, 00499 const FlatMatrix<Complex> & elx, 00500 FlatVector<Complex> & ely, 00501 LocalHeap & lh) const; 00502 00503 00504 00505 virtual void ApplyDMat (const FiniteElement & bfel, 00506 const BaseMappedIntegrationPoint & bmip, 00507 const FlatVector<double> & elx, 00508 FlatVector<double> & eldx, 00509 LocalHeap & lh) const; 00510 00511 virtual void ApplyDMat (const FiniteElement & bfel, 00512 const BaseMappedIntegrationPoint & bmip, 00513 const FlatVector<Complex> & elx, 00514 FlatVector<Complex> & eldx, 00515 LocalHeap & lh) const; 00516 00517 virtual void ApplyDMat (const FiniteElement & bfel, 00518 const BaseMappedIntegrationRule & mir, 00519 const FlatMatrix<double> & elx, 00520 FlatMatrix<double> & eldx, 00521 LocalHeap & lh) const; 00522 00523 virtual void ApplyDMat (const FiniteElement & bfel, 00524 const BaseMappedIntegrationRule & mir, 00525 const FlatMatrix<Complex> & elx, 00526 FlatMatrix<Complex> & eldx, 00527 LocalHeap & lh) const; 00528 00529 virtual void ApplyDMatInv (const FiniteElement & bfel, 00530 const BaseMappedIntegrationPoint & bmip, 00531 const FlatVector<double> & elx, 00532 FlatVector<double> & eldx, 00533 LocalHeap & lh) const; 00534 00535 virtual void ApplyDMatInv (const FiniteElement & bfel, 00536 const BaseMappedIntegrationPoint & bmip, 00537 const FlatVector<Complex> & elx, 00538 FlatVector<Complex> & eldx, 00539 LocalHeap & lh) const; 00540 00541 virtual void ApplyDMatInv (const FiniteElement & bfel, 00542 const BaseMappedIntegrationRule & mir, 00543 const FlatMatrix<double> & elx, 00544 FlatMatrix<double> & eldx, 00545 LocalHeap & lh) const; 00546 00547 virtual void ApplyDMatInv (const FiniteElement & bfel, 00548 const BaseMappedIntegrationRule & mir, 00549 const FlatMatrix<Complex> & elx, 00550 FlatMatrix<Complex> & eldx, 00551 LocalHeap & lh) const; 00552 00553 00554 00555 /* 00556 virtual const IntegrationRule & GetIntegrationRule (const FiniteElement & fel, 00557 const bool use_higher_integration_order = false) const; 00558 */ 00559 }; 00560 00561 00562 /* 00563 class FacetNeighbourElInfo{ 00564 public: 00565 //finite Element of the neighbour element 00566 const FiniteElement & volumefel; 00567 //local Facet Number from volumeElements view 00568 int LocalFacetNr; 00569 //Transformation of the neighbouring element 00570 const ElementTransformation & eltrans; 00571 //Vertices of the Element 00572 FlatArray<int> & ElVertices; 00573 bool nonempty; 00574 FacetNeighbourElInfo(const FiniteElement & vvolumefel, int lLocalFacetNr, 00575 const ElementTransformation & eeltrans, 00576 FlatArray<int> & eElVertices):volumefel(vvolumefel), 00577 LocalFacetNr(lLocalFacetNr), eltrans(eeltrans), 00578 ElVertices(eElVertices), nonempty(true) 00579 {;} 00580 FacetNeighbourElInfo():nonempty(false){;}; 00581 bool IsEmpty(){return !nonempty;} 00582 }; 00583 */ 00584 00585 class NGS_DLL_HEADER FacetBilinearFormIntegrator : public BilinearFormIntegrator 00586 { 00587 public: 00588 00589 FacetBilinearFormIntegrator(Array<CoefficientFunction*> & coeffs) 00590 : BilinearFormIntegrator() { ; } 00591 00592 ~FacetBilinearFormIntegrator() { ; } 00593 00594 00595 00596 00597 virtual bool BoundaryForm () const 00598 { return 0; } 00599 00600 virtual bool SkeletonForm () const 00601 { return 1; } 00602 00603 virtual void CalcElementMatrix (const FiniteElement & fel, 00604 const ElementTransformation & eltrans, 00605 FlatMatrix<double> & elmat, 00606 LocalHeap & lh) { 00607 throw Exception ("FacetBilinearFormIntegrator can not assemble volumetric element matrices!"); 00608 } 00609 00610 00611 virtual void 00612 CalcFacetMatrix (const FiniteElement & volumefel1, int LocalFacetNr1, 00613 const ElementTransformation & eltrans1, FlatArray<int> & ElVertices1, 00614 const FiniteElement & volumefel2, int LocalFacetNr2, 00615 const ElementTransformation & eltrans2, FlatArray<int> & ElVertices2, 00616 FlatMatrix<double> & elmat, 00617 LocalHeap & lh) const{ 00618 throw Exception ("FacetBilinearFormIntegrator::CalcFacetMatrix for inner facets not implemented!"); 00619 } 00620 virtual void 00621 CalcFacetMatrix (const FiniteElement & volumefel1, int LocalFacetNr1, 00622 const ElementTransformation & eltrans1, FlatArray<int> & ElVertices1, 00623 const FiniteElement & volumefel2, int LocalFacetNr2, 00624 const ElementTransformation & eltrans2, FlatArray<int> & ElVertices2, 00625 FlatMatrix<Complex> & elmat, 00626 LocalHeap & lh) const{ 00627 throw Exception ("FacetBilinearFormIntegrator::CalcFacetMatrix<Complex> for inner facets not implemented!"); 00628 } 00629 00630 00631 virtual void 00632 CalcFacetMatrix (const FiniteElement & volumefel, int LocalFacetNr, 00633 const ElementTransformation & eltrans, FlatArray<int> & ElVertices, 00634 const ElementTransformation & seltrans, 00635 FlatMatrix<double> & elmat, 00636 LocalHeap & lh) const{ 00637 throw Exception ("FacetBilinearFormIntegrator::CalcFacetMatrix for boundary facets not implemented!"); 00638 } 00639 virtual void 00640 CalcFacetMatrix (const FiniteElement & volumefel, int LocalFacetNr, 00641 const ElementTransformation & eltrans, FlatArray<int> & ElVertices, 00642 const ElementTransformation & seltrans, 00643 FlatMatrix<Complex> & elmat, 00644 LocalHeap & lh) const{ 00645 throw Exception ("FacetBilinearFormIntegrator::CalcFacetMatrix<Complex> for boundary facets not implemented!"); 00646 } 00647 00648 00649 }; 00650 00651 00652 00653 00654 /* 00655 // returns coef * Identity-matrix 00656 class RegularizationBilinearFormIntegrator : public BilinearFormIntegrator 00657 { 00658 CoefficientFunction * coef; 00659 public: 00660 RegularizationBilinearFormIntegrator (CoefficientFunction * acoef) 00661 : coef(acoef) 00662 { ; } 00663 00664 virtual bool BoundaryForm () const { return 0; } 00665 virtual int DimFlux () const { return 1; } 00666 virtual int DimElement () const { return 3; } 00667 virtual int DimSpace () const { return 3; } 00668 00669 virtual void 00670 CalcElementMatrix (const FiniteElement & fel, 00671 const ElementTransformation & eltrans, 00672 FlatMatrix<double> & elmat, 00673 LocalHeap & lh) const 00674 { 00675 int ndof = fel.GetNDof(); 00676 00677 elmat = 0; 00678 MappedIntegrationPoint<3,3> mip (GetIntegrationRules().SelectIntegrationRule (fel.ElementType(), 0)[0], 00679 eltrans); 00680 double val = Evaluate (*coef, mip); 00681 elmat = 0; 00682 for (int i = 0; i < ndof; i++) 00683 elmat(i,i) = val; 00684 } 00685 00686 virtual string Name () const { return "Regularization"; } 00687 }; 00688 */ 00689 00690 00691 00692 00693 class NGS_DLL_HEADER BlockBilinearFormIntegrator : public BilinearFormIntegrator 00694 { 00695 BilinearFormIntegrator & bfi; 00696 int dim; 00697 int comp; 00698 public: 00699 BlockBilinearFormIntegrator (BilinearFormIntegrator & abfi, int adim, int acomp); 00700 BlockBilinearFormIntegrator (BilinearFormIntegrator & abfi, int adim); 00701 virtual ~BlockBilinearFormIntegrator (); 00702 00703 virtual bool BoundaryForm () const 00704 { return bfi.BoundaryForm(); } 00705 00706 virtual int DimFlux () const 00707 { return (comp == -1) ? dim * bfi.DimFlux() : bfi.DimFlux(); } 00708 00709 const BilinearFormIntegrator & Block () const { return bfi; } 00710 00711 virtual void 00712 CalcElementMatrix (const FiniteElement & bfel, 00713 const ElementTransformation & eltrans, 00714 FlatMatrix<double> & elmat, 00715 LocalHeap & lh) const; 00716 00717 virtual void 00718 CalcElementMatrix (const FiniteElement & bfel, 00719 const ElementTransformation & eltrans, 00720 FlatMatrix<Complex> & elmat, 00721 LocalHeap & lh) const; 00722 00723 virtual void 00724 ApplyElementMatrix (const FiniteElement & bfel, 00725 const ElementTransformation & eltrans, 00726 const FlatVector<double> & elx, 00727 FlatVector<double> & ely, 00728 void * precomputed, 00729 LocalHeap & lh) const; 00730 00731 virtual void 00732 ApplyElementMatrix (const FiniteElement & bfel, 00733 const ElementTransformation & eltrans, 00734 const FlatVector<Complex> & elx, 00735 FlatVector<Complex> & ely, 00736 void * precomputed, 00737 LocalHeap & lh) const; 00738 00739 virtual void 00740 CalcLinearizedElementMatrix (const FiniteElement & bfel, 00741 const ElementTransformation & eltrans, 00742 FlatVector<double> & elveclin, 00743 FlatMatrix<double> & elmat, 00744 LocalHeap & lh) const; 00745 virtual void 00746 CalcLinearizedElementMatrix (const FiniteElement & bfel, 00747 const ElementTransformation & eltrans, 00748 FlatVector<Complex> & elveclin, 00749 FlatMatrix<Complex> & elmat, 00750 LocalHeap & lh) const; 00751 00752 00753 virtual void 00754 CalcFlux (const FiniteElement & fel, 00755 const ElementTransformation & eltrans, 00756 const IntegrationPoint & ip, 00757 const FlatVector<double> & elx, 00758 FlatVector<double> & flux, 00759 bool applyd, 00760 LocalHeap & lh) const; 00761 00762 virtual void 00763 CalcFlux (const FiniteElement & fel, 00764 const ElementTransformation & eltrans, 00765 const IntegrationPoint & ip, 00766 const FlatVector<Complex> & elx, 00767 FlatVector<Complex> & flux, 00768 bool applyd, 00769 LocalHeap & lh) const; 00770 00771 00772 00773 virtual void 00774 CalcFlux (const FiniteElement & fel, 00775 const BaseMappedIntegrationPoint & bmip, 00776 const FlatVector<double> & elx, 00777 FlatVector<double> & flux, 00778 bool applyd, 00779 LocalHeap & lh) const; 00780 00781 virtual void 00782 CalcFlux (const FiniteElement & fel, 00783 const BaseMappedIntegrationPoint & bmip, 00784 const FlatVector<Complex> & elx, 00785 FlatVector<Complex> & flux, 00786 bool applyd, 00787 LocalHeap & lh) const; 00788 00789 virtual void 00790 CalcFlux (const FiniteElement & fel, 00791 const BaseMappedIntegrationRule & mir, 00792 const FlatVector<double> & elx, 00793 FlatMatrix<double> & flux, 00794 bool applyd, 00795 LocalHeap & lh) const; 00796 00797 00798 00799 virtual void 00800 ApplyBTrans (const FiniteElement & bfel, 00801 const BaseMappedIntegrationPoint & bmip, 00802 const FlatVector<double> & elx, 00803 FlatVector<double> & ely, 00804 LocalHeap & lh) const; 00805 00806 virtual void 00807 ApplyBTrans (const FiniteElement & bfel, 00808 const BaseMappedIntegrationPoint & bmip, 00809 const FlatVector<Complex> & elx, 00810 FlatVector<Complex> & ely, 00811 LocalHeap & lh) const; 00812 00813 virtual void 00814 ApplyBTrans (const FiniteElement & fel, 00815 const BaseMappedIntegrationRule & mir, 00816 const FlatMatrix<double> & elx, 00817 FlatVector<double> & ely, 00818 LocalHeap & lh) const; 00819 00820 virtual double Energy (const FiniteElement & fel, 00821 const ElementTransformation & eltrans, 00822 const FlatVector<double> & elx, 00823 LocalHeap & lh) const; 00824 00825 virtual string Name () const; 00826 }; 00827 00828 00829 00830 00831 00832 00833 00834 class NGS_DLL_HEADER ComplexBilinearFormIntegrator : public BilinearFormIntegrator 00835 { 00836 const BilinearFormIntegrator & bfi; 00837 Complex factor; 00838 public: 00839 ComplexBilinearFormIntegrator (const BilinearFormIntegrator & abfi, 00840 Complex afactor); 00841 00842 virtual bool BoundaryForm () const 00843 { return bfi.BoundaryForm(); } 00844 00845 virtual int DimFlux () const 00846 { return bfi.DimFlux(); } 00847 virtual int DimElement () const 00848 { return bfi.DimElement(); } 00849 virtual int DimSpace () const 00850 { return bfi.DimSpace(); } 00851 00852 virtual void GetFactor(Complex & fac) const {fac = factor;} 00853 virtual void GetFactor(double & fac) const {fac = factor.real();} 00854 00855 virtual const BilinearFormIntegrator & GetBFI(void) const {return bfi;} 00856 00857 virtual void CheckElement (const FiniteElement & el) const { bfi.CheckElement(el); } 00858 00859 00860 virtual void 00861 CalcElementMatrix (const FiniteElement & fel, 00862 const ElementTransformation & eltrans, 00863 FlatMatrix<double> & elmat, 00864 LocalHeap & lh) const; 00865 00866 virtual void 00867 CalcElementMatrix (const FiniteElement & fel, 00868 const ElementTransformation & eltrans, 00869 FlatMatrix<Complex> & elmat, 00870 LocalHeap & lh) const; 00871 00872 virtual void 00873 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00874 const FiniteElement & bfel_master_element, 00875 const FiniteElement & bfel_slave, 00876 const ElementTransformation & eltrans_master, 00877 const ElementTransformation & eltrans_master_element, 00878 const ElementTransformation & eltrans_slave, 00879 const IntegrationPoint & ip_master, 00880 const IntegrationPoint & ip_master_element, 00881 const IntegrationPoint & ip_slave, 00882 FlatMatrix<double> & elmat, 00883 LocalHeap & lh) const; 00884 00885 virtual void 00886 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00887 const FiniteElement & bfel_master_element, 00888 const FiniteElement & bfel_slave, 00889 const ElementTransformation & eltrans_master, 00890 const ElementTransformation & eltrans_master_element, 00891 const ElementTransformation & eltrans_slave, 00892 const IntegrationPoint & ip_master, 00893 const IntegrationPoint & ip_master_element, 00894 const IntegrationPoint & ip_slave, 00895 FlatMatrix<Complex> & elmat, 00896 LocalHeap & lh) const; 00897 00898 virtual void 00899 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00900 const FiniteElement & bfel_slave, 00901 const ElementTransformation & eltrans_master, 00902 const ElementTransformation & eltrans_slave, 00903 const IntegrationPoint & ip_master, 00904 const IntegrationPoint & ip_slave, 00905 FlatMatrix<double> & elmat, 00906 LocalHeap & lh) const; 00907 00908 virtual void 00909 AssembleElementMatrixIndependent (const FiniteElement & bfel_master, 00910 const FiniteElement & bfel_slave, 00911 const ElementTransformation & eltrans_master, 00912 const ElementTransformation & eltrans_slave, 00913 const IntegrationPoint & ip_master, 00914 const IntegrationPoint & ip_slave, 00915 FlatMatrix<Complex> & elmat, 00916 LocalHeap & lh) const; 00917 00918 00919 00920 virtual void 00921 ApplyElementMatrix (const FiniteElement & fel, 00922 const ElementTransformation & eltrans, 00923 const FlatVector<Complex> & elx, 00924 FlatVector<Complex> & ely, 00925 void * precomputed, 00926 LocalHeap & lh) const; 00927 00928 virtual void 00929 CalcFlux (const FiniteElement & fel, 00930 const ElementTransformation & eltrans, 00931 const IntegrationPoint & ip, 00932 const FlatVector<Complex> & elx, 00933 FlatVector<Complex> & flux, 00934 bool applyd, 00935 LocalHeap & lh) const; 00936 00937 00938 virtual void 00939 CalcFlux (const FiniteElement & fel, 00940 const BaseMappedIntegrationPoint & bmip, 00941 const FlatVector<Complex> & elx, 00942 FlatVector<Complex> & flux, 00943 bool applyd, 00944 LocalHeap & lh) const; 00945 00946 virtual string Name () const; 00947 00948 /* 00949 virtual const IntegrationRule & GetIntegrationRule (const FiniteElement & fel, 00950 const bool use_higher_integration_order = false) const; 00951 */ 00952 }; 00953 00954 00955 00956 00957 00958 00959 class NGS_DLL_HEADER CompoundBilinearFormIntegrator : public BilinearFormIntegrator 00960 { 00961 const BilinearFormIntegrator & bfi; 00962 int comp; 00963 public: 00964 CompoundBilinearFormIntegrator (const BilinearFormIntegrator & abfi, int acomp); 00965 00966 const BilinearFormIntegrator * GetBFI(void) const {return &bfi;} 00967 00968 virtual bool BoundaryForm () const 00969 { return bfi.BoundaryForm(); } 00970 00971 virtual int DimFlux () const 00972 { return bfi.DimFlux(); } 00973 virtual int DimElement () const 00974 { return bfi.DimElement(); } 00975 virtual int DimSpace () const 00976 { return bfi.DimSpace(); } 00977 00978 00979 virtual void 00980 CalcElementMatrix (const FiniteElement & bfel, 00981 const ElementTransformation & eltrans, 00982 FlatMatrix<double> & elmat, 00983 LocalHeap & lh) const; 00984 00985 virtual void 00986 CalcElementMatrix (const FiniteElement & bfel, 00987 const ElementTransformation & eltrans, 00988 FlatMatrix<Complex> & elmat, 00989 LocalHeap & lh) const; 00990 00991 00992 virtual void 00993 CalcLinearizedElementMatrix (const FiniteElement & fel, 00994 const ElementTransformation & eltrans, 00995 FlatVector<double> & elveclin, 00996 FlatMatrix<double> & elmat, 00997 LocalHeap & lh) const; 00998 00999 virtual void 01000 CalcLinearizedElementMatrix (const FiniteElement & fel, 01001 const ElementTransformation & eltrans, 01002 FlatVector<Complex> & elveclin, 01003 FlatMatrix<Complex> & elmat, 01004 LocalHeap & lh) const; 01005 01006 virtual void 01007 ApplyElementMatrix (const FiniteElement & bfel, 01008 const ElementTransformation & eltrans, 01009 const FlatVector<double> & elx, 01010 FlatVector<double> & ely, 01011 void * precomputed, 01012 LocalHeap & lh) const; 01013 01014 virtual void 01015 ApplyElementMatrix (const FiniteElement & bfel, 01016 const ElementTransformation & eltrans, 01017 const FlatVector<Complex> & elx, 01018 FlatVector<Complex> & ely, 01019 void * precomputed, 01020 LocalHeap & lh) const; 01021 01022 virtual void 01023 ApplyLinearizedElementMatrix (const FiniteElement & bfel, 01024 const ElementTransformation & eltrans, 01025 const FlatVector<double> & ellin, 01026 const FlatVector<double> & elx, 01027 FlatVector<double> & ely, 01028 LocalHeap & lh) const; 01029 01030 virtual void 01031 ApplyLinearizedElementMatrix (const FiniteElement & bfel, 01032 const ElementTransformation & eltrans, 01033 const FlatVector<Complex> & ellin, 01034 const FlatVector<Complex> & elx, 01035 FlatVector<Complex> & ely, 01036 LocalHeap & lh) const; 01037 01038 01039 virtual void 01040 CalcFlux (const FiniteElement & fel, 01041 const BaseMappedIntegrationPoint & bmip, 01042 const FlatVector<double> & elx, 01043 FlatVector<double> & flux, 01044 bool applyd, 01045 LocalHeap & lh) const; 01046 01047 01048 virtual void 01049 CalcFlux (const FiniteElement & fel, 01050 const BaseMappedIntegrationPoint & bmip, 01051 const FlatVector<Complex> & elx, 01052 FlatVector<Complex> & flux, 01053 bool applyd, 01054 LocalHeap & lh) const; 01055 01056 01057 virtual void 01058 CalcFlux (const FiniteElement & fel, 01059 const BaseMappedIntegrationRule & mir, 01060 const FlatVector<double> & elx, 01061 FlatMatrix<double> & flux, 01062 bool applyd, 01063 LocalHeap & lh) const; 01064 01065 01066 virtual void 01067 CalcFlux (const FiniteElement & fel, 01068 const BaseMappedIntegrationRule & mir, 01069 const FlatVector<Complex> & elx, 01070 FlatMatrix<Complex> & flux, 01071 bool applyd, 01072 LocalHeap & lh) const; 01073 01074 01075 virtual void 01076 ApplyBTrans (const FiniteElement & fel, 01077 const BaseMappedIntegrationPoint & bmip, 01078 const FlatVector<double> & elx, 01079 FlatVector<double> & ely, 01080 LocalHeap & lh) const; 01081 01082 virtual void 01083 ApplyBTrans (const FiniteElement & fel, 01084 const BaseMappedIntegrationPoint & bmip, 01085 const FlatVector<Complex> & elx, 01086 FlatVector<Complex> & ely, 01087 LocalHeap & lh) const; 01088 01089 virtual string Name () const; 01090 }; 01091 01092 01093 01094 01095 01096 01097 01098 01099 01103 class NGS_DLL_HEADER LinearFormIntegrator : public Integrator 01104 01105 { 01106 public: 01108 LinearFormIntegrator () { ; } 01110 virtual ~LinearFormIntegrator () { ; } 01111 01112 01116 virtual void 01117 CalcElementVector (const FiniteElement & fel, 01118 const ElementTransformation & eltrans, 01119 FlatVector<double> & elvec, 01120 LocalHeap & lh) const 01121 { 01122 AssembleElementVector (fel, eltrans, elvec, lh); 01123 } 01124 01125 01126 virtual void 01127 CalcElementVector (const FiniteElement & fel, 01128 const ElementTransformation & eltrans, 01129 FlatVector<Complex> & elvec, 01130 LocalHeap & lh) const; 01131 01132 01133 // old version: 01134 virtual void 01135 AssembleElementVector (const FiniteElement & fel, 01136 const ElementTransformation & eltrans, 01137 FlatVector<double> & elvec, 01138 LocalHeap & lh) const; 01139 01140 virtual void 01141 AssembleElementVector (const FiniteElement & fel, 01142 const ElementTransformation & eltrans, 01143 FlatVector<Complex> & elvec, 01144 LocalHeap & lh) const; 01145 01146 virtual void 01147 AssembleElementVectorIndependent (const FiniteElement & gfel, 01148 const BaseMappedIntegrationPoint & s_mip, 01149 const BaseMappedIntegrationPoint & g_mip, 01150 FlatVector<double> & elvec, 01151 LocalHeap & lh, 01152 const bool curveint = false) const 01153 { 01154 cerr << "AssembleElementVectorIndependent called for base-class!" << endl; 01155 exit(10); 01156 } 01157 01158 virtual void 01159 AssembleElementVectorIndependent (const FiniteElement & gfel, 01160 const BaseMappedIntegrationPoint & s_mip, 01161 const BaseMappedIntegrationPoint & g_mip, 01162 FlatVector<Complex> & elvec, 01163 LocalHeap & lh, 01164 const bool curveint = false) const 01165 { 01166 FlatVector<double> rvec(elvec.Size(), lh); 01167 AssembleElementVectorIndependent (gfel, s_mip, g_mip, rvec, lh,curveint); 01168 elvec = rvec; 01169 } 01170 01171 01172 }; 01173 01174 01175 01176 01177 01178 01179 01180 01181 01182 01183 01184 01185 01186 01187 01188 class NGS_DLL_HEADER FacetLinearFormIntegrator : public LinearFormIntegrator 01189 { 01190 public: 01191 01192 FacetLinearFormIntegrator(Array<CoefficientFunction*> & coeffs) 01193 : LinearFormIntegrator() { ; } 01194 01195 ~FacetLinearFormIntegrator() { ; } 01196 01197 virtual bool BoundaryForm () const 01198 { return 1; } 01199 01200 virtual bool SkeletonForm () const 01201 { return 1; } 01202 01203 virtual void 01204 CalcElementVector (const FiniteElement & bfel, 01205 const ElementTransformation & eltrans, 01206 FlatVector<double> & elvec, 01207 LocalHeap & lh) const{ 01208 throw Exception ("FacetLinearFormIntegrator can not assemble volumetric element matrices!"); 01209 } 01210 01211 01212 virtual void 01213 CalcFacetVector (const FiniteElement & volumefel, int LocalFacetNr, 01214 const ElementTransformation & eltrans, FlatArray<int> & ElVertices, 01215 const ElementTransformation & seltrans, 01216 FlatVector<double> & elvec, 01217 LocalHeap & lh) const{ 01218 throw Exception ("FacetLinearFormIntegrator::CalcFacetVector not implemented!"); 01219 } 01220 01221 virtual void 01222 CalcFacetVector (const FiniteElement & volumefel, int LocalFacetNr, 01223 const ElementTransformation & eltrans, FlatArray<int> & ElVertices, 01224 const ElementTransformation & seltrans, 01225 FlatVector<Complex> & elvec, 01226 LocalHeap & lh) const { 01227 throw Exception ("FacetLinearFormIntegrator::CalcFacetVector<complex> not implemented!"); 01228 } 01229 01230 }; 01231 01232 01233 01234 01235 01236 class NGS_DLL_HEADER BlockLinearFormIntegrator : public LinearFormIntegrator 01237 { 01238 const LinearFormIntegrator & lfi; 01239 int dim; 01240 int comp; 01241 public: 01242 BlockLinearFormIntegrator (const LinearFormIntegrator & alfi, int adim, int acomp); 01243 01244 virtual bool BoundaryForm () const 01245 { return lfi.BoundaryForm(); } 01246 01247 01248 virtual void 01249 CalcElementVector (const FiniteElement & bfel, 01250 const ElementTransformation & eltrans, 01251 FlatVector<double> & elvec, 01252 LocalHeap & lh) const; 01253 }; 01254 01255 01256 01257 01258 01259 class NGS_DLL_HEADER ComplexLinearFormIntegrator : public LinearFormIntegrator 01260 { 01261 const LinearFormIntegrator & lfi; 01262 Complex factor; 01263 public: 01264 ComplexLinearFormIntegrator (const LinearFormIntegrator & alfi, 01265 Complex afactor) 01266 : lfi(alfi), factor(afactor) 01267 { ; } 01268 01269 virtual bool BoundaryForm () const { return lfi.BoundaryForm(); } 01270 01271 virtual void CheckElement (const FiniteElement & el) const { lfi.CheckElement(el); } 01272 01273 01274 virtual void 01275 CalcElementVector (const FiniteElement & fel, 01276 const ElementTransformation & eltrans, 01277 FlatVector<double> & elvec, 01278 LocalHeap & lh) const 01279 { 01280 throw Exception ("ComplexLinearFormIntegrator: cannot assemble double vector"); 01281 } 01282 01283 virtual void 01284 CalcElementVector (const FiniteElement & fel, 01285 const ElementTransformation & eltrans, 01286 FlatVector<Complex> & elvec, 01287 LocalHeap & lh) const 01288 { 01289 FlatVector<Complex> rvec(elvec.Size(), lh); 01290 lfi.CalcElementVector (fel, eltrans, rvec, lh); 01291 elvec = factor * rvec; 01292 } 01293 01294 01295 virtual void 01296 AssembleElementVectorIndependent (const FiniteElement & gfel, 01297 const BaseMappedIntegrationPoint & s_mip, 01298 const BaseMappedIntegrationPoint & g_mip, 01299 FlatVector<double> & elvec, 01300 LocalHeap & lh, 01301 const bool curveint = false) const 01302 { 01303 throw Exception ("ComplexLinearFormIntegrator: cannot assemble double vector"); 01304 } 01305 01306 01307 virtual void 01308 AssembleElementVectorIndependent (const FiniteElement & gfel, 01309 const BaseMappedIntegrationPoint & s_mip, 01310 const BaseMappedIntegrationPoint & g_mip, 01311 FlatVector<Complex> & elvec, 01312 LocalHeap & lh, 01313 const bool curveint = false) const 01314 { 01315 FlatVector<double> rvec; 01316 01317 lfi.AssembleElementVectorIndependent (gfel, s_mip, g_mip, 01318 rvec, lh, curveint); 01319 elvec.AssignMemory (rvec.Size(), lh); 01320 elvec = factor * rvec; 01321 } 01322 01323 01324 01325 01326 virtual string Name () const 01327 { 01328 return string ("ComplexIntegrator (") + lfi.Name() + ")"; 01329 } 01330 01331 }; 01332 01333 01334 01335 01336 class NGS_DLL_HEADER CompoundLinearFormIntegrator : public LinearFormIntegrator 01337 { 01338 const LinearFormIntegrator & lfi; 01339 int comp; 01340 public: 01341 CompoundLinearFormIntegrator (const LinearFormIntegrator & alfi, int acomp) 01342 : lfi(alfi), comp(acomp) { ; } 01343 01344 virtual bool BoundaryForm () const 01345 { return lfi.BoundaryForm(); } 01346 01347 01348 virtual void 01349 CalcElementVector (const FiniteElement & bfel, 01350 const ElementTransformation & eltrans, 01351 FlatVector<double> & elvec, 01352 LocalHeap & lh) const; 01353 01354 virtual void 01355 CalcElementVector (const FiniteElement & bfel, 01356 const ElementTransformation & eltrans, 01357 FlatVector<Complex> & elvec, 01358 LocalHeap & lh) const; 01359 01360 virtual void 01361 AssembleElementVectorIndependent (const FiniteElement & gfel, 01362 const BaseMappedIntegrationPoint & s_mip, 01363 const BaseMappedIntegrationPoint & g_mip, 01364 FlatVector<double> & elvec, 01365 LocalHeap & lh, 01366 const bool curveint = false) const; 01367 01368 virtual void 01369 AssembleElementVectorIndependent (const FiniteElement & gfel, 01370 const BaseMappedIntegrationPoint & s_mip, 01371 const BaseMappedIntegrationPoint & g_mip, 01372 FlatVector<Complex> & elvec, 01373 LocalHeap & lh, 01374 const bool curveint = false) const; 01375 01376 virtual string Name () const 01377 { 01378 return string ("CompoundIntegrator (") + lfi.Name() + ")"; 01379 } 01380 }; 01381 01382 01383 01384 01385 01386 01387 01388 01389 01391 class NGS_DLL_HEADER Integrators 01392 { 01393 public: 01394 01396 class IntegratorInfo 01397 { 01398 public: 01399 string name; 01400 int spacedim; 01401 int numcoeffs; 01402 Integrator* (*creator)(Array<CoefficientFunction*> &); 01403 01404 IntegratorInfo (const string & aname, 01405 int aspacedim, 01406 int anumcoffs, 01407 Integrator* (*acreator)(Array<CoefficientFunction*>&)); 01408 }; 01409 01410 01411 Array<IntegratorInfo*> bfis; 01412 Array<IntegratorInfo*> lfis; 01413 01414 public: 01416 Integrators(); 01418 ~Integrators(); 01420 void AddBFIntegrator (const string & aname, int aspacedim, int anumcoeffs, 01421 Integrator* (*acreator)(Array<CoefficientFunction*>&)); 01423 void AddLFIntegrator (const string & aname, int aspacedim, int anumcoeffs, 01424 Integrator* (*acreator)(Array<CoefficientFunction*>&)); 01425 01427 const Array<IntegratorInfo*> & GetBFIs() const { return bfis; } 01429 const IntegratorInfo * GetBFI(const string & name, int dim) const; 01431 BilinearFormIntegrator * CreateBFI(const string & name, int dim, 01432 Array<CoefficientFunction*> & coeffs) const; 01434 BilinearFormIntegrator * CreateBFI(const string & name, int dim, 01435 CoefficientFunction * coef) const; 01436 01438 const Array<IntegratorInfo*> & GetLFIs() const { return lfis; } 01440 const IntegratorInfo * GetLFI(const string & name, int dim) const; 01442 LinearFormIntegrator * CreateLFI(const string & name, int dim, 01443 Array<CoefficientFunction*> & coeffs) const; 01444 01445 LinearFormIntegrator * CreateLFI(const string & name, int dim, 01446 CoefficientFunction * coef) const; 01447 01449 void Print (ostream & ost) const; 01450 }; 01451 01453 extern NGS_DLL_HEADER Integrators & GetIntegrators (); 01454 01455 01456 01457 01458 01459 01460 template <typename BFI> 01461 class RegisterBilinearFormIntegrator 01462 { 01463 public: 01464 RegisterBilinearFormIntegrator (string label, int dim, int numcoeffs) 01465 { 01466 GetIntegrators().AddBFIntegrator (label, dim, numcoeffs, Create); 01467 // cout << "register bf-integrator '" << label << "'" << endl; 01468 } 01469 01470 static Integrator * Create (Array<CoefficientFunction*> & coefs) 01471 { 01472 return new BFI (coefs); 01473 } 01474 }; 01475 01476 01477 01478 template <typename LFI> 01479 class RegisterLinearFormIntegrator 01480 { 01481 public: 01482 RegisterLinearFormIntegrator (string label, int dim, int numcoeffs) 01483 { 01484 GetIntegrators().AddLFIntegrator (label, dim, numcoeffs, Create); 01485 // cout << "register lf-integrator '" << label << "'" << endl; 01486 } 01487 01488 static Integrator * Create (Array<CoefficientFunction*> & coefs) 01489 { 01490 return new LFI (coefs); 01491 } 01492 }; 01493 01494 01495 01496 01497 } 01498 01499 #endif