MFEM v2.0
bilinearform.hpp
Go to the documentation of this file.
00001 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
00002 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
00003 // reserved. See file COPYRIGHT for details.
00004 //
00005 // This file is part of the MFEM library. For more information and source code
00006 // availability see http://mfem.googlecode.com.
00007 //
00008 // MFEM is free software; you can redistribute it and/or modify it under the
00009 // terms of the GNU Lesser General Public License (as published by the Free
00010 // Software Foundation) version 2.1 dated February 1999.
00011 
00012 #ifndef MFEM_BILINEARFORM
00013 #define MFEM_BILINEARFORM
00014 
00015 
00018 class BilinearForm : public Matrix
00019 {
00020 protected:
00022    SparseMatrix *mat;
00023 
00024    // Matrix used to eliminate b.c.
00025    SparseMatrix *mat_e;
00026 
00028    FiniteElementSpace *fes;
00029 
00030    int extern_bfs;
00031 
00033    Array<BilinearFormIntegrator*> dbfi;
00034 
00036    Array<BilinearFormIntegrator*> bbfi;
00037 
00039    Array<BilinearFormIntegrator*> fbfi;
00040 
00042    Array<BilinearFormIntegrator*> bfbfi;
00043 
00044    DenseMatrix elemmat;
00045    Array<int>  vdofs;
00046 
00047    DenseTensor *element_matrices;
00048 
00049    // may be used in the construction of derived classes
00050    BilinearForm() : Matrix (0)
00051    { fes = NULL; mat = mat_e = NULL; extern_bfs = 0; element_matrices = NULL; }
00052 
00053 public:
00055    BilinearForm (FiniteElementSpace * f);
00056 
00057    BilinearForm (FiniteElementSpace * f, BilinearForm * bf);
00058 
00059    Array<BilinearFormIntegrator*> *GetDBFI() { return &dbfi; }
00060 
00061    Array<BilinearFormIntegrator*> *GetBBFI() { return &bbfi; }
00062 
00063    Array<BilinearFormIntegrator*> *GetFBFI() { return &fbfi; }
00064 
00065    Array<BilinearFormIntegrator*> *GetBFBFI() { return &bfbfi; }
00066 
00067    const double &operator()(int i, int j) { return (*mat)(i,j); }
00068 
00070    virtual double &Elem(int i, int j);
00071 
00073    virtual const double &Elem(int i, int j) const;
00074 
00076    virtual void Mult(const Vector &x, Vector &y) const;
00077 
00078    void FullMult(const Vector &x, Vector &y) const
00079    { mat->Mult(x, y); mat_e->AddMult(x, y); }
00080 
00081    virtual void AddMult(const Vector &x, Vector &y, const double a = 1.0) const
00082    { mat -> AddMult (x, y, a); }
00083 
00084    void FullAddMult(const Vector &x, Vector &y) const
00085    { mat->AddMult(x, y); mat_e->AddMult(x, y); }
00086 
00087    double InnerProduct(const Vector &x, const Vector &y) const
00088    { return mat->InnerProduct (x, y); }
00089 
00091    virtual MatrixInverse *Inverse() const;
00092 
00094    virtual void Finalize(int skip_zeros = 1);
00095 
00097    const SparseMatrix &SpMat() const { return *mat; }
00098    SparseMatrix &SpMat() { return *mat; }
00099 
00101    void AddDomainIntegrator(BilinearFormIntegrator *bfi);
00102 
00104    void AddBoundaryIntegrator(BilinearFormIntegrator *bfi);
00105 
00107    void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi);
00108 
00110    void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi);
00111 
00112    void operator=(const double a)
00113    { if (mat != NULL) *mat = a; if (mat_e != NULL) *mat_e = a; }
00114 
00116    void Assemble(int skip_zeros = 1);
00117 
00119    void ComputeElementMatrices();
00120 
00122    void FreeElementMatrices()
00123    { delete element_matrices; element_matrices = NULL; }
00124 
00125    void ComputeElementMatrix(int i, DenseMatrix &elmat);
00126    void AssembleElementMatrix(int i, const DenseMatrix &elmat,
00127                               Array<int> &vdofs, int skip_zeros = 1);
00128 
00131    void EliminateEssentialBC(Array<int> &bdr_attr_is_ess,
00132                              Vector &sol, Vector &rhs, int d = 0);
00133 
00135    void EliminateVDofs(Array<int> &vdofs, Vector &sol, Vector &rhs, int d = 0);
00136 
00139    void EliminateVDofs(Array<int> &vdofs, int d = 0);
00140 
00143    void EliminateVDofsInRHS(Array<int> &vdofs, const Vector &x, Vector &b);
00144 
00145    double FullInnerProduct(const Vector &x, const Vector &y) const
00146    { return mat->InnerProduct(x, y) + mat_e->InnerProduct(x, y); }
00147 
00148    void EliminateEssentialBC(Array<int> &bdr_attr_is_ess, int d = 0);
00149 
00152    void EliminateEssentialBCFromDofs(Array<int> &ess_dofs, Vector &sol,
00153                                      Vector &rhs, int d = 0);
00154 
00157    void EliminateEssentialBCFromDofs(Array<int> &ess_dofs, int d = 0);
00158 
00159    void Update(FiniteElementSpace *nfes = NULL);
00160 
00161    FiniteElementSpace *GetFES() { return fes; }
00162 
00164    virtual ~BilinearForm();
00165 };
00166 
00182 class MixedBilinearForm : public Matrix
00183 {
00184 protected:
00185    SparseMatrix *mat;
00186 
00187    FiniteElementSpace *trial_fes, *test_fes;
00188 
00189    Array<BilinearFormIntegrator*> dom;
00190    Array<BilinearFormIntegrator*> bdr;
00191 
00192    int width;
00193 
00194 public:
00195    MixedBilinearForm (FiniteElementSpace *tr_fes,
00196                       FiniteElementSpace *te_fes);
00197 
00198    int Width() const { return width; };
00199 
00200    virtual double& Elem (int i, int j);
00201 
00202    virtual const double& Elem (int i, int j) const;
00203 
00204    virtual void Mult (const Vector & x, Vector & y) const;
00205 
00206    virtual void AddMult (const Vector & x, Vector & y,
00207                          const double a = 1.0) const;
00208 
00209    virtual void AddMultTranspose (const Vector & x, Vector & y,
00210                                   const double a = 1.0) const;
00211 
00212    virtual void MultTranspose (const Vector & x, Vector & y) const
00213    { y = 0.0; AddMultTranspose (x, y); };
00214 
00215    virtual MatrixInverse * Inverse() const;
00216 
00217    virtual void Finalize (int skip_zeros = 1);
00218 
00222    void GetBlocks(Array2D<SparseMatrix *> &blocks) const;
00223 
00224    const SparseMatrix &SpMat() const { return *mat; };
00225    SparseMatrix &SpMat() { return *mat; };
00226 
00227    void AddDomainIntegrator (BilinearFormIntegrator * bfi);
00228 
00229    void AddBoundaryIntegrator (BilinearFormIntegrator * bfi);
00230 
00231    void operator= (const double a) { *mat = a; }
00232 
00233    virtual void Assemble (int skip_zeros = 1);
00234 
00235    void EliminateTrialDofs (Array<int> &bdr_attr_is_ess,
00236                             Vector &sol, Vector &rhs);
00237 
00238    virtual void EliminateTestDofs (Array<int> &bdr_attr_is_ess);
00239 
00240    void Update();
00241 
00242    virtual ~MixedBilinearForm();
00243 };
00244 
00245 
00276 class DiscreteLinearOperator : public MixedBilinearForm
00277 {
00278 public:
00279    DiscreteLinearOperator(FiniteElementSpace *domain_fes,
00280                           FiniteElementSpace *range_fes)
00281       : MixedBilinearForm(domain_fes, range_fes) { }
00282 
00283    void AddDomainInterpolator(DiscreteInterpolator *di)
00284    { AddDomainIntegrator(di); }
00285 
00286    virtual void Assemble(int skip_zeros = 1);
00287 };
00288 
00289 #endif
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Defines