MFEM  v3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
bilinearform.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_BILINEARFORM
13 #define MFEM_BILINEARFORM
14 
15 #include "../config/config.hpp"
16 #include "../linalg/linalg.hpp"
17 #include "fespace.hpp"
18 #include "gridfunc.hpp"
19 #include "linearform.hpp"
20 #include "bilininteg.hpp"
21 #include "staticcond.hpp"
22 #include "hybridization.hpp"
23 
24 namespace mfem
25 {
26 
29 class BilinearForm : public Matrix
30 {
31 protected:
34 
37 
40 
43  long sequence;
44 
46 
49 
52 
55 
58 
61 
63 
66 
68  // Allocate appropriate SparseMatrix and assign it to mat
69  void AllocMat();
70 
71  void ConformingAssemble();
72 
73  // may be used in the construction of derived classes
75  {
76  fes = NULL; sequence = -1;
77  mat = mat_e = NULL; extern_bfs = 0; element_matrices = NULL;
78  static_cond = NULL; hybridization = NULL;
80  }
81 
82 public:
85 
86  BilinearForm(FiniteElementSpace *f, BilinearForm *bf, int ps = 0);
87 
89  int Size() const { return height; }
90 
96 
99  bool StaticCondensationIsEnabled() const { return static_cond; }
100 
103  { return static_cond ? static_cond->GetTraceFESpace() : NULL; }
104 
108  void EnableHybridization(FiniteElementSpace *constr_space,
109  BilinearFormIntegrator *constr_integ,
110  const Array<int> &ess_tdof_list);
111 
115  void UsePrecomputedSparsity(int ps = 1) { precompute_sparsity = ps; }
116 
120  void AllocateMatrix() { if (mat == NULL) { AllocMat(); } }
121 
123 
125 
127 
129 
130  const double &operator()(int i, int j) { return (*mat)(i,j); }
131 
133  virtual double &Elem(int i, int j);
134 
136  virtual const double &Elem(int i, int j) const;
137 
139  virtual void Mult(const Vector &x, Vector &y) const;
140 
141  void FullMult(const Vector &x, Vector &y) const
142  { mat->Mult(x, y); mat_e->AddMult(x, y); }
143 
144  virtual void AddMult(const Vector &x, Vector &y, const double a = 1.0) const
145  { mat -> AddMult (x, y, a); }
146 
147  void FullAddMult(const Vector &x, Vector &y) const
148  { mat->AddMult(x, y); mat_e->AddMult(x, y); }
149 
150  double InnerProduct(const Vector &x, const Vector &y) const
151  { return mat->InnerProduct (x, y); }
152 
154  virtual MatrixInverse *Inverse() const;
155 
157  virtual void Finalize(int skip_zeros = 1);
158 
160  const SparseMatrix &SpMat() const
161  {
162  MFEM_VERIFY(mat, "mat is NULL and can't be dereferenced");
163  return *mat;
164  }
166  {
167  MFEM_VERIFY(mat, "mat is NULL and can't be dereferenced");
168  return *mat;
169  }
170  SparseMatrix *LoseMat() { SparseMatrix *tmp = mat; mat = NULL; return tmp; }
171 
173  const SparseMatrix &SpMatElim() const
174  {
175  MFEM_VERIFY(mat_e, "mat_e is NULL and can't be dereferenced");
176  return *mat_e;
177  }
179  {
180  MFEM_VERIFY(mat_e, "mat_e is NULL and can't be dereferenced");
181  return *mat_e;
182  }
183 
186 
189 
192 
195 
196  void operator=(const double a)
197  {
198  if (mat != NULL) { *mat = a; }
199  if (mat_e != NULL) { *mat_e = a; }
200  }
201 
203  void Assemble(int skip_zeros = 1);
204 
227  void FormLinearSystem(Array<int> &ess_tdof_list, Vector &x, Vector &b,
228  SparseMatrix &A, Vector &X, Vector &B,
229  int copy_interior = 0);
230 
234  void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
235 
237  void ComputeElementMatrices();
238 
241  { delete element_matrices; element_matrices = NULL; }
242 
243  void ComputeElementMatrix(int i, DenseMatrix &elmat);
244  void AssembleElementMatrix(int i, const DenseMatrix &elmat,
245  Array<int> &vdofs, int skip_zeros = 1);
246  void AssembleBdrElementMatrix(int i, const DenseMatrix &elmat,
247  Array<int> &vdofs, int skip_zeros = 1);
248 
253  void EliminateEssentialBC(Array<int> &bdr_attr_is_ess,
254  Vector &sol, Vector &rhs, int d = 0);
255 
256  void EliminateEssentialBC(Array<int> &bdr_attr_is_ess, int d = 0);
258  void EliminateEssentialBCDiag(Array<int> &bdr_attr_is_ess, double value);
259 
261  void EliminateVDofs(Array<int> &vdofs, Vector &sol, Vector &rhs, int d = 0);
262 
267  void EliminateVDofs(Array<int> &vdofs, int d = 0);
268 
271  void EliminateEssentialBCFromDofs(Array<int> &ess_dofs, Vector &sol,
272  Vector &rhs, int d = 0);
273 
276  void EliminateEssentialBCFromDofs(Array<int> &ess_dofs, int d = 0);
278  void EliminateEssentialBCFromDofsDiag(Array<int> &ess_dofs, double value);
279 
282  void EliminateVDofsInRHS(Array<int> &vdofs, const Vector &x, Vector &b);
283 
284  double FullInnerProduct(const Vector &x, const Vector &y) const
285  { return mat->InnerProduct(x, y) + mat_e->InnerProduct(x, y); }
286 
287  virtual void Update(FiniteElementSpace *nfes = NULL);
288 
291 
293  virtual ~BilinearForm();
294 };
295 
311 class MixedBilinearForm : public Matrix
312 {
313 protected:
315 
317 
320  Array<BilinearFormIntegrator*> skt; // trace face integrators
321 
322 public:
324  FiniteElementSpace *te_fes);
325 
326  virtual double& Elem (int i, int j);
327 
328  virtual const double& Elem (int i, int j) const;
329 
330  virtual void Mult (const Vector & x, Vector & y) const;
331 
332  virtual void AddMult (const Vector & x, Vector & y,
333  const double a = 1.0) const;
334 
335  virtual void AddMultTranspose (const Vector & x, Vector & y,
336  const double a = 1.0) const;
337 
338  virtual void MultTranspose (const Vector & x, Vector & y) const
339  { y = 0.0; AddMultTranspose (x, y); }
340 
341  virtual MatrixInverse * Inverse() const;
342 
343  virtual void Finalize (int skip_zeros = 1);
344 
348  void GetBlocks(Array2D<SparseMatrix *> &blocks) const;
349 
350  const SparseMatrix &SpMat() const { return *mat; }
351  SparseMatrix &SpMat() { return *mat; }
352  SparseMatrix *LoseMat() { SparseMatrix *tmp = mat; mat = NULL; return tmp; }
353 
355 
357 
362 
364 
366 
368 
369  void operator= (const double a) { *mat = a; }
370 
371  void Assemble (int skip_zeros = 1);
372 
378  void ConformingAssemble();
379 
380  void EliminateTrialDofs(Array<int> &bdr_attr_is_ess,
381  Vector &sol, Vector &rhs);
382 
384  Vector &sol, Vector &rhs);
385 
386  virtual void EliminateTestDofs(Array<int> &bdr_attr_is_ess);
387 
388  void Update();
389 
390  virtual ~MixedBilinearForm();
391 };
392 
393 
425 {
426 public:
428  FiniteElementSpace *range_fes)
429  : MixedBilinearForm(domain_fes, range_fes) { }
430 
432  { AddDomainIntegrator(di); }
433 
435  { AddTraceFaceIntegrator(di); }
436 
438 
439  virtual void Assemble(int skip_zeros = 1);
440 };
441 
442 }
443 
444 #endif
void EliminateEssentialBCFromDofs(Array< int > &ess_dofs, Vector &sol, Vector &rhs, int d=0)
void FullAddMult(const Vector &x, Vector &y) const
SparseMatrix & SpMat()
Array< BilinearFormIntegrator * > * GetBBFI()
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
virtual void AddMult(const Vector &x, Vector &y, const double a=1.0) const
Array< BilinearFormIntegrator * > skt
virtual MatrixInverse * Inverse() const
Returns a pointer to (approximation) of the matrix inverse.
Array< BilinearFormIntegrator * > fbfi
Set of interior face Integrators to be applied.
void Assemble(int skip_zeros=1)
Assembles the form i.e. sums over all domain/bdr integrators.
Array< BilinearFormIntegrator * > * GetDBFI()
Array< BilinearFormIntegrator * > dbfi
Set of Domain Integrators to be applied.
virtual void AddMult(const Vector &x, Vector &y, const double a=1.0) const
Array< BilinearFormIntegrator * > * GetBFBFI()
void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
Data type dense matrix using column-major storage.
Definition: densemat.hpp:22
void AddMult(const Vector &x, Vector &y, const double a=1.0) const
y += A * x (default) or y += a * A * x
Definition: sparsemat.cpp:445
Abstract data type for matrix inverse.
Definition: matrix.hpp:58
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
Array< BilinearFormIntegrator * > * GetFBFI()
MixedBilinearForm(FiniteElementSpace *tr_fes, FiniteElementSpace *te_fes)
void FormLinearSystem(Array< int > &ess_tdof_list, Vector &x, Vector &b, SparseMatrix &A, Vector &X, Vector &B, int copy_interior=0)
double FullInnerProduct(const Vector &x, const Vector &y) const
void AddDomainInterpolator(DiscreteInterpolator *di)
void EnableHybridization(FiniteElementSpace *constr_space, BilinearFormIntegrator *constr_integ, const Array< int > &ess_tdof_list)
void EliminateEssentialBCFromTrialDofs(Array< int > &marked_vdofs, Vector &sol, Vector &rhs)
Array< BilinearFormIntegrator * > bdr
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
SparseMatrix & SpMat()
void EliminateVDofs(Array< int > &vdofs, Vector &sol, Vector &rhs, int d=0)
Eliminate the given vdofs. NOTE: here, vdofs is a list of DOFs.
DiscreteLinearOperator(FiniteElementSpace *domain_fes, FiniteElementSpace *range_fes)
virtual void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const
Data type sparse matrix.
Definition: sparsemat.hpp:38
StaticCondensation * static_cond
SparseMatrix * mat
Sparse matrix to be associated with the form.
int Size() const
Get the size of the BilinearForm as a square matrix.
Abstract data type matrix.
Definition: matrix.hpp:27
void Assemble(int skip_zeros=1)
void FullMult(const Vector &x, Vector &y) const
virtual void Update(FiniteElementSpace *nfes=NULL)
void EliminateTrialDofs(Array< int > &bdr_attr_is_ess, Vector &sol, Vector &rhs)
FiniteElementSpace * GetFES()
Return the FE space associated with the BilinearForm.
void AddTraceFaceIntegrator(BilinearFormIntegrator *bfi)
void GetBlocks(Array2D< SparseMatrix * > &blocks) const
Array< BilinearFormIntegrator * > * GetTFBFI()
Array< BilinearFormIntegrator * > * GetBBFI()
void ComputeElementMatrix(int i, DenseMatrix &elmat)
void EliminateEssentialBCDiag(Array< int > &bdr_attr_is_ess, double value)
Perform elimination and set the diagonal entry to the given value.
SparseMatrix & SpMatElim()
SparseMatrix * mat_e
Matrix used to eliminate b.c.
bool StaticCondensationIsEnabled() const
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:22
void AssembleElementMatrix(int i, const DenseMatrix &elmat, Array< int > &vdofs, int skip_zeros=1)
void UsePrecomputedSparsity(int ps=1)
SparseMatrix * LoseMat()
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
Definition: sparsemat.cpp:439
Array< BilinearFormIntegrator * > bfbfi
Set of boundary face Integrators to be applied.
Array< BilinearFormIntegrator * > bbfi
Set of Boundary Integrators to be applied.
DenseTensor * element_matrices
SparseMatrix * LoseMat()
void operator=(const double a)
void AddTraceFaceInterpolator(DiscreteInterpolator *di)
virtual double & Elem(int i, int j)
Returns reference to a_{ij}.
Array< BilinearFormIntegrator * > dom
double InnerProduct(const Vector &x, const Vector &y) const
const SparseMatrix & SpMat() const
FiniteElementSpace * test_fes
const double & operator()(int i, int j)
FiniteElementSpace * GetTraceFESpace()
Return a pointer to the reduced/trace FE space.
Definition: staticcond.hpp:107
virtual void EliminateTestDofs(Array< int > &bdr_attr_is_ess)
void EliminateVDofsInRHS(Array< int > &vdofs, const Vector &x, Vector &b)
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
const SparseMatrix & SpMatElim() const
Returns a reference to the sparse matrix of eliminated b.c.
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi)
Adds new Boundary Integrator.
void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new interior Face Integrator.
FiniteElementSpace * fes
FE space on which the form lives.
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator.
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator.
void ComputeElementMatrices()
Compute and store internally all element matrices.
FiniteElementSpace * SCFESpace() const
Return the trace FE space associated with static condensation.
void AssembleBdrElementMatrix(int i, const DenseMatrix &elmat, Array< int > &vdofs, int skip_zeros=1)
DenseMatrix elemmat
FiniteElementSpace * trial_fes
virtual double & Elem(int i, int j)
Returns reference to a_{ij}.
virtual ~BilinearForm()
Destroys bilinear form.
Vector data type.
Definition: vector.hpp:33
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi)
virtual MatrixInverse * Inverse() const
Returns a pointer to (an approximation) of the matrix inverse.
Hybridization * hybridization
void operator=(const double a)
Array< BilinearFormIntegrator * > * GetDBFI()
const SparseMatrix & SpMat() const
Returns a reference to the sparse matrix.
Array< int > vdofs
void FreeElementMatrices()
Free the memory used by the element matrices.
void EliminateEssentialBCFromDofsDiag(Array< int > &ess_dofs, double value)
Perform elimination and set the diagonal entry to the given value.
virtual void Assemble(int skip_zeros=1)
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:569
Array< BilinearFormIntegrator * > * GetDI()
double InnerProduct(const Vector &x, const Vector &y) const
Compute y^t A x.
Definition: sparsemat.cpp:641
void EliminateEssentialBC(Array< int > &bdr_attr_is_ess, Vector &sol, Vector &rhs, int d=0)
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
void EnableStaticCondensation()
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new boundary Face Integrator.