MFEM  v3.3
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 
27 /** Class for bilinear form - "Matrix" with associated FE space and
28  BLFIntegrators. */
29 class BilinearForm : public Matrix
30 {
31 protected:
32  /// Sparse matrix to be associated with the form.
34 
35  /// Matrix used to eliminate b.c.
37 
38  /// FE space on which the form lives.
40 
41  /// Indicates the Mesh::sequence corresponding to the current state of the
42  /// BilinearForm.
43  long sequence;
44 
46 
47  /// Set of Domain Integrators to be applied.
49 
50  /// Set of Boundary Integrators to be applied.
52 
53  /// Set of interior face Integrators to be applied.
55 
56  /// Set of boundary face Integrators to be applied.
59 
62 
64 
67 
69  // Allocate appropriate SparseMatrix and assign it to mat
70  void AllocMat();
71 
72  void ConformingAssemble();
73 
74  // may be used in the construction of derived classes
76  {
77  fes = NULL; sequence = -1;
78  mat = mat_e = NULL; extern_bfs = 0; element_matrices = NULL;
79  static_cond = NULL; hybridization = NULL;
81  }
82 
83 public:
84  /// Creates bilinear form associated with FE space *f.
86 
87  BilinearForm(FiniteElementSpace *f, BilinearForm *bf, int ps = 0);
88 
89  /// Get the size of the BilinearForm as a square matrix.
90  int Size() const { return height; }
91 
92  /** Enable the use of static condensation. For details see the description
93  for class StaticCondensation in fem/staticcond.hpp This method should be
94  called before assembly. If the number of unknowns after static
95  condensation is not reduced, it is not enabled. */
97 
98  /** Check if static condensation was actually enabled by a previous call to
99  EnableStaticCondensation. */
100  bool StaticCondensationIsEnabled() const { return static_cond; }
101 
102  /// Return the trace FE space associated with static condensation.
104  { return static_cond ? static_cond->GetTraceFESpace() : NULL; }
105 
106  /** Enable hybridization; for details see the description for class
107  Hybridization in fem/hybridization.hpp. This method should be called
108  before assembly. */
109  void EnableHybridization(FiniteElementSpace *constr_space,
110  BilinearFormIntegrator *constr_integ,
111  const Array<int> &ess_tdof_list);
112 
113  /** For scalar FE spaces, precompute the sparsity pattern of the matrix
114  (assuming dense element matrices) based on the types of integrators
115  present in the bilinear form. */
116  void UsePrecomputedSparsity(int ps = 1) { precompute_sparsity = ps; }
117 
118  /** @brief Use the given CSR sparsity pattern to allocate the internal
119  SparseMatrix.
120 
121  - The @a I and @a J arrays must define a square graph with size equal to
122  GetVSize() of the associated FiniteElementSpace.
123  - This method should be called after enabling static condensation or
124  hybridization, if used.
125  - In the case of static condensation, @a I and @a J are not used.
126  - The ownership of the arrays @a I and @a J remains with the caller. */
127  void UseSparsity(int *I, int *J, bool isSorted);
128 
129  /// Use the sparsity of @a A to allocate the internal SparseMatrix.
130  void UseSparsity(SparseMatrix &A);
131 
132  /** Pre-allocate the internal SparseMatrix before assembly. If the flag
133  'precompute sparsity' is set, the matrix is allocated in CSR format (i.e.
134  finalized) and the entries are initialized with zeros. */
135  void AllocateMatrix() { if (mat == NULL) { AllocMat(); } }
136 
138 
140 
142 
144 
145  const double &operator()(int i, int j) { return (*mat)(i,j); }
146 
147  /// Returns reference to a_{ij}.
148  virtual double &Elem(int i, int j);
149 
150  /// Returns constant reference to a_{ij}.
151  virtual const double &Elem(int i, int j) const;
152 
153  /// Matrix vector multiplication.
154  virtual void Mult(const Vector &x, Vector &y) const;
155 
156  void FullMult(const Vector &x, Vector &y) const
157  { mat->Mult(x, y); mat_e->AddMult(x, y); }
158 
159  virtual void AddMult(const Vector &x, Vector &y, const double a = 1.0) const
160  { mat -> AddMult (x, y, a); }
161 
162  void FullAddMult(const Vector &x, Vector &y) const
163  { mat->AddMult(x, y); mat_e->AddMult(x, y); }
164 
165  double InnerProduct(const Vector &x, const Vector &y) const
166  { return mat->InnerProduct (x, y); }
167 
168  /// Returns a pointer to (approximation) of the matrix inverse.
169  virtual MatrixInverse *Inverse() const;
170 
171  /// Finalizes the matrix initialization.
172  virtual void Finalize(int skip_zeros = 1);
173 
174  /// Returns a reference to the sparse matrix
175  const SparseMatrix &SpMat() const
176  {
177  MFEM_VERIFY(mat, "mat is NULL and can't be dereferenced");
178  return *mat;
179  }
181  {
182  MFEM_VERIFY(mat, "mat is NULL and can't be dereferenced");
183  return *mat;
184  }
185  SparseMatrix *LoseMat() { SparseMatrix *tmp = mat; mat = NULL; return tmp; }
186 
187  /// Returns a reference to the sparse matrix of eliminated b.c.
188  const SparseMatrix &SpMatElim() const
189  {
190  MFEM_VERIFY(mat_e, "mat_e is NULL and can't be dereferenced");
191  return *mat_e;
192  }
194  {
195  MFEM_VERIFY(mat_e, "mat_e is NULL and can't be dereferenced");
196  return *mat_e;
197  }
198 
199  /// Adds new Domain Integrator.
201 
202  /// Adds new Boundary Integrator.
204 
205  /// Adds new interior Face Integrator.
207 
208  /// Adds new boundary Face Integrator.
210 
211  /** @brief Adds new boundary Face Integrator, restricted to specific boundary
212  attributes. */
214  Array<int> &bdr_marker);
215 
216  void operator=(const double a)
217  {
218  if (mat != NULL) { *mat = a; }
219  if (mat_e != NULL) { *mat_e = a; }
220  }
221 
222  /// Assembles the form i.e. sums over all domain/bdr integrators.
223  void Assemble(int skip_zeros = 1);
224 
225  /// Get the finite element space prolongation matrix
226  virtual const Operator *GetProlongation() const
227  { return fes->GetConformingProlongation(); }
228  /// Get the finite element space restriction matrix
229  virtual const Operator *GetRestriction() const
230  { return fes->GetConformingRestriction(); }
231 
232  /** Form the linear system A X = B, corresponding to the current bilinear
233  form and b(.), by applying any necessary transformations such as:
234  eliminating boundary conditions; applying conforming constraints for
235  non-conforming AMR; static condensation; hybridization.
236 
237  The GridFunction-size vector x must contain the essential b.c. The
238  BilinearForm and the LinearForm-size vector b must be assembled.
239 
240  The vector X is initialized with a suitable initial guess: when using
241  hybridization, the vector X is set to zero; otherwise, the essential
242  entries of X are set to the corresponding b.c. and all other entries are
243  set to zero (copy_interior == 0) or copied from x (copy_interior != 0).
244 
245  This method can be called multiple times (with the same ess_tdof_list
246  array) to initialize different right-hand sides and boundary condition
247  values.
248 
249  After solving the linear system, the finite element solution x can be
250  recovered by calling RecoverFEMSolution (with the same vectors X, b, and
251  x).
252 
253  NOTE: If there are no transformations, X simply reuses the data of x. */
254  void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x, Vector &b,
255  SparseMatrix &A, Vector &X, Vector &B,
256  int copy_interior = 0);
257 
258  /// Form the linear system matrix A, see FormLinearSystem for details.
259  void FormSystemMatrix(const Array<int> &ess_tdof_list, SparseMatrix &A);
260 
261  /** Call this method after solving a linear system constructed using the
262  FormLinearSystem method to recover the solution as a GridFunction-size
263  vector in x. Use the same arguments as in the FormLinearSystem call. */
264  virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
265 
266  /// Compute and store internally all element matrices.
267  void ComputeElementMatrices();
268 
269  /// Free the memory used by the element matrices.
271  { delete element_matrices; element_matrices = NULL; }
272 
273  void ComputeElementMatrix(int i, DenseMatrix &elmat);
274  void AssembleElementMatrix(int i, const DenseMatrix &elmat,
275  Array<int> &vdofs, int skip_zeros = 1);
276  void AssembleBdrElementMatrix(int i, const DenseMatrix &elmat,
277  Array<int> &vdofs, int skip_zeros = 1);
278 
279  /** Eliminate essential boundary DOFs from the system. The array
280  'bdr_attr_is_ess' marks boundary attributes that constitute the essential
281  part of the boundary. If d == 0, the diagonal at the essential DOFs is
282  set to 1.0, otherwise it is left the same. */
283  void EliminateEssentialBC(const Array<int> &bdr_attr_is_ess,
284  Vector &sol, Vector &rhs, int d = 0);
285 
286  void EliminateEssentialBC(const Array<int> &bdr_attr_is_ess, int d = 0);
287  /// Perform elimination and set the diagonal entry to the given value
288  void EliminateEssentialBCDiag(const Array<int> &bdr_attr_is_ess,
289  double value);
290 
291  /// Eliminate the given vdofs. NOTE: here, vdofs is a list of DOFs.
292  void EliminateVDofs(const Array<int> &vdofs, Vector &sol, Vector &rhs,
293  int d = 0);
294 
295  /** Eliminate the given vdofs storing the eliminated part internally; this
296  method works in conjunction with EliminateVDofsInRHS and allows
297  elimination of boundary conditions in multiple right-hand sides. In this
298  method, vdofs is a list of DOFs. */
299  void EliminateVDofs(const Array<int> &vdofs, int d = 0);
300 
301  /** Similar to EliminateVDofs but here ess_dofs is a marker
302  (boolean) array on all vdofs (ess_dofs[i] < 0 is true). */
303  void EliminateEssentialBCFromDofs(const Array<int> &ess_dofs, Vector &sol,
304  Vector &rhs, int d = 0);
305 
306  /** Similar to EliminateVDofs but here ess_dofs is a marker
307  (boolean) array on all vdofs (ess_dofs[i] < 0 is true). */
308  void EliminateEssentialBCFromDofs(const Array<int> &ess_dofs, int d = 0);
309  /// Perform elimination and set the diagonal entry to the given value
310  void EliminateEssentialBCFromDofsDiag(const Array<int> &ess_dofs,
311  double value);
312 
313  /** Use the stored eliminated part of the matrix (see EliminateVDofs) to
314  modify r.h.s.; vdofs is a list of DOFs (non-directional, i.e. >= 0). */
315  void EliminateVDofsInRHS(const Array<int> &vdofs, const Vector &x,
316  Vector &b);
317 
318  double FullInnerProduct(const Vector &x, const Vector &y) const
319  { return mat->InnerProduct(x, y) + mat_e->InnerProduct(x, y); }
320 
321  virtual void Update(FiniteElementSpace *nfes = NULL);
322 
323  /// Return the FE space associated with the BilinearForm.
325 
326  /// Destroys bilinear form.
327  virtual ~BilinearForm();
328 };
329 
330 /**
331  Class for assembling of bilinear forms `a(u,v)` defined on different
332  trial and test spaces. The assembled matrix `A` is such that
333 
334  a(u,v) = V^t A U
335 
336  where `U` and `V` are the vectors representing the functions `u` and `v`,
337  respectively. The first argument, `u`, of `a(,)` is in the trial space
338  and the second argument, `v`, is in the test space. Thus,
339 
340  # of rows of A = dimension of the test space and
341  # of cols of A = dimension of the trial space.
342 
343  Both trial and test spaces should be defined on the same mesh.
344 */
345 class MixedBilinearForm : public Matrix
346 {
347 protected:
349 
351 
354  Array<BilinearFormIntegrator*> skt; // trace face integrators
355 
356 public:
358  FiniteElementSpace *te_fes);
359 
360  virtual double& Elem (int i, int j);
361 
362  virtual const double& Elem (int i, int j) const;
363 
364  virtual void Mult (const Vector & x, Vector & y) const;
365 
366  virtual void AddMult (const Vector & x, Vector & y,
367  const double a = 1.0) const;
368 
369  virtual void AddMultTranspose (const Vector & x, Vector & y,
370  const double a = 1.0) const;
371 
372  virtual void MultTranspose (const Vector & x, Vector & y) const
373  { y = 0.0; AddMultTranspose (x, y); }
374 
375  virtual MatrixInverse * Inverse() const;
376 
377  virtual void Finalize (int skip_zeros = 1);
378 
379  /** Extract the associated matrix as SparseMatrix blocks. The number of
380  block rows and columns is given by the vector dimensions (vdim) of the
381  test and trial spaces, respectively. */
382  void GetBlocks(Array2D<SparseMatrix *> &blocks) const;
383 
384  const SparseMatrix &SpMat() const { return *mat; }
385  SparseMatrix &SpMat() { return *mat; }
386  SparseMatrix *LoseMat() { SparseMatrix *tmp = mat; mat = NULL; return tmp; }
387 
389 
391 
392  /** Add a trace face integrator. This type of integrator assembles terms
393  over all faces of the mesh using the face FE from the trial space and the
394  two adjacent volume FEs from the test space. */
396 
398 
400 
402 
403  void operator= (const double a) { *mat = a; }
404 
405  void Assemble (int skip_zeros = 1);
406 
407  /** For partially conforming trial and/or test FE spaces, complete the
408  assembly process by performing A := P2^t A P1 where A is the internal
409  sparse matrix; P1 and P2 are the conforming prolongation matrices of the
410  trial and test FE spaces, respectively. After this call the
411  MixedBilinearForm becomes an operator on the conforming FE spaces. */
412  void ConformingAssemble();
413 
414  void EliminateTrialDofs(Array<int> &bdr_attr_is_ess,
415  Vector &sol, Vector &rhs);
416 
418  Vector &sol, Vector &rhs);
419 
420  virtual void EliminateTestDofs(Array<int> &bdr_attr_is_ess);
421 
422  void Update();
423 
424  virtual ~MixedBilinearForm();
425 };
426 
427 
428 /**
429  Class for constructing the matrix representation of a linear operator,
430  `v = L u`, from one FiniteElementSpace (domain) to another FiniteElementSpace
431  (range). The constructed matrix `A` is such that
432 
433  V = A U
434 
435  where `U` and `V` are the vectors of degrees of freedom representing the
436  functions `u` and `v`, respectively. The dimensions of `A` are
437 
438  number of rows of A = dimension of the range space and
439  number of cols of A = dimension of the domain space.
440 
441  This class is very similar to MixedBilinearForm. One difference is that
442  the linear operator `L` is defined using a special kind of
443  BilinearFormIntegrator (we reuse its functionality instead of defining a
444  new class). The other difference with the MixedBilinearForm class is that
445  the "assembly" process overwrites the global matrix entries using the
446  local element matrices instead of adding them.
447 
448  Note that if we define the bilinear form `b(u,v) := (Lu,v)` using an inner
449  product in the range space, then its matrix representation, `B`, is
450 
451  B = M A, (since V^t B U = b(u,v) = (Lu,v) = V^t M A U)
452 
453  where `M` denotes the mass matrix for the inner product in the range space:
454  `V1^t M V2 = (v1,v2)`. Similarly, if `c(u,w) := (Lu,Lw)` then
455 
456  C = A^t M A.
457 */
459 {
460 public:
462  FiniteElementSpace *range_fes)
463  : MixedBilinearForm(domain_fes, range_fes) { }
464 
466  { AddDomainIntegrator(di); }
467 
469  { AddTraceFaceIntegrator(di); }
470 
472 
473  virtual void Assemble(int skip_zeros = 1);
474 };
475 
476 }
477 
478 #endif
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()
virtual 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:481
Abstract data type for matrix inverse.
Definition: matrix.hpp:58
void EliminateVDofs(const Array< int > &vdofs, Vector &sol, Vector &rhs, int d=0)
Eliminate the given vdofs. NOTE: here, vdofs is a list of DOFs.
void FormSystemMatrix(const Array< int > &ess_tdof_list, SparseMatrix &A)
Form the linear system matrix A, see FormLinearSystem for details.
virtual const Operator * GetProlongation() const
Get the finite element space prolongation matrix.
void EliminateEssentialBCFromDofs(const Array< int > &ess_dofs, Vector &sol, Vector &rhs, int d=0)
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
Array< BilinearFormIntegrator * > * GetFBFI()
MixedBilinearForm(FiniteElementSpace *tr_fes, FiniteElementSpace *te_fes)
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 EliminateEssentialBCFromDofsDiag(const Array< int > &ess_dofs, double value)
Perform elimination and set the diagonal entry to the given value.
void EliminateEssentialBCFromTrialDofs(Array< int > &marked_vdofs, Vector &sol, Vector &rhs)
Array< BilinearFormIntegrator * > bdr
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
SparseMatrix & SpMat()
const SparseMatrix * GetConformingRestriction()
Definition: fespace.cpp:697
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
void EliminateEssentialBCDiag(const Array< int > &bdr_attr_is_ess, double value)
Perform elimination and set the diagonal entry to the given value.
virtual void Update(FiniteElementSpace *nfes=NULL)
void EliminateTrialDofs(Array< int > &bdr_attr_is_ess, Vector &sol, Vector &rhs)
const SparseMatrix * GetConformingProlongation()
Definition: fespace.cpp:690
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)
SparseMatrix & SpMatElim()
SparseMatrix * mat_e
Matrix used to eliminate b.c.
bool StaticCondensationIsEnabled() const
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:22
virtual const Operator * GetRestriction() const
Get the finite element space restriction matrix.
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:475
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
int height
Dimension of the output / number of rows in the matrix.
Definition: operator.hpp:24
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:109
virtual void EliminateTestDofs(Array< int > &bdr_attr_is_ess)
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
void EliminateVDofsInRHS(const Array< int > &vdofs, const Vector &x, Vector &b)
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: y=A^t(x). The default behavior in class Operator is to generate an ...
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:36
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< Array< int > * > bfbfi_marker
Array< int > vdofs
void FreeElementMatrices()
Free the memory used by the element matrices.
Abstract operator.
Definition: operator.hpp:21
void EliminateEssentialBC(const Array< int > &bdr_attr_is_ess, Vector &sol, Vector &rhs, int d=0)
virtual void Assemble(int skip_zeros=1)
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, SparseMatrix &A, Vector &X, Vector &B, int copy_interior=0)
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:593
Array< BilinearFormIntegrator * > * GetDI()
void UseSparsity(int *I, int *J, bool isSorted)
Use the given CSR sparsity pattern to allocate the internal SparseMatrix.
double InnerProduct(const Vector &x, const Vector &y) const
Compute y^t A x.
Definition: sparsemat.cpp:679
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
void EnableStaticCondensation()
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new boundary Face Integrator.