MFEM  v3.3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
bilininteg.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_BILININTEG
13 #define MFEM_BILININTEG
14 
15 #include "../config/config.hpp"
16 #include "nonlininteg.hpp"
17 
18 namespace mfem
19 {
20 
21 /// Abstract base class BilinearFormIntegrator
23 {
24 protected:
27 
28 public:
29  /// Given a particular Finite Element computes the element matrix elmat.
30  virtual void AssembleElementMatrix(const FiniteElement &el,
32  DenseMatrix &elmat);
33 
34  /** Compute the local matrix representation of a bilinear form
35  a(u,v) defined on different trial (given by u) and test
36  (given by v) spaces. The rows in the local matrix correspond
37  to the test dofs and the columns -- to the trial dofs. */
38  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
39  const FiniteElement &test_fe,
41  DenseMatrix &elmat);
42 
43  virtual void AssembleFaceMatrix(const FiniteElement &el1,
44  const FiniteElement &el2,
46  DenseMatrix &elmat);
47 
48  /** Abstract method used for assembling TraceFaceIntegrators in a
49  MixedBilinearForm. */
50  virtual void AssembleFaceMatrix(const FiniteElement &trial_face_fe,
51  const FiniteElement &test_fe1,
52  const FiniteElement &test_fe2,
54  DenseMatrix &elmat);
55 
56  /// Perform the local action of the BilinearFormIntegrator
57  virtual void AssembleElementVector(const FiniteElement &el,
59  const Vector &elfun, Vector &elvect);
60 
61  virtual void AssembleElementGrad(const FiniteElement &el,
63  const Vector &elfun, DenseMatrix &elmat)
64  { AssembleElementMatrix(el, Tr, elmat); }
65 
66  virtual void AssembleFaceGrad(const FiniteElement &el1,
67  const FiniteElement &el2,
69  const Vector &elfun, DenseMatrix &elmat)
70  { AssembleFaceMatrix(el1, el2, Tr, elmat); }
71 
72  virtual void ComputeElementFlux(const FiniteElement &el,
74  Vector &u,
75  const FiniteElement &fluxelem,
76  Vector &flux, int with_coef = 1) { }
77 
78  virtual double ComputeFluxEnergy(const FiniteElement &fluxelem,
80  Vector &flux, Vector *d_energy = NULL)
81  { return 0.0; }
82 
84 };
85 
87 {
88 private:
89  int own_bfi;
91 
92  DenseMatrix bfi_elmat;
93 
94 public:
95  TransposeIntegrator (BilinearFormIntegrator *_bfi, int _own_bfi = 1)
96  { bfi = _bfi; own_bfi = _own_bfi; }
97 
98  virtual void AssembleElementMatrix(const FiniteElement &el,
100  DenseMatrix &elmat);
101 
102  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
103  const FiniteElement &test_fe,
105  DenseMatrix &elmat);
106 
108  virtual void AssembleFaceMatrix(const FiniteElement &el1,
109  const FiniteElement &el2,
111  DenseMatrix &elmat);
112 
113  virtual ~TransposeIntegrator() { if (own_bfi) { delete bfi; } }
114 };
115 
117 {
118 private:
119  int own_bfi;
121 
122 public:
123  LumpedIntegrator (BilinearFormIntegrator *_bfi, int _own_bfi = 1)
124  { bfi = _bfi; own_bfi = _own_bfi; }
125 
126  virtual void AssembleElementMatrix(const FiniteElement &el,
128  DenseMatrix &elmat);
129 
130  virtual ~LumpedIntegrator() { if (own_bfi) { delete bfi; } }
131 };
132 
133 /// Integrator that inverts the matrix assembled by another integrator.
135 {
136 private:
137  int own_integrator;
138  BilinearFormIntegrator *integrator;
139 
140 public:
141  InverseIntegrator(BilinearFormIntegrator *integ, int own_integ = 1)
142  { integrator = integ; own_integrator = own_integ; }
143 
144  virtual void AssembleElementMatrix(const FiniteElement &el,
146  DenseMatrix &elmat);
147 
148  virtual ~InverseIntegrator() { if (own_integrator) { delete integrator; } }
149 };
150 
151 /// Integrator defining a sum of multiple Integrators.
153 {
154 private:
155  int own_integrators;
156  DenseMatrix elem_mat;
157  Array<BilinearFormIntegrator*> integrators;
158 
159 public:
160  SumIntegrator(int own_integs = 1) { own_integrators = own_integs; }
161 
163  { integrators.Append(integ); }
164 
165  virtual void AssembleElementMatrix(const FiniteElement &el,
167  DenseMatrix &elmat);
168 
169  virtual ~SumIntegrator();
170 };
171 
172 /** An abstract class for integrating the product of two scalar basis functions
173  with an optional scalar coefficient. */
175 {
176 public:
177 
178  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
179  const FiniteElement &test_fe,
181  DenseMatrix &elmat);
182 
183  /// Support for use in BilinearForm. Can be used only when appropriate.
184  virtual void AssembleElementMatrix(const FiniteElement &fe,
185  ElementTransformation &Trans,
186  DenseMatrix &elmat)
187  { AssembleElementMatrix2(fe, fe, Trans, elmat); }
188 
189 protected:
190  /// This parameter can be set by derived methods to enable single shape
191  /// evaluation in case CalcTestShape() and CalcTrialShape() return the same
192  /// result if given the same FiniteElement. The default is false.
194 
195  MixedScalarIntegrator() : same_calc_shape(false), Q(NULL) {}
197 
198  inline virtual bool VerifyFiniteElementTypes(
199  const FiniteElement & trial_fe,
200  const FiniteElement & test_fe) const
201  {
202  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
204  }
205 
206  inline virtual const char * FiniteElementTypeFailureMessage() const
207  {
208  return "MixedScalarIntegrator: "
209  "Trial and test spaces must both be scalar fields.";
210  }
211 
212  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
213  const FiniteElement & test_fe,
214  ElementTransformation &Trans)
215  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW(); }
216 
217 
218  inline virtual void CalcTestShape(const FiniteElement & test_fe,
219  ElementTransformation &Trans,
220  Vector & shape)
221  { test_fe.CalcPhysShape(Trans, shape); }
222 
223  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
224  ElementTransformation &Trans,
225  Vector & shape)
226  { trial_fe.CalcPhysShape(Trans, shape); }
227 
228 private:
229 
230  Coefficient *Q;
231 
232 #ifndef MFEM_THREAD_SAFE
233  Vector test_shape;
234  Vector trial_shape;
235 #endif
236 
237 };
238 
239 /** An abstract class for integrating the inner product of two vector basis
240  functions with an optional scalar, vector, or matrix coefficient. */
242 {
243 public:
244 
245  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
246  const FiniteElement &test_fe,
248  DenseMatrix &elmat);
249 
250  /// Support for use in BilinearForm. Can be used only when appropriate.
251  virtual void AssembleElementMatrix(const FiniteElement &fe,
252  ElementTransformation &Trans,
253  DenseMatrix &elmat)
254  { AssembleElementMatrix2(fe, fe, Trans, elmat); }
255 
256 protected:
257  /// This parameter can be set by derived methods to enable single shape
258  /// evaluation in case CalcTestShape() and CalcTrialShape() return the same
259  /// result if given the same FiniteElement. The default is false.
261 
263  : same_calc_shape(false), Q(NULL), VQ(NULL), DQ(NULL), MQ(NULL) {}
265  : same_calc_shape(false), Q(&q), VQ(NULL), DQ(NULL), MQ(NULL) {}
267  : same_calc_shape(false), Q(NULL), VQ(diag?NULL:&dq), DQ(diag?&dq:NULL),
268  MQ(NULL) {}
270  : same_calc_shape(false), Q(NULL), VQ(NULL), DQ(NULL), MQ(&mq) {}
271 
272  inline virtual bool VerifyFiniteElementTypes(
273  const FiniteElement & trial_fe,
274  const FiniteElement & test_fe) const
275  {
276  return (trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
278  }
279 
280  inline virtual const char * FiniteElementTypeFailureMessage() const
281  {
282  return "MixedVectorIntegrator: "
283  "Trial and test spaces must both be vector fields";
284  }
285 
286  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
287  const FiniteElement & test_fe,
288  ElementTransformation &Trans)
289  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW(); }
290 
291 
292  inline virtual void CalcTestShape(const FiniteElement & test_fe,
293  ElementTransformation &Trans,
294  DenseMatrix & shape)
295  { test_fe.CalcVShape(Trans, shape); }
296 
297  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
298  ElementTransformation &Trans,
299  DenseMatrix & shape)
300  { trial_fe.CalcVShape(Trans, shape); }
301 
302 private:
303 
304  Coefficient *Q;
305  VectorCoefficient *VQ;
306  VectorCoefficient *DQ;
307  MatrixCoefficient *MQ;
308 
309 #ifndef MFEM_THREAD_SAFE
310  Vector V;
311  Vector D;
312  DenseMatrix M;
313  DenseMatrix test_shape;
314  DenseMatrix trial_shape;
315  DenseMatrix test_shape_tmp;
316 #endif
317 
318 };
319 
320 /** An abstract class for integrating the product of a scalar basis function and
321  the inner product of a vector basis function with a vector coefficient. In
322  2D the inner product can be replaced with a cross product. */
324 {
325 public:
326 
327  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
328  const FiniteElement &test_fe,
330  DenseMatrix &elmat);
331 
332 protected:
333 
334  MixedScalarVectorIntegrator(VectorCoefficient &vq, bool _transpose = false,
335  bool _cross_2d = false)
336  : VQ(&vq), transpose(_transpose) , cross_2d(_cross_2d) {}
337 
338  inline virtual bool VerifyFiniteElementTypes(
339  const FiniteElement & trial_fe,
340  const FiniteElement & test_fe) const
341  {
342  return ((transpose &&
344  test_fe.GetRangeType() == mfem::FiniteElement::SCALAR ) ||
345  (!transpose &&
348  );
349  }
350 
351  inline virtual const char * FiniteElementTypeFailureMessage() const
352  {
353  if ( transpose )
354  {
355  return "MixedScalarVectorIntegrator: "
356  "Trial space must be a vector field "
357  "and the test space must be a scalar field";
358  }
359  else
360  {
361  return "MixedScalarVectorIntegrator: "
362  "Trial space must be a scalar field "
363  "and the test space must be a vector field";
364  }
365  }
366 
367  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
368  const FiniteElement & test_fe,
369  ElementTransformation &Trans)
370  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW(); }
371 
372 
373  inline virtual void CalcVShape(const FiniteElement & vector_fe,
374  ElementTransformation &Trans,
375  DenseMatrix & shape)
376  { vector_fe.CalcVShape(Trans, shape); }
377 
378  inline virtual void CalcShape(const FiniteElement & scalar_fe,
379  ElementTransformation &Trans,
380  Vector & shape)
381  { scalar_fe.CalcPhysShape(Trans, shape); }
382 
383 private:
384 
385  VectorCoefficient *VQ;
386  bool transpose;
387  bool cross_2d; // In 2D use a cross product rather than a dot product
388 
389 #ifndef MFEM_THREAD_SAFE
390  Vector V;
391  DenseMatrix vshape;
392  Vector shape;
393  Vector vshape_tmp;
394 #endif
395 
396 };
397 
398 /** Class for integrating the bilinear form a(u,v) := (Q u, v) in either 1D, 2D,
399  or 3D and where Q is an optional scalar coefficient, u and v are each in H1
400  or L2. */
402 {
403 public:
406  : MixedScalarIntegrator(q) { same_calc_shape = true; }
407 };
408 
409 /** Class for integrating the bilinear form a(u,v) := (Q u, v) in either 2D, or
410  3D and where Q is a vector coefficient, u is in H1 or L2 and v is in H(Curl)
411  or H(Div). */
413 {
414 public:
417 };
418 
419 /** Class for integrating the bilinear form a(u,v) := (Q D u, v) in 1D where Q
420  is an optional scalar coefficient, u is in H1, and v is in L2. */
422 {
423 public:
426  : MixedScalarIntegrator(q) {}
427 
428 protected:
429  inline virtual bool VerifyFiniteElementTypes(
430  const FiniteElement & trial_fe,
431  const FiniteElement & test_fe) const
432  {
433  return (trial_fe.GetDim() == 1 && test_fe.GetDim() == 1 &&
434  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
436  }
437 
438  inline virtual const char * FiniteElementTypeFailureMessage() const
439  {
440  return "MixedScalarDerivativeIntegrator: "
441  "Trial and test spaces must both be scalar fields in 1D "
442  "and the trial space must implement CaldDShape.";
443  }
444 
445  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
447  Vector & shape)
448  {
449  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
450  trial_fe.CalcPhysDShape(Trans, dshape);
451  }
452 };
453 
454 /** Class for integrating the bilinear form a(u,v) := -(Q u, D v) in 1D where Q
455  is an optional scalar coefficient, u is in L2, and v is in H1. */
457 {
458 public:
461  : MixedScalarIntegrator(q) {}
462 
463 protected:
464  inline virtual bool VerifyFiniteElementTypes(
465  const FiniteElement & trial_fe,
466  const FiniteElement & test_fe) const
467  {
468  return (trial_fe.GetDim() == 1 && test_fe.GetDim() == 1 &&
471  }
472 
473  inline virtual const char * FiniteElementTypeFailureMessage() const
474  {
475  return "MixedScalarWeakDerivativeIntegrator: "
476  "Trial and test spaces must both be scalar fields in 1D "
477  "and the test space must implement CalcDShape with "
478  "map type \"VALUE\".";
479  }
480 
481  inline virtual void CalcTestShape(const FiniteElement & test_fe,
483  Vector & shape)
484  {
485  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
486  test_fe.CalcPhysDShape(Trans, dshape);
487  shape *= -1.0;
488  }
489 };
490 
491 /** Class for integrating the bilinear form a(u,v) := (Q div u, v) in either 2D
492  or 3D where Q is an optional scalar coefficient, u is in H(Div), and v is a
493  scalar field. */
495 {
496 public:
499  : MixedScalarIntegrator(q) {}
500 
501 protected:
502  inline virtual bool VerifyFiniteElementTypes(
503  const FiniteElement & trial_fe,
504  const FiniteElement & test_fe) const
505  {
506  return (trial_fe.GetDerivType() == mfem::FiniteElement::DIV &&
508  }
509 
510  inline virtual const char * FiniteElementTypeFailureMessage() const
511  {
512  return "MixedScalarDivergenceIntegrator: "
513  "Trial must be H(Div) and the test space must be a "
514  "scalar field";
515  }
516 
517  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
518  const FiniteElement & test_fe,
520  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW() - 1; }
521 
522  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
524  Vector & shape)
525  { trial_fe.CalcPhysDivShape(Trans, shape); }
526 };
527 
528 /** Class for integrating the bilinear form a(u,v) := (V div u, v) in either 2D
529  or 3D where V is a vector coefficient, u is in H(Div), and v is a vector
530  field. */
532 {
533 public:
536 
537 protected:
538  inline virtual bool VerifyFiniteElementTypes(
539  const FiniteElement & trial_fe,
540  const FiniteElement & test_fe) const
541  {
542  return (trial_fe.GetDerivType() == mfem::FiniteElement::DIV &&
544  }
545 
546  inline virtual const char * FiniteElementTypeFailureMessage() const
547  {
548  return "MixedVectorDivergenceIntegrator: "
549  "Trial must be H(Div) and the test space must be a "
550  "vector field";
551  }
552 
553  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
554  const FiniteElement & test_fe,
556  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW() - 1; }
557 
558  inline virtual void CalcShape(const FiniteElement & scalar_fe,
560  Vector & shape)
561  { scalar_fe.CalcPhysDivShape(Trans, shape); }
562 };
563 
564 /** Class for integrating the bilinear form a(u,v) := -(Q u, div v) in either 2D
565  or 3D where Q is an optional scalar coefficient, u is in L2 or H1, and v is
566  in H(Div). */
568 {
569 public:
572  : MixedScalarIntegrator(q) {}
573 
574 protected:
575  inline virtual bool VerifyFiniteElementTypes(
576  const FiniteElement & trial_fe,
577  const FiniteElement & test_fe) const
578  {
579  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
580  test_fe.GetDerivType() == mfem::FiniteElement::DIV );
581  }
582 
583  inline virtual const char * FiniteElementTypeFailureMessage() const
584  {
585  return "MixedScalarWeakGradientIntegrator: "
586  "Trial space must be a scalar field "
587  "and the test space must be H(Div)";
588  }
589 
590  virtual void CalcTestShape(const FiniteElement & test_fe,
592  Vector & shape)
593  {
594  test_fe.CalcPhysDivShape(Trans, shape);
595  shape *= -1.0;
596  }
597 };
598 
599 /** Class for integrating the bilinear form a(u,v) := (Q curl u, v) in 2D where
600  Q is an optional scalar coefficient, u is in H(Curl), and v is in L2 or
601  H1. */
603 {
604 public:
607  : MixedScalarIntegrator(q) {}
608 
609 protected:
610  inline virtual bool VerifyFiniteElementTypes(
611  const FiniteElement & trial_fe,
612  const FiniteElement & test_fe) const
613  {
614  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
615  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
617  }
618 
619  inline virtual const char * FiniteElementTypeFailureMessage() const
620  {
621  return "MixedScalarCurlIntegrator: "
622  "Trial must be H(Curl) and the test space must be a "
623  "scalar field";
624  }
625 
626  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
627  const FiniteElement & test_fe,
629  { return trial_fe.GetOrder() + test_fe.GetOrder() + Trans.OrderW() - 1; }
630 
631  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
633  Vector & shape)
634  {
635  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
636  trial_fe.CalcPhysCurlShape(Trans, dshape);
637  }
638 };
639 
640 /** Class for integrating the bilinear form a(u,v) := (Q u, curl v) in 2D where
641  Q is an optional scalar coefficient, u is in L2 or H1, and v is in
642  H(Curl). */
644 {
645 public:
648  : MixedScalarIntegrator(q) {}
649 
650 protected:
651  inline virtual bool VerifyFiniteElementTypes(
652  const FiniteElement & trial_fe,
653  const FiniteElement & test_fe) const
654  {
655  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
658  }
659 
660  inline virtual const char * FiniteElementTypeFailureMessage() const
661  {
662  return "MixedScalarWeakCurlIntegrator: "
663  "Trial space must be a scalar field "
664  "and the test space must be H(Curl)";
665  }
666 
667  inline virtual void CalcTestShape(const FiniteElement & test_fe,
669  Vector & shape)
670  {
671  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
672  test_fe.CalcPhysCurlShape(Trans, dshape);
673  }
674 };
675 
676 /** Class for integrating the bilinear form a(u,v) := (Q u, v) in either 2D or
677  3D and where Q is an optional coefficient (of type scalar, matrix, or
678  diagonal matrix) u and v are each in H(Curl) or H(Div). */
680 {
681 public:
684  : MixedVectorIntegrator(q) { same_calc_shape = true; }
686  : MixedVectorIntegrator(dq, true) { same_calc_shape = true; }
688  : MixedVectorIntegrator(mq) { same_calc_shape = true; }
689 };
690 
691 /** Class for integrating the bilinear form a(u,v) := (V x u, v) in 3D and where
692  V is a vector coefficient u and v are each in H(Curl) or H(Div). */
694 {
695 public:
697  : MixedVectorIntegrator(vq, false) { same_calc_shape = true; }
698 };
699 
700 /** Class for integrating the bilinear form a(u,v) := (V . u, v) in 2D or 3D and
701  where V is a vector coefficient u is in H(Curl) or H(Div) and v is in H1 or
702  L2. */
704 {
705 public:
707  : MixedScalarVectorIntegrator(vq, true) {}
708 
709  inline virtual bool VerifyFiniteElementTypes(
710  const FiniteElement & trial_fe,
711  const FiniteElement & test_fe) const
712  {
713  return (trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
715  }
716 
717  inline virtual const char * FiniteElementTypeFailureMessage() const
718  {
719  return "MixedDotProductIntegrator: "
720  "Trial space must be a vector field "
721  "and the test space must be a scalar field";
722  }
723 };
724 
725 /** Class for integrating the bilinear form a(u,v) := (-V . u, Div v) in 2D or
726  3D and where V is a vector coefficient u is in H(Curl) or H(Div) and v is in
727  RT. */
729 {
730 public:
732  : MixedScalarVectorIntegrator(vq, true) {}
733 
734  inline virtual bool VerifyFiniteElementTypes(
735  const FiniteElement & trial_fe,
736  const FiniteElement & test_fe) const
737  {
738  return (trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
740  test_fe.GetDerivType() == mfem::FiniteElement::DIV );
741  }
742 
743  inline virtual const char * FiniteElementTypeFailureMessage() const
744  {
745  return "MixedWeakGradDotIntegrator: "
746  "Trial space must be a vector field "
747  "and the test space must be a vector field with a divergence";
748  }
749 
750  inline virtual void CalcShape(const FiniteElement & scalar_fe,
752  Vector & shape)
753  { scalar_fe.CalcPhysDivShape(Trans, shape); shape *= -1.0; }
754 };
755 
756 /** Class for integrating the bilinear form a(u,v) := (V x u, Grad v) in 3D and
757  where V is a vector coefficient u is in H(Curl) or H(Div) and v is in H1. */
759 {
760 public:
762  : MixedVectorIntegrator(vq, false) {}
763 
764  inline virtual bool VerifyFiniteElementTypes(
765  const FiniteElement & trial_fe,
766  const FiniteElement & test_fe) const
767  {
768  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
772  }
773 
774  inline virtual const char * FiniteElementTypeFailureMessage() const
775  {
776  return "MixedWeakDivCrossIntegrator: "
777  "Trial space must be a vector field in 3D "
778  "and the test space must be a scalar field with a gradient";
779  }
780 
781  inline virtual void CalcTestShape(const FiniteElement & test_fe,
783  DenseMatrix & shape)
784  { test_fe.CalcPhysDShape(Trans, shape); shape *= -1.0; }
785 };
786 
787 /** Class for integrating the bilinear form a(u,v) := (Q Grad u, Grad v) in 3D
788  or in 2D and where Q is a scalar or matrix coefficient u and v are both in
789  H1. */
791 {
792 public:
795  : MixedVectorIntegrator(q) { same_calc_shape = true; }
797  : MixedVectorIntegrator(dq, true) { same_calc_shape = true; }
799  : MixedVectorIntegrator(mq) { same_calc_shape = true; }
800 
801  inline virtual bool VerifyFiniteElementTypes(
802  const FiniteElement & trial_fe,
803  const FiniteElement & test_fe) const
804  {
805  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
806  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
809  }
810 
811  inline virtual const char * FiniteElementTypeFailureMessage() const
812  {
813  return "MixedGradGradIntegrator: "
814  "Trial and test spaces must both be scalar fields "
815  "with a gradient operator.";
816  }
817 
818  inline virtual int GetIntegrationOrder(const FiniteElement & trial_fe,
819  const FiniteElement & test_fe,
821  {
822  // Same as DiffusionIntegrator
823  return test_fe.Space() == FunctionSpace::Pk ?
824  trial_fe.GetOrder() + test_fe.GetOrder() - 2 :
825  trial_fe.GetOrder() + test_fe.GetOrder() + test_fe.GetDim() - 1;
826  }
827 
828  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
830  DenseMatrix & shape)
831  { trial_fe.CalcPhysDShape(Trans, shape); }
832 
833  inline virtual void CalcTestShape(const FiniteElement & test_fe,
835  DenseMatrix & shape)
836  { test_fe.CalcPhysDShape(Trans, shape); }
837 };
838 
839 /** Class for integrating the bilinear form a(u,v) := (V x Grad u, Grad v) in 3D
840  or in 2D and where V is a vector coefficient u and v are both in H1. */
842 {
843 public:
845  : MixedVectorIntegrator(vq, false) { same_calc_shape = true; }
846 
847  inline virtual bool VerifyFiniteElementTypes(
848  const FiniteElement & trial_fe,
849  const FiniteElement & test_fe) const
850  {
851  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
852  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
855  }
856 
857  inline virtual const char * FiniteElementTypeFailureMessage() const
858  {
859  return "MixedCrossGradGradIntegrator: "
860  "Trial and test spaces must both be scalar fields "
861  "with a gradient operator.";
862  }
863 
864  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
866  DenseMatrix & shape)
867  { trial_fe.CalcPhysDShape(Trans, shape); }
868 
869  inline virtual void CalcTestShape(const FiniteElement & test_fe,
871  DenseMatrix & shape)
872  { test_fe.CalcPhysDShape(Trans, shape); }
873 };
874 
875 /** Class for integrating the bilinear form a(u,v) := (Q Curl u, Curl v) in 3D
876  and where Q is a scalar or matrix coefficient u and v are both in
877  H(Curl). */
879 {
880 public:
883  : MixedVectorIntegrator(q) { same_calc_shape = true; }
885  : MixedVectorIntegrator(dq, true) { same_calc_shape = true; }
887  : MixedVectorIntegrator(mq) { same_calc_shape = true; }
888 
889  inline virtual bool VerifyFiniteElementTypes(
890  const FiniteElement & trial_fe,
891  const FiniteElement & test_fe) const
892  {
893  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
895  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
898  }
899 
900  inline virtual const char * FiniteElementTypeFailureMessage() const
901  {
902  return "MixedCurlCurlIntegrator"
903  "Trial and test spaces must both be vector fields in 3D "
904  "with a curl.";
905  }
906 
907  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
909  DenseMatrix & shape)
910  { trial_fe.CalcPhysCurlShape(Trans, shape); }
911 
912  inline virtual void CalcTestShape(const FiniteElement & test_fe,
914  DenseMatrix & shape)
915  { test_fe.CalcPhysCurlShape(Trans, shape); }
916 };
917 
918 /** Class for integrating the bilinear form a(u,v) := (V x Curl u, Curl v) in 3D
919  and where V is a vector coefficient u and v are both in H(Curl). */
921 {
922 public:
924  : MixedVectorIntegrator(vq, false) { same_calc_shape = true; }
925 
926  inline virtual bool VerifyFiniteElementTypes(
927  const FiniteElement & trial_fe,
928  const FiniteElement & test_fe) const
929  {
930  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
932  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
935  }
936 
937  inline virtual const char * FiniteElementTypeFailureMessage() const
938  {
939  return "MixedCrossCurlCurlIntegrator: "
940  "Trial and test spaces must both be vector fields in 3D "
941  "with a curl.";
942  }
943 
944  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
946  DenseMatrix & shape)
947  { trial_fe.CalcPhysCurlShape(Trans, shape); }
948 
949  inline virtual void CalcTestShape(const FiniteElement & test_fe,
951  DenseMatrix & shape)
952  { test_fe.CalcPhysCurlShape(Trans, shape); }
953 };
954 
955 /** Class for integrating the bilinear form a(u,v) := (V x Curl u, Grad v) in 3D
956  and where V is a vector coefficient u is in H(Curl) and v is in H1. */
958 {
959 public:
961  : MixedVectorIntegrator(vq, false) {}
962 
963  inline virtual bool VerifyFiniteElementTypes(
964  const FiniteElement & trial_fe,
965  const FiniteElement & test_fe) const
966  {
967  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
969  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
972  }
973 
974  inline virtual const char * FiniteElementTypeFailureMessage() const
975  {
976  return "MixedCrossCurlGradIntegrator"
977  "Trial space must be a vector field in 3D with a curl"
978  "and the test space must be a scalar field with a gradient";
979  }
980 
981  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
983  DenseMatrix & shape)
984  { trial_fe.CalcPhysCurlShape(Trans, shape); }
985 
986  inline virtual void CalcTestShape(const FiniteElement & test_fe,
988  DenseMatrix & shape)
989  { test_fe.CalcPhysDShape(Trans, shape); }
990 };
991 
992 /** Class for integrating the bilinear form a(u,v) := (V x Grad u, Curl v) in 3D
993  and where V is a scalar coefficient u is in H1 and v is in H(Curl). */
995 {
996 public:
998  : MixedVectorIntegrator(vq, false) {}
999 
1000  inline virtual bool VerifyFiniteElementTypes(
1001  const FiniteElement & trial_fe,
1002  const FiniteElement & test_fe) const
1003  {
1004  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1005  trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1006  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1008  test_fe.GetDerivType() == mfem::FiniteElement::CURL );
1009  }
1010 
1011  inline virtual const char * FiniteElementTypeFailureMessage() const
1012  {
1013  return "MixedCrossGradCurlIntegrator"
1014  "Trial space must be a scalar field in 3D with a gradient"
1015  "and the test space must be a vector field with a curl";
1016  }
1017 
1018  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
1020  DenseMatrix & shape)
1021  { trial_fe.CalcPhysDShape(Trans, shape); }
1022 
1023  inline virtual void CalcTestShape(const FiniteElement & test_fe,
1025  DenseMatrix & shape)
1026  { test_fe.CalcPhysCurlShape(Trans, shape); }
1027 };
1028 
1029 /** Class for integrating the bilinear form a(u,v) := (V x u, Curl v) in 3D and
1030  where V is a vector coefficient u is in H(Curl) or H(Div) and v is in
1031  H(Curl). */
1033 {
1034 public:
1036  : MixedVectorIntegrator(vq, false) {}
1037 
1038  inline virtual bool VerifyFiniteElementTypes(
1039  const FiniteElement & trial_fe,
1040  const FiniteElement & test_fe) const
1041  {
1042  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1043  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1045  test_fe.GetDerivType() == mfem::FiniteElement::CURL );
1046  }
1047 
1048  inline virtual const char * FiniteElementTypeFailureMessage() const
1049  {
1050  return "MixedWeakCurlCrossIntegrator: "
1051  "Trial space must be a vector field in 3D "
1052  "and the test space must be a vector field with a curl";
1053  }
1054 
1055  inline virtual void CalcTestShape(const FiniteElement & test_fe,
1057  DenseMatrix & shape)
1058  { test_fe.CalcPhysCurlShape(Trans, shape); }
1059 };
1060 
1061 /** Class for integrating the bilinear form a(u,v) := (V x u, Curl v) in 2D and
1062  where V is a vector coefficient u is in H(Curl) or H(Div) and v is in
1063  H(Curl). */
1065 {
1066 public:
1068  : MixedScalarVectorIntegrator(vq, true, true) {}
1069 
1070  inline virtual bool VerifyFiniteElementTypes(
1071  const FiniteElement & trial_fe,
1072  const FiniteElement & test_fe) const
1073  {
1074  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
1075  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1077  test_fe.GetDerivType() == mfem::FiniteElement::CURL );
1078  }
1079 
1080  inline virtual const char * FiniteElementTypeFailureMessage() const
1081  {
1082  return "MixedScalarWeakCurlCrossIntegrator: "
1083  "Trial space must be a vector field in 2D "
1084  "and the test space must be a vector field with a curl";
1085  }
1086 
1087  inline virtual void CalcShape(const FiniteElement & scalar_fe,
1089  Vector & shape)
1090  {
1091  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
1092  scalar_fe.CalcPhysCurlShape(Trans, dshape);
1093  }
1094 };
1095 
1096 /** Class for integrating the bilinear form a(u,v) := (V x Grad u, v) in 3D or
1097  in 2D and where V is a vector coefficient u is in H1 and v is in H(Curl) or
1098  H(Div). */
1100 {
1101 public:
1103  : MixedVectorIntegrator(vq, false) {}
1104 
1105  inline virtual bool VerifyFiniteElementTypes(
1106  const FiniteElement & trial_fe,
1107  const FiniteElement & test_fe) const
1108  {
1109  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1110  trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1111  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1113  }
1114 
1115  inline virtual const char * FiniteElementTypeFailureMessage() const
1116  {
1117  return "MixedCrossGradIntegrator: "
1118  "Trial space must be a scalar field with a gradient operator"
1119  " and the test space must be a vector field both in 3D.";
1120  }
1121 
1122  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
1124  DenseMatrix & shape)
1125  { trial_fe.CalcPhysDShape(Trans, shape); }
1126 
1127  inline virtual void CalcTestShape(const FiniteElement & test_fe,
1129  DenseMatrix & shape)
1130  { test_fe.CalcVShape(Trans, shape); }
1131 };
1132 
1133 /** Class for integrating the bilinear form a(u,v) := (V x Curl u, v) in 3D and
1134  where V is a vector coefficient u is in H(Curl) and v is in H(Curl) or
1135  H(Div). */
1137 {
1138 public:
1140  : MixedVectorIntegrator(vq, false) {}
1141 
1142  inline virtual bool VerifyFiniteElementTypes(
1143  const FiniteElement & trial_fe,
1144  const FiniteElement & test_fe) const
1145  {
1146  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1147  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1148  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
1150  }
1151 
1152  inline virtual const char * FiniteElementTypeFailureMessage() const
1153  {
1154  return "MixedCrossCurlIntegrator: "
1155  "Trial space must be a vector field in 3D with a curl "
1156  "and the test space must be a vector field";
1157  }
1158 
1159  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
1161  DenseMatrix & shape)
1162  { trial_fe.CalcPhysCurlShape(Trans, shape); }
1163 };
1164 
1165 /** Class for integrating the bilinear form a(u,v) := (V x Curl u, v) in 2D and
1166  where V is a vector coefficient u is in H(Curl) and v is in H(Curl) or
1167  H(Div). */
1169 {
1170 public:
1172  : MixedScalarVectorIntegrator(vq, false, true) {}
1173 
1174  inline virtual bool VerifyFiniteElementTypes(
1175  const FiniteElement & trial_fe,
1176  const FiniteElement & test_fe) const
1177  {
1178  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
1179  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1180  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
1182  }
1183 
1184  inline virtual const char * FiniteElementTypeFailureMessage() const
1185  {
1186  return "MixedCrossCurlIntegrator: "
1187  "Trial space must be a vector field in 2D with a curl "
1188  "and the test space must be a vector field";
1189  }
1190 
1191  inline virtual void CalcShape(const FiniteElement & scalar_fe,
1193  Vector & shape)
1194  {
1195  DenseMatrix dshape(shape.GetData(), shape.Size(), 1);
1196  scalar_fe.CalcPhysCurlShape(Trans, dshape); shape *= -1.0;
1197  }
1198 };
1199 
1200 /** Class for integrating the bilinear form a(u,v) := (V x Grad u, v) in 2D and
1201  where V is a vector coefficient u is in H1 and v is in H1 or L2. */
1203 {
1204 public:
1206  : MixedScalarVectorIntegrator(vq, true, true) {}
1207 
1208  inline virtual bool VerifyFiniteElementTypes(
1209  const FiniteElement & trial_fe,
1210  const FiniteElement & test_fe) const
1211  {
1212  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
1213  trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1214  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1216  }
1217 
1218  inline virtual const char * FiniteElementTypeFailureMessage() const
1219  {
1220  return "MixedScalarCrossGradIntegrator: "
1221  "Trial space must be a scalar field in 2D with a gradient "
1222  "and the test space must be a scalar field";
1223  }
1224 
1225  inline virtual void CalcVShape(const FiniteElement & vector_fe,
1227  DenseMatrix & shape)
1228  { vector_fe.CalcPhysDShape(Trans, shape); }
1229 };
1230 
1231 /** Class for integrating the bilinear form a(u,v) := (V x u, v) in 2D and where
1232  V is a vector coefficient u is in ND or RT and v is in H1 or L2. */
1234 {
1235 public:
1237  : MixedScalarVectorIntegrator(vq, true, true) {}
1238 
1239  inline virtual bool VerifyFiniteElementTypes(
1240  const FiniteElement & trial_fe,
1241  const FiniteElement & test_fe) const
1242  {
1243  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
1244  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1246  }
1247 
1248  inline virtual const char * FiniteElementTypeFailureMessage() const
1249  {
1250  return "MixedScalarCrossProductIntegrator: "
1251  "Trial space must be a vector field in 2D "
1252  "and the test space must be a scalar field";
1253  }
1254 };
1255 
1256 /** Class for integrating the bilinear form a(u,v) := (V x z u, v) in 2D and
1257  where V is a vector coefficient u is in H1 or L2 and v is in ND or RT. */
1259 {
1260 public:
1262  : MixedScalarVectorIntegrator(vq, false, true) {}
1263 
1264  inline virtual bool VerifyFiniteElementTypes(
1265  const FiniteElement & trial_fe,
1266  const FiniteElement & test_fe) const
1267  {
1268  return (trial_fe.GetDim() == 2 && test_fe.GetDim() == 2 &&
1269  trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1271  }
1272 
1273  inline virtual const char * FiniteElementTypeFailureMessage() const
1274  {
1275  return "MixedScalarWeakCrossProductIntegrator: "
1276  "Trial space must be a scalar field in 2D "
1277  "and the test space must be a vector field";
1278  }
1279 
1280  inline virtual void CalcShape(const FiniteElement & scalar_fe,
1282  Vector & shape)
1283  { scalar_fe.CalcPhysShape(Trans, shape); shape *= -1.0; }
1284 };
1285 
1286 /** Class for integrating the bilinear form a(u,v) := (V . Grad u, v) in 2D or
1287  3D and where V is a vector coefficient, u is in H1 and v is in H1 or L2. */
1289 {
1290 public:
1292  : MixedScalarVectorIntegrator(vq, true) {}
1293 
1294  inline virtual bool VerifyFiniteElementTypes(
1295  const FiniteElement & trial_fe,
1296  const FiniteElement & test_fe) const
1297  {
1298  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1299  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1301  }
1302 
1303  inline virtual const char * FiniteElementTypeFailureMessage() const
1304  {
1305  return "MixedDirectionalDerivativeIntegrator: "
1306  "Trial space must be a scalar field with a gradient "
1307  "and the test space must be a scalar field";
1308  }
1309 
1310  inline virtual void CalcVShape(const FiniteElement & vector_fe,
1312  DenseMatrix & shape)
1313  { vector_fe.CalcPhysDShape(Trans, shape); }
1314 };
1315 
1316 /** Class for integrating the bilinear form a(u,v) := (-V . Grad u, Div v) in 2D
1317  or 3D and where V is a vector coefficient, u is in H1 and v is in RT. */
1319 {
1320 public:
1322  : MixedScalarVectorIntegrator(vq, true) {}
1323 
1324  inline virtual bool VerifyFiniteElementTypes(
1325  const FiniteElement & trial_fe,
1326  const FiniteElement & test_fe) const
1327  {
1328  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1329  trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1331  test_fe.GetDerivType() == mfem::FiniteElement::DIV );
1332  }
1333 
1334  inline virtual const char * FiniteElementTypeFailureMessage() const
1335  {
1336  return "MixedGradDivIntegrator: "
1337  "Trial space must be a scalar field with a gradient"
1338  "and the test space must be a vector field with a divergence";
1339  }
1340 
1341  inline virtual void CalcVShape(const FiniteElement & vector_fe,
1343  DenseMatrix & shape)
1344  { vector_fe.CalcPhysDShape(Trans, shape); shape *= -1.0; }
1345 
1346  inline virtual void CalcShape(const FiniteElement & scalar_fe,
1348  Vector & shape)
1349  { scalar_fe.CalcPhysDivShape(Trans, shape); }
1350 };
1351 
1352 /** Class for integrating the bilinear form a(u,v) := (-V Div u, Grad v) in 2D
1353  or 3D and where V is a vector coefficient, u is in RT and v is in H1. */
1355 {
1356 public:
1358  : MixedScalarVectorIntegrator(vq, false) {}
1359 
1360  inline virtual bool VerifyFiniteElementTypes(
1361  const FiniteElement & trial_fe,
1362  const FiniteElement & test_fe) const
1363  {
1364  return (trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1365  trial_fe.GetDerivType() == mfem::FiniteElement::DIV &&
1368  );
1369  }
1370 
1371  inline virtual const char * FiniteElementTypeFailureMessage() const
1372  {
1373  return "MixedDivGradIntegrator: "
1374  "Trial space must be a vector field with a divergence"
1375  "and the test space must be a scalar field with a gradient";
1376  }
1377 
1378  inline virtual void CalcVShape(const FiniteElement & vector_fe,
1380  DenseMatrix & shape)
1381  { vector_fe.CalcPhysDShape(Trans, shape); shape *= -1.0; }
1382 
1383  inline virtual void CalcShape(const FiniteElement & scalar_fe,
1385  Vector & shape)
1386  { scalar_fe.CalcPhysDivShape(Trans, shape); }
1387 };
1388 
1389 /** Class for integrating the bilinear form a(u,v) := (-V u, Grad v) in 2D or 3D
1390  and where V is a vector coefficient, u is in H1 and v is in H1. */
1392 {
1393 public:
1395  : MixedScalarVectorIntegrator(vq, false) {}
1396 
1397  inline virtual bool VerifyFiniteElementTypes(
1398  const FiniteElement & trial_fe,
1399  const FiniteElement & test_fe) const
1400  {
1401  return (trial_fe.GetRangeType() == mfem::FiniteElement::SCALAR &&
1403  test_fe.GetDerivType() == mfem::FiniteElement::GRAD );
1404  }
1405 
1406  inline virtual const char * FiniteElementTypeFailureMessage() const
1407  {
1408  return "MixedScalarWeakDivergenceIntegrator: "
1409  "Trial space must be a scalar field "
1410  "and the test space must be a scalar field with a gradient";
1411  }
1412 
1413  inline virtual void CalcVShape(const FiniteElement & vector_fe,
1415  DenseMatrix & shape)
1416  { vector_fe.CalcPhysDShape(Trans, shape); shape *= -1.0; }
1417 };
1418 
1419 /** Class for integrating the bilinear form a(u,v) := (Q grad u, v) in either 2D
1420  or 3D and where Q is an optional coefficient (of type scalar, matrix, or
1421  diagonal matrix) u is in H1 and v is in H(Curl) or H(Div). */
1423 {
1424 public:
1427  : MixedVectorIntegrator(q) {}
1429  : MixedVectorIntegrator(dq, true) {}
1431  : MixedVectorIntegrator(mq) {}
1432 
1433 protected:
1434  inline virtual bool VerifyFiniteElementTypes(
1435  const FiniteElement & trial_fe,
1436  const FiniteElement & test_fe) const
1437  {
1438  return (trial_fe.GetDerivType() == mfem::FiniteElement::GRAD &&
1440  }
1441 
1442  inline virtual const char * FiniteElementTypeFailureMessage() const
1443  {
1444  return "MixedVectorGradientIntegrator: "
1445  "Trial spaces must be H1 and the test space must be a "
1446  "vector field in 2D or 3D";
1447  }
1448 
1449  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
1451  DenseMatrix & shape)
1452  {
1453  trial_fe.CalcPhysDShape(Trans, shape);
1454  }
1455 
1456 private:
1457  DenseMatrix Jinv;
1458 };
1459 
1460 /** Class for integrating the bilinear form a(u,v) := (Q curl u, v) in 3D and
1461  where Q is an optional coefficient (of type scalar, matrix, or diagonal
1462  matrix) u is in H(Curl) and v is in H(Div) or H(Curl). */
1464 {
1465 public:
1468  : MixedVectorIntegrator(q) {}
1470  : MixedVectorIntegrator(dq, true) {}
1472  : MixedVectorIntegrator(mq) {}
1473 
1474 protected:
1475  inline virtual bool VerifyFiniteElementTypes(
1476  const FiniteElement & trial_fe,
1477  const FiniteElement & test_fe) const
1478  {
1479  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1480  trial_fe.GetDerivType() == mfem::FiniteElement::CURL &&
1482  }
1483 
1484  inline virtual const char * FiniteElementTypeFailureMessage() const
1485  {
1486  return "MixedVectorCurlIntegrator: "
1487  "Trial space must be H(Curl) and the test space must be a "
1488  "vector field in 3D";
1489  }
1490 
1491  inline virtual void CalcTrialShape(const FiniteElement & trial_fe,
1493  DenseMatrix & shape)
1494  {
1495  trial_fe.CalcPhysCurlShape(Trans, shape);
1496  }
1497 };
1498 
1499 /** Class for integrating the bilinear form a(u,v) := (Q u, curl v) in 3D and
1500  where Q is an optional coefficient (of type scalar, matrix, or diagonal
1501  matrix) u is in H(Div) or H(Curl) and v is in H(Curl). */
1503 {
1504 public:
1507  : MixedVectorIntegrator(q) {}
1509  : MixedVectorIntegrator(dq, true) {}
1511  : MixedVectorIntegrator(mq) {}
1512 
1513 protected:
1514  inline virtual bool VerifyFiniteElementTypes(
1515  const FiniteElement & trial_fe,
1516  const FiniteElement & test_fe) const
1517  {
1518  return (trial_fe.GetDim() == 3 && test_fe.GetDim() == 3 &&
1519  trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1520  test_fe.GetDerivType() == mfem::FiniteElement::CURL );
1521  }
1522 
1523  inline virtual const char * FiniteElementTypeFailureMessage() const
1524  {
1525  return "MixedVectorWeakCurlIntegrator: "
1526  "Trial space must be vector field in 3D and the "
1527  "test space must be H(Curl)";
1528  }
1529 
1530  inline virtual void CalcTestShape(const FiniteElement & test_fe,
1532  DenseMatrix & shape)
1533  {
1534  test_fe.CalcPhysCurlShape(Trans, shape);
1535  }
1536 };
1537 
1538 /** Class for integrating the bilinear form a(u,v) := - (Q u, grad v) in either
1539  2D or 3D and where Q is an optional coefficient (of type scalar, matrix, or
1540  diagonal matrix) u is in H(Div) or H(Curl) and v is in H1. */
1542 {
1543 public:
1546  : MixedVectorIntegrator(q) {}
1548  : MixedVectorIntegrator(dq, true) {}
1550  : MixedVectorIntegrator(mq) {}
1551 
1552 protected:
1553  inline virtual bool VerifyFiniteElementTypes(
1554  const FiniteElement & trial_fe,
1555  const FiniteElement & test_fe) const
1556  {
1557  return (trial_fe.GetRangeType() == mfem::FiniteElement::VECTOR &&
1558  test_fe.GetDerivType() == mfem::FiniteElement::GRAD );
1559  }
1560 
1561  inline virtual const char * FiniteElementTypeFailureMessage() const
1562  {
1563  return "MixedVectorWeakDivergenceIntegrator: "
1564  "Trial space must be vector field and the "
1565  "test space must be H1";
1566  }
1567 
1568  inline virtual void CalcTestShape(const FiniteElement & test_fe,
1570  DenseMatrix & shape)
1571  {
1572  test_fe.CalcPhysDShape(Trans, shape);
1573  shape *= -1.0;
1574  }
1575 };
1576 
1577 /** Class for integrating the bilinear form a(u,v) := (Q grad u, grad v) where Q
1578  can be a scalar or a matrix coefficient. */
1580 {
1581 private:
1582  Vector vec, pointflux, shape;
1583 #ifndef MFEM_THREAD_SAFE
1584  DenseMatrix dshape, dshapedxt, invdfdx, mq;
1585  DenseMatrix te_dshape, te_dshapedxt;
1586 #endif
1587  Coefficient *Q;
1588  MatrixCoefficient *MQ;
1589 
1590 public:
1591  /// Construct a diffusion integrator with coefficient Q = 1
1592  DiffusionIntegrator() { Q = NULL; MQ = NULL; }
1593 
1594  /// Construct a diffusion integrator with a scalar coefficient q
1595  DiffusionIntegrator (Coefficient &q) : Q(&q) { MQ = NULL; }
1596 
1597  /// Construct a diffusion integrator with a matrix coefficient q
1598  DiffusionIntegrator (MatrixCoefficient &q) : MQ(&q) { Q = NULL; }
1599 
1600  /** Given a particular Finite Element
1601  computes the element stiffness matrix elmat. */
1602  virtual void AssembleElementMatrix(const FiniteElement &el,
1604  DenseMatrix &elmat);
1605  /** Given a trial and test Finite Element computes the element stiffness
1606  matrix elmat. */
1607  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1608  const FiniteElement &test_fe,
1610  DenseMatrix &elmat);
1611 
1612  /// Perform the local action of the BilinearFormIntegrator
1613  virtual void AssembleElementVector(const FiniteElement &el,
1615  const Vector &elfun, Vector &elvect);
1616 
1617  virtual void ComputeElementFlux(const FiniteElement &el,
1619  Vector &u, const FiniteElement &fluxelem,
1620  Vector &flux, int with_coef = 1);
1621 
1622  virtual double ComputeFluxEnergy(const FiniteElement &fluxelem,
1624  Vector &flux, Vector *d_energy = NULL);
1625 };
1626 
1627 /** Class for local mass matrix assembling a(u,v) := (Q u, v) */
1629 {
1630 protected:
1631 #ifndef MFEM_THREAD_SAFE
1633 #endif
1635 
1636 public:
1638  : BilinearFormIntegrator(ir) { Q = NULL; }
1639  /// Construct a mass integrator with coefficient q
1641  : BilinearFormIntegrator(ir), Q(&q) { }
1642 
1643  /** Given a particular Finite Element
1644  computes the element mass matrix elmat. */
1645  virtual void AssembleElementMatrix(const FiniteElement &el,
1647  DenseMatrix &elmat);
1648  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1649  const FiniteElement &test_fe,
1651  DenseMatrix &elmat);
1652 };
1653 
1655 {
1656 public:
1658 
1660 
1661  virtual void AssembleFaceMatrix(const FiniteElement &el1,
1662  const FiniteElement &el2,
1664  DenseMatrix &elmat);
1665 };
1666 
1667 /// alpha (q . grad u, v)
1669 {
1670 private:
1671 #ifndef MFEM_THREAD_SAFE
1672  DenseMatrix dshape, adjJ, Q_ir;
1673  Vector shape, vec2, BdFidxT;
1674 #endif
1675  VectorCoefficient &Q;
1676  double alpha;
1677 
1678 public:
1680  : Q(q) { alpha = a; }
1681  virtual void AssembleElementMatrix(const FiniteElement &,
1683  DenseMatrix &);
1684 };
1685 
1686 /// alpha (q . grad u, v) using the "group" FE discretization
1688 {
1689 private:
1690  DenseMatrix dshape, adjJ, Q_nodal, grad;
1691  Vector shape;
1692  VectorCoefficient &Q;
1693  double alpha;
1694 
1695 public:
1697  : Q(q) { alpha = a; }
1698  virtual void AssembleElementMatrix(const FiniteElement &,
1700  DenseMatrix &);
1701 };
1702 
1703 /** Class for integrating the bilinear form a(u,v) := (Q u, v),
1704  where u=(u1,...,un) and v=(v1,...,vn); ui and vi are defined
1705  by scalar FE through standard transformation. */
1707 {
1708 private:
1709  Vector shape, te_shape, vec;
1710  DenseMatrix partelmat;
1711  DenseMatrix mcoeff;
1712  Coefficient *Q;
1713  VectorCoefficient *VQ;
1714  MatrixCoefficient *MQ;
1715 
1716  int Q_order;
1717 
1718 public:
1719  /// Construct an integrator with coefficient 1.0
1721  { Q = NULL; VQ = NULL; MQ = NULL; Q_order = 0; }
1722  /** Construct an integrator with scalar coefficient q.
1723  If possible, save memory by using a scalar integrator since
1724  the resulting matrix is block diagonal with the same diagonal
1725  block repeated. */
1727  : Q(&q) { VQ = NULL; MQ = NULL; Q_order = qo; }
1729  : BilinearFormIntegrator(ir), Q(&q)
1730  { VQ = NULL; MQ = NULL; Q_order = 0; }
1731  /// Construct an integrator with diagonal coefficient q
1733  : VQ(&q) { Q = NULL; MQ = NULL; Q_order = qo; }
1734  /// Construct an integrator with matrix coefficient q
1736  : MQ(&q) { Q = NULL; VQ = NULL; Q_order = qo; }
1737 
1738  virtual void AssembleElementMatrix(const FiniteElement &el,
1740  DenseMatrix &elmat);
1741  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1742  const FiniteElement &test_fe,
1744  DenseMatrix &elmat);
1745 };
1746 
1747 
1748 /** Class for integrating (div u, p) where u is a vector field given
1749  by VectorFiniteElement through Piola transformation (for RT
1750  elements); p is scalar function given by FiniteElement through
1751  standard transformation. Here, u is the trial function and p is
1752  the test function.
1753  Note: the element matrix returned by AssembleElementMatrix2
1754  does NOT depend on the ElementTransformation Trans. */
1756 {
1757 private:
1758  Coefficient *Q;
1759 #ifndef MFEM_THREAD_SAFE
1760  Vector divshape, shape;
1761 #endif
1762 public:
1765  virtual void AssembleElementMatrix(const FiniteElement &el,
1767  DenseMatrix &elmat) { }
1768  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1769  const FiniteElement &test_fe,
1771  DenseMatrix &elmat);
1772 };
1773 
1774 
1775 /** Integrator for `(-Q u, grad v)` for Nedelec (`u`) and H1 (`v`) elements.
1776  This is equivalent to a weak divergence of the Nedelec basis functions. */
1778 {
1779 private:
1780  Coefficient *Q;
1781 #ifndef MFEM_THREAD_SAFE
1782  DenseMatrix dshape;
1783  DenseMatrix dshapedxt;
1784  DenseMatrix vshape;
1785  DenseMatrix invdfdx;
1786 #endif
1787 public:
1790  virtual void AssembleElementMatrix(const FiniteElement &el,
1792  DenseMatrix &elmat) { }
1793  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1794  const FiniteElement &test_fe,
1796  DenseMatrix &elmat);
1797 };
1798 
1799 /** Integrator for (curl u, v) for Nedelec and RT elements. If the trial and
1800  test spaces are switched, assembles the form (u, curl v). */
1802 {
1803 private:
1804  Coefficient *Q;
1805 #ifndef MFEM_THREAD_SAFE
1806  DenseMatrix curlshapeTrial;
1807  DenseMatrix vshapeTest;
1808  DenseMatrix curlshapeTrial_dFT;
1809 #endif
1810 public:
1811  VectorFECurlIntegrator() { Q = NULL; }
1813  virtual void AssembleElementMatrix(const FiniteElement &el,
1815  DenseMatrix &elmat) { }
1816  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1817  const FiniteElement &test_fe,
1819  DenseMatrix &elmat);
1820 };
1821 
1822 
1823 /// Class for integrating (Q D_i(u), v); u and v are scalars
1825 {
1826 private:
1827  Coefficient & Q;
1828  int xi;
1829  DenseMatrix dshape, dshapedxt, invdfdx;
1830  Vector shape, dshapedxi;
1831 public:
1832  DerivativeIntegrator(Coefficient &q, int i) : Q(q), xi(i) { }
1833  virtual void AssembleElementMatrix(const FiniteElement &el,
1835  DenseMatrix &elmat)
1836  { AssembleElementMatrix2(el,el,Trans,elmat); }
1837  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1838  const FiniteElement &test_fe,
1840  DenseMatrix &elmat);
1841 };
1842 
1843 /// Integrator for (curl u, curl v) for Nedelec elements
1845 {
1846 private:
1847  Vector vec, pointflux;
1848 #ifndef MFEM_THREAD_SAFE
1849  DenseMatrix curlshape, curlshape_dFt, M;
1850  DenseMatrix vshape, projcurl;
1851 #endif
1852  Coefficient *Q;
1853  MatrixCoefficient *MQ;
1854 
1855 public:
1856  CurlCurlIntegrator() { Q = NULL; MQ = NULL; }
1857  /// Construct a bilinear form integrator for Nedelec elements
1858  CurlCurlIntegrator(Coefficient &q) : Q(&q) { MQ = NULL; }
1859  CurlCurlIntegrator(MatrixCoefficient &m) : MQ(&m) { Q = NULL; }
1860 
1861  /* Given a particular Finite Element, compute the
1862  element curl-curl matrix elmat */
1863  virtual void AssembleElementMatrix(const FiniteElement &el,
1865  DenseMatrix &elmat);
1866 
1867  virtual void ComputeElementFlux(const FiniteElement &el,
1869  Vector &u, const FiniteElement &fluxelem,
1870  Vector &flux, int with_coef);
1871 
1872  virtual double ComputeFluxEnergy(const FiniteElement &fluxelem,
1874  Vector &flux, Vector *d_energy = NULL);
1875 };
1876 
1877 /** Integrator for (curl u, curl v) for FE spaces defined by 'dim' copies of a
1878  scalar FE space. */
1880 {
1881 private:
1882 #ifndef MFEM_THREAD_SAFE
1883  DenseMatrix dshape_hat, dshape, curlshape, Jadj, grad_hat, grad;
1884 #endif
1885  Coefficient *Q;
1886 
1887 public:
1889 
1891 
1892  /// Assemble an element matrix
1893  virtual void AssembleElementMatrix(const FiniteElement &el,
1895  DenseMatrix &elmat);
1896  /// Compute element energy: (1/2) (curl u, curl u)_E
1897  virtual double GetElementEnergy(const FiniteElement &el,
1899  const Vector &elfun);
1900 };
1901 
1902 /// Integrator for (Q u, v) for VectorFiniteElements
1904 {
1905 private:
1906  Coefficient *Q;
1907  VectorCoefficient *VQ;
1908  MatrixCoefficient *MQ;
1909  void Init(Coefficient *q, VectorCoefficient *vq, MatrixCoefficient *mq)
1910  { Q = q; VQ = vq; MQ = mq; }
1911 
1912 #ifndef MFEM_THREAD_SAFE
1913  Vector shape;
1914  Vector D;
1915  DenseMatrix K;
1916  DenseMatrix test_vshape;
1917  DenseMatrix trial_vshape;
1918 #endif
1919 
1920 public:
1921  VectorFEMassIntegrator() { Init(NULL, NULL, NULL); }
1922  VectorFEMassIntegrator(Coefficient *_q) { Init(_q, NULL, NULL); }
1923  VectorFEMassIntegrator(Coefficient &q) { Init(&q, NULL, NULL); }
1924  VectorFEMassIntegrator(VectorCoefficient *_vq) { Init(NULL, _vq, NULL); }
1925  VectorFEMassIntegrator(VectorCoefficient &vq) { Init(NULL, &vq, NULL); }
1926  VectorFEMassIntegrator(MatrixCoefficient *_mq) { Init(NULL, NULL, _mq); }
1927  VectorFEMassIntegrator(MatrixCoefficient &mq) { Init(NULL, NULL, &mq); }
1928 
1929  virtual void AssembleElementMatrix(const FiniteElement &el,
1931  DenseMatrix &elmat);
1932  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1933  const FiniteElement &test_fe,
1935  DenseMatrix &elmat);
1936 };
1937 
1938 /** Integrator for (Q div u, p) where u=(v1,...,vn) and all vi are in the same
1939  scalar FE space; p is also in a (different) scalar FE space. */
1941 {
1942 private:
1943  Coefficient *Q;
1944 
1945  Vector shape;
1946  Vector divshape;
1947  DenseMatrix dshape;
1948  DenseMatrix gshape;
1949  DenseMatrix Jadj;
1950 
1951 public:
1955 
1956  virtual void AssembleElementMatrix2(const FiniteElement &trial_fe,
1957  const FiniteElement &test_fe,
1959  DenseMatrix &elmat);
1960 };
1961 
1962 /// (Q div u, div v) for RT elements
1964 {
1965 private:
1966  Coefficient *Q;
1967 
1968 #ifndef MFEM_THREAD_SAFE
1969  Vector divshape;
1970 #endif
1971 
1972 public:
1973  DivDivIntegrator() { Q = NULL; }
1975 
1976  virtual void AssembleElementMatrix(const FiniteElement &el,
1978  DenseMatrix &elmat);
1979 };
1980 
1981 /** Integrator for
1982  (Q grad u, grad v) = sum_i (Q grad u_i, grad v_i) e_i e_i^T
1983  for FE spaces defined by 'dim' copies of a scalar FE space. Where e_i
1984  is the unit vector in the i-th direction. The resulting local element
1985  matrix is a block-diagonal matrix consisting of 'dim' copies of a scalar
1986  diffusion matrix in each diagonal block. */
1988 {
1989 private:
1990  Coefficient *Q;
1991 
1992  DenseMatrix Jinv;
1993  DenseMatrix dshape;
1994  DenseMatrix gshape;
1995  DenseMatrix pelmat;
1996 
1997 public:
2000 
2001  virtual void AssembleElementMatrix(const FiniteElement &el,
2003  DenseMatrix &elmat);
2004  virtual void AssembleElementVector(const FiniteElement &el,
2006  const Vector &elfun, Vector &elvect);
2007 };
2008 
2009 /** Integrator for the linear elasticity form:
2010  a(u,v) = (lambda div(u), div(v)) + (2 mu e(u), e(v)),
2011  where e(v) = (1/2) (grad(v) + grad(v)^T).
2012  This is a 'Vector' integrator, i.e. defined for FE spaces
2013  using multiple copies of a scalar FE space. */
2015 {
2016 private:
2017  double q_lambda, q_mu;
2018  Coefficient *lambda, *mu;
2019 
2020 #ifndef MFEM_THREAD_SAFE
2021  DenseMatrix dshape, Jinv, gshape, pelmat;
2022  Vector divshape;
2023 #endif
2024 
2025 public:
2027  { lambda = &l; mu = &m; }
2028  /** With this constructor lambda = q_l * m and mu = q_m * m;
2029  if dim * q_l + 2 * q_m = 0 then trace(sigma) = 0. */
2030  ElasticityIntegrator(Coefficient &m, double q_l, double q_m)
2031  { lambda = NULL; mu = &m; q_lambda = q_l; q_mu = q_m; }
2032 
2033  virtual void AssembleElementMatrix(const FiniteElement &,
2035  DenseMatrix &);
2036 };
2037 
2038 /** Integrator for the DG form:
2039  alpha < rho_u (u.n) {v},[w] > + beta < rho_u |u.n| [v],[w] >,
2040  where v and w are the trial and test variables, respectively, and rho/u are
2041  given scalar/vector coefficients. The vector coefficient, u, is assumed to
2042  be continuous across the faces and when given the scalar coefficient, rho,
2043  is assumed to be discontinuous. The integrator uses the upwind value of rho,
2044  rho_u, which is value from the side into which the vector coefficient, u,
2045  points. */
2047 {
2048 private:
2049  Coefficient *rho;
2050  VectorCoefficient *u;
2051  double alpha, beta;
2052 
2053  Vector shape1, shape2;
2054 
2055 public:
2056  /// Construct integrator with rho = 1.
2057  DGTraceIntegrator(VectorCoefficient &_u, double a, double b)
2058  { rho = NULL; u = &_u; alpha = a; beta = b; }
2059 
2061  double a, double b)
2062  { rho = &_rho; u = &_u; alpha = a; beta = b; }
2063 
2065  virtual void AssembleFaceMatrix(const FiniteElement &el1,
2066  const FiniteElement &el2,
2068  DenseMatrix &elmat);
2069 };
2070 
2071 /** Integrator for the DG form:
2072 
2073  - < {(Q grad(u)).n}, [v] > + sigma < [u], {(Q grad(v)).n} >
2074  + kappa < {h^{-1} Q} [u], [v] >,
2075 
2076  where Q is a scalar or matrix diffusion coefficient and u, v are the trial
2077  and test spaces, respectively. The parameters sigma and kappa determine the
2078  DG method to be used (when this integrator is added to the "broken"
2079  DiffusionIntegrator):
2080  * sigma = -1, kappa >= kappa0: symm. interior penalty (IP or SIPG) method,
2081  * sigma = +1, kappa > 0: non-symmetric interior penalty (NIPG) method,
2082  * sigma = +1, kappa = 0: the method of Baumann and Oden. */
2084 {
2085 protected:
2088  double sigma, kappa;
2089 
2090  // these are not thread-safe!
2093 
2094 public:
2095  DGDiffusionIntegrator(const double s, const double k)
2096  : Q(NULL), MQ(NULL), sigma(s), kappa(k) { }
2097  DGDiffusionIntegrator(Coefficient &q, const double s, const double k)
2098  : Q(&q), MQ(NULL), sigma(s), kappa(k) { }
2099  DGDiffusionIntegrator(MatrixCoefficient &q, const double s, const double k)
2100  : Q(NULL), MQ(&q), sigma(s), kappa(k) { }
2102  virtual void AssembleFaceMatrix(const FiniteElement &el1,
2103  const FiniteElement &el2,
2105  DenseMatrix &elmat);
2106 };
2107 
2108 /** Integrator for the DG elasticity form, for the formulations see:
2109  - PhD Thesis of Jonas De Basabe, High-Order Finite %Element Methods for
2110  Seismic Wave Propagation, UT Austin, 2009, p. 23, and references therein
2111  - Peter Hansbo and Mats G. Larson, Discontinuous Galerkin and the
2112  Crouzeix-Raviart %Element: Application to Elasticity, PREPRINT 2000-09,
2113  p.3
2114 
2115  \f[
2116  - \left< \{ \tau(u) \}, [v] \right> + \alpha \left< \{ \tau(v) \}, [u]
2117  \right> + \kappa \left< h^{-1} \{ \lambda + 2 \mu \} [u], [v] \right>
2118  \f]
2119 
2120  where \f$ \left<u, v\right> = \int_{F} u \cdot v \f$, and \f$ F \f$ is a
2121  face which is either a boundary face \f$ F_b \f$ of an element \f$ K \f$ or
2122  an interior face \f$ F_i \f$ separating elements \f$ K_1 \f$ and \f$ K_2 \f$.
2123 
2124  In the bilinear form above \f$ \tau(u) \f$ is traction, and it's also
2125  \f$ \tau(u) = \sigma(u) \cdot \vec{n} \f$, where \f$ \sigma(u) \f$ is
2126  stress, and \f$ \vec{n} \f$ is the unit normal vector w.r.t. to \f$ F \f$.
2127 
2128  In other words, we have
2129  \f[
2130  - \left< \{ \sigma(u) \cdot \vec{n} \}, [v] \right> + \alpha \left< \{
2131  \sigma(v) \cdot \vec{n} \}, [u] \right> + \kappa \left< h^{-1} \{
2132  \lambda + 2 \mu \} [u], [v] \right>
2133  \f]
2134 
2135  For isotropic media
2136  \f[
2137  \begin{split}
2138  \sigma(u) &= \lambda \nabla \cdot u I + 2 \mu \varepsilon(u) \\
2139  &= \lambda \nabla \cdot u I + 2 \mu \frac{1}{2} (\nabla u + \nabla
2140  u^T) \\
2141  &= \lambda \nabla \cdot u I + \mu (\nabla u + \nabla u^T)
2142  \end{split}
2143  \f]
2144 
2145  where \f$ I \f$ is identity matrix, \f$ \lambda \f$ and \f$ \mu \f$ are Lame
2146  coefficients (see ElasticityIntegrator), \f$ u, v \f$ are the trial and test
2147  functions, respectively.
2148 
2149  The parameters \f$ \alpha \f$ and \f$ \kappa \f$ determine the DG method to
2150  use (when this integrator is added to the "broken" ElasticityIntegrator):
2151 
2152  - IIPG, \f$\alpha = 0\f$,
2153  C. Dawson, S. Sun, M. Wheeler, Compatible algorithms for coupled flow and
2154  transport, Comp. Meth. Appl. Mech. Eng., 193(23-26), 2565-2580, 2004.
2155 
2156  - SIPG, \f$\alpha = -1\f$,
2157  M. Grote, A. Schneebeli, D. Schotzau, Discontinuous Galerkin Finite
2158  %Element Method for the Wave Equation, SINUM, 44(6), 2408-2431, 2006.
2159 
2160  - NIPG, \f$\alpha = 1\f$,
2161  B. Riviere, M. Wheeler, V. Girault, A Priori Error Estimates for Finite
2162  %Element Methods Based on Discontinuous Approximation Spaces for Elliptic
2163  Problems, SINUM, 39(3), 902-931, 2001.
2164 
2165  This is a '%Vector' integrator, i.e. defined for FE spaces using multiple
2166  copies of a scalar FE space.
2167  */
2169 {
2170 public:
2171  DGElasticityIntegrator(double alpha_, double kappa_)
2172  : lambda(NULL), mu(NULL), alpha(alpha_), kappa(kappa_) { }
2173 
2175  double alpha_, double kappa_)
2176  : lambda(&lambda_), mu(&mu_), alpha(alpha_), kappa(kappa_) { }
2177 
2179  virtual void AssembleFaceMatrix(const FiniteElement &el1,
2180  const FiniteElement &el2,
2182  DenseMatrix &elmat);
2183 
2184 protected:
2186  double alpha, kappa;
2187 
2188 #ifndef MFEM_THREAD_SAFE
2189  // values of all scalar basis functions for one component of u (which is a
2190  // vector) at the integration point in the reference space
2192  // values of derivatives of all scalar basis functions for one component
2193  // of u (which is a vector) at the integration point in the reference space
2195  // Adjugate of the Jacobian of the transformation: adjJ = det(J) J^{-1}
2197  // gradient of shape functions in the real (physical, not reference)
2198  // coordinates, scaled by det(J):
2199  // dshape_ps(jdof,jm) = sum_{t} adjJ(t,jm)*dshape(jdof,t)
2201  Vector nor; // nor = |weight(J_face)| n
2202  Vector nL1, nL2; // nL1 = (lambda1 * ip.weight / detJ1) nor
2203  Vector nM1, nM2; // nM1 = (mu1 * ip.weight / detJ1) nor
2204  Vector dshape1_dnM, dshape2_dnM; // dshape1_dnM = dshape1_ps . nM1
2205  // 'jmat' corresponds to the term: kappa <h^{-1} {lambda + 2 mu} [u], [v]>
2207 #endif
2208 
2209  static void AssembleBlock(
2210  const int dim, const int row_ndofs, const int col_ndofs,
2211  const int row_offset, const int col_offset,
2212  const double jmatcoef, const Vector &col_nL, const Vector &col_nM,
2213  const Vector &row_shape, const Vector &col_shape,
2214  const Vector &col_dshape_dnM, const DenseMatrix &col_dshape,
2215  DenseMatrix &elmat, DenseMatrix &jmat);
2216 };
2217 
2218 /** Integrator for the DPG form: < v, [w] > over all faces (the interface) where
2219  the trial variable v is defined on the interface and the test variable w is
2220  defined inside the elements, generally in a DG space. */
2222 {
2223 private:
2224  Vector face_shape, shape1, shape2;
2225 
2226 public:
2229  virtual void AssembleFaceMatrix(const FiniteElement &trial_face_fe,
2230  const FiniteElement &test_fe1,
2231  const FiniteElement &test_fe2,
2233  DenseMatrix &elmat);
2234 };
2235 
2236 /** Integrator for the form: < v, [w.n] > over all faces (the interface) where
2237  the trial variable v is defined on the interface and the test variable w is
2238  in an H(div)-conforming space. */
2240 {
2241 private:
2242  Vector face_shape, normal, shape1_n, shape2_n;
2243  DenseMatrix shape1, shape2;
2244 
2245 public:
2248  virtual void AssembleFaceMatrix(const FiniteElement &trial_face_fe,
2249  const FiniteElement &test_fe1,
2250  const FiniteElement &test_fe2,
2252  DenseMatrix &elmat);
2253 };
2254 
2255 /** Abstract class to serve as a base for local interpolators to be used in the
2256  DiscreteLinearOperator class. */
2258 
2259 
2260 /** Class for constructing the gradient as a DiscreteLinearOperator from an
2261  H1-conforming space to an H(curl)-conforming space. The range space can be
2262  vector L2 space as well. */
2264 {
2265 public:
2266  virtual void AssembleElementMatrix2(const FiniteElement &h1_fe,
2267  const FiniteElement &nd_fe,
2269  DenseMatrix &elmat)
2270  { nd_fe.ProjectGrad(h1_fe, Trans, elmat); }
2271 };
2272 
2273 
2274 /** Class for constructing the identity map as a DiscreteLinearOperator. This
2275  is the discrete embedding matrix when the domain space is a subspace of
2276  the range space. Otherwise, a dof projection matrix is constructed. */
2278 {
2279 public:
2280  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2281  const FiniteElement &ran_fe,
2283  DenseMatrix &elmat)
2284  { ran_fe.Project(dom_fe, Trans, elmat); }
2285 };
2286 
2287 
2288 /** Class for constructing the (local) discrete curl matrix which can be used
2289  as an integrator in a DiscreteLinearOperator object to assemble the global
2290  discrete curl matrix. */
2292 {
2293 public:
2294  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2295  const FiniteElement &ran_fe,
2297  DenseMatrix &elmat)
2298  { ran_fe.ProjectCurl(dom_fe, Trans, elmat); }
2299 };
2300 
2301 
2302 /** Class for constructing the (local) discrete divergence matrix which can
2303  be used as an integrator in a DiscreteLinearOperator object to assemble
2304  the global discrete divergence matrix.
2305 
2306  Note: Since the dofs in the L2_FECollection are nodal values, the local
2307  discrete divergence matrix (with an RT-type domain space) will depend on
2308  the transformation. On the other hand, the local matrix returned by
2309  VectorFEDivergenceIntegrator is independent of the transformation. */
2311 {
2312 public:
2313  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2314  const FiniteElement &ran_fe,
2316  DenseMatrix &elmat)
2317  { ran_fe.ProjectDiv(dom_fe, Trans, elmat); }
2318 };
2319 
2320 
2321 /** A trace face interpolator class for interpolating the normal component of
2322  the domain space, e.g. vector H1, into the range space, e.g. the trace of
2323  RT which uses FiniteElement::INTEGRAL map type. */
2325 {
2326 public:
2327  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2328  const FiniteElement &ran_fe,
2330  DenseMatrix &elmat);
2331 };
2332 
2333 /** Interpolator of a scalar coefficient multiplied by a scalar field onto
2334  another scalar field. Note that this can produce inaccurate fields unless
2335  the target is sufficiently high order. */
2337 {
2338 public:
2340 
2341  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2342  const FiniteElement &ran_fe,
2344  DenseMatrix &elmat);
2345 
2346 protected:
2348 };
2349 
2350 /** Interpolator of a scalar coefficient multiplied by a vector field onto
2351  another vector field. Note that this can produce inaccurate fields unless
2352  the target is sufficiently high order. */
2354 {
2355 public:
2357  : Q(sc) { }
2358 
2359  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2360  const FiniteElement &ran_fe,
2362  DenseMatrix &elmat);
2363 protected:
2365 };
2366 
2367 /** Interpolator of a vector coefficient multiplied by a scalar field onto
2368  another vector field. Note that this can produce inaccurate fields unless
2369  the target is sufficiently high order. */
2371 {
2372 public:
2374  : VQ(vc) { }
2375 
2376  virtual void AssembleElementMatrix2(const FiniteElement &dom_fe,
2377  const FiniteElement &ran_fe,
2379  DenseMatrix &elmat);
2380 protected:
2382 };
2383 
2384 /** Interpolator of the cross product between a vector coefficient and an
2385  H(curl)-conforming field onto an H(div)-conforming field. The range space
2386  can also be vector L2. */
2388 {
2389 public:
2391  : VQ(vc) { }
2392 
2393  virtual void AssembleElementMatrix2(const FiniteElement &nd_fe,
2394  const FiniteElement &rt_fe,
2396  DenseMatrix &elmat);
2397 protected:
2399 };
2400 
2401 /** Interpolator of the inner product between a vector coefficient and an
2402  H(div)-conforming field onto an L2-conforming field. The range space can
2403  also be H1. */
2405 {
2406 public:
2408 
2409  virtual void AssembleElementMatrix2(const FiniteElement &rt_fe,
2410  const FiniteElement &l2_fe,
2412  DenseMatrix &elmat);
2413 protected:
2415 };
2416 
2417 }
2418 
2419 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:47
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:464
virtual void AssembleElementMatrix(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &elmat)
Support for use in BilinearForm. Can be used only when appropriate.
Definition: bilininteg.hpp:184
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:857
virtual void AssembleElementMatrix(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &elmat)
Support for use in BilinearForm. Can be used only when appropriate.
Definition: bilininteg.hpp:251
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
VectorFEMassIntegrator(Coefficient &q)
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:65
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:963
int GetDim() const
Returns the reference space dimension for the finite element.
Definition: fe.hpp:116
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:619
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:83
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:590
MixedCrossCurlGradIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:960
DGTraceIntegrator(VectorCoefficient &_u, double a, double b)
Construct integrator with rho = 1.
MixedGradDivIntegrator(VectorCoefficient &vq)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
DGDiffusionIntegrator(Coefficient &q, const double s, const double k)
Integrator defining a sum of multiple Integrators.
Definition: bilininteg.hpp:152
SumIntegrator(int own_integs=1)
Definition: bilininteg.hpp:160
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
Definition: fe.cpp:38
virtual const char * FiniteElementTypeFailureMessage() const
virtual const char * FiniteElementTypeFailureMessage() const
DiffusionIntegrator(Coefficient &q)
Construct a diffusion integrator with a scalar coefficient q.
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:272
MixedScalarCrossProductIntegrator(VectorCoefficient &vq)
CurlCurlIntegrator(MatrixCoefficient &m)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:864
virtual void ProjectDiv(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const
Definition: fe.cpp:162
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:981
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:367
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:115
MixedVectorWeakDivergenceIntegrator(VectorCoefficient &dq)
virtual double ComputeFluxEnergy(const FiniteElement &fluxelem, ElementTransformation &Trans, Vector &flux, Vector *d_energy=NULL)
void CalcPhysDivShape(ElementTransformation &Trans, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in physical space at the po...
Definition: fe.cpp:59
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:750
ElasticityIntegrator(Coefficient &l, Coefficient &m)
VectorFEMassIntegrator(VectorCoefficient *_vq)
VectorDiffusionIntegrator(Coefficient &q)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:218
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:445
Integrator for (curl u, curl v) for Nedelec elements.
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
VectorCurlCurlIntegrator(Coefficient &q)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:481
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
MixedScalarMassIntegrator(Coefficient &q)
Definition: bilininteg.hpp:405
VectorMassIntegrator(Coefficient &q, const IntegrationRule *ir)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:833
alpha (q . grad u, v) using the &quot;group&quot; FE discretization
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:957
int GetOrder() const
Returns the order of the finite element.
Definition: fe.hpp:125
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:575
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
DivDivIntegrator(Coefficient &q)
int Size() const
Returns the size of the vector.
Definition: vector.hpp:113
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
MixedVectorWeakCurlIntegrator(MatrixCoefficient &mq)
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
Polynomials of order k.
Definition: fe.hpp:34
virtual const char * FiniteElementTypeFailureMessage() const
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:912
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:82
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:91
BilinearFormIntegrator(const IntegrationRule *ir=NULL)
Definition: bilininteg.hpp:25
int Space() const
Returns the type of space on each element.
Definition: fe.hpp:128
void AddIntegrator(BilinearFormIntegrator *integ)
Definition: bilininteg.hpp:162
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:811
MixedVectorCurlIntegrator(VectorCoefficient &dq)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:223
virtual double ComputeFluxEnergy(const FiniteElement &fluxelem, ElementTransformation &Trans, Vector &flux, Vector *d_energy=NULL)
Definition: bilininteg.hpp:78
MixedCrossProductIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:696
MixedScalarVectorIntegrator(VectorCoefficient &vq, bool _transpose=false, bool _cross_2d=false)
Definition: bilininteg.hpp:334
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.cpp:146
double * GetData() const
Definition: vector.hpp:121
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
MixedVectorGradientIntegrator(Coefficient &q)
(Q div u, div v) for RT elements
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:734
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:351
MixedScalarCrossGradIntegrator(VectorCoefficient &vq)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:828
virtual void ComputeElementFlux(const FiniteElement &el, ElementTransformation &Trans, Vector &u, const FiniteElement &fluxelem, Vector &flux, int with_coef=1)
Definition: bilininteg.cpp:595
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:811
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:801
void CalcPhysCurlShape(ElementTransformation &Trans, DenseMatrix &curl_shape) const
Evaluate the curl of all shape functions of a vector finite element in physical space at the point de...
Definition: fe.cpp:73
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:129
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Implements CalcDivShape methods.
Definition: fe.hpp:102
MixedCurlCurlIntegrator(MatrixCoefficient &mq)
Definition: bilininteg.hpp:886
MatrixCoefficient * MQ
VectorFEMassIntegrator(VectorCoefficient &vq)
VectorMassIntegrator()
Construct an integrator with coefficient 1.0.
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
VectorMassIntegrator(VectorCoefficient &q, int qo=0)
Construct an integrator with diagonal coefficient q.
VectorCrossProductInterpolator(VectorCoefficient &vc)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:198
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:378
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:297
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:73
MixedGradGradIntegrator(VectorCoefficient &dq)
Definition: bilininteg.hpp:796
Integrator that inverts the matrix assembled by another integrator.
Definition: bilininteg.hpp:134
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.cpp:154
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
ElasticityIntegrator(Coefficient &m, double q_l, double q_m)
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleFaceMatrix(const FiniteElement &trial_face_fe, const FiniteElement &test_fe1, const FiniteElement &test_fe2, FaceElementTransformations &Trans, DenseMatrix &elmat)
virtual const char * FiniteElementTypeFailureMessage() const
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:721
VectorInnerProductInterpolator(VectorCoefficient &vc)
MixedWeakDivCrossIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:761
virtual const char * FiniteElementTypeFailureMessage() const
virtual const char * FiniteElementTypeFailureMessage() const
int dim
Definition: ex3.cpp:47
VectorDivergenceIntegrator(Coefficient &q)
MixedCrossCurlCurlIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:923
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:98
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:986
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
MixedVectorWeakCurlIntegrator(VectorCoefficient &dq)
MixedScalarWeakCurlIntegrator(Coefficient &q)
Definition: bilininteg.hpp:647
MixedCurlCurlIntegrator(Coefficient &q)
Definition: bilininteg.hpp:882
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:502
Implements CalcDShape methods.
Definition: fe.hpp:101
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:105
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: bilininteg.hpp:61
VectorFEDivergenceIntegrator(Coefficient &q)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:907
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:764
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
VectorFECurlIntegrator(Coefficient &q)
virtual void AssembleElementMatrix2(const FiniteElement &h1_fe, const FiniteElement &nd_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual void ComputeElementFlux(const FiniteElement &el, ElementTransformation &Trans, Vector &u, const FiniteElement &fluxelem, Vector &flux, int with_coef=1)
Definition: bilininteg.hpp:72
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:280
DiffusionIntegrator()
Construct a diffusion integrator with coefficient Q = 1.
MixedVectorMassIntegrator(MatrixCoefficient &mq)
Definition: bilininteg.hpp:687
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
MixedScalarDivergenceIntegrator(Coefficient &q)
Definition: bilininteg.hpp:498
MixedVectorIntegrator(MatrixCoefficient &mq)
Definition: bilininteg.hpp:269
Implements CalcCurlShape methods.
Definition: fe.hpp:103
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:292
MixedVectorMassIntegrator(Coefficient &q)
Definition: bilininteg.hpp:683
void CalcPhysDShape(ElementTransformation &Trans, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in physical space at the poi...
Definition: fe.cpp:180
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
MixedDirectionalDerivativeIntegrator(VectorCoefficient &vq)
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
MixedCrossGradIntegrator(VectorCoefficient &vq)
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:553
ScalarProductInterpolator(Coefficient &sc)
virtual const char * FiniteElementTypeFailureMessage() const
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:651
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:944
MixedVectorCurlIntegrator(MatrixCoefficient &mq)
MixedCurlCurlIntegrator(VectorCoefficient &dq)
Definition: bilininteg.hpp:884
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:667
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:510
BoundaryMassIntegrator(Coefficient &q)
virtual void AssembleElementMatrix(const FiniteElement &, ElementTransformation &, DenseMatrix &)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:903
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:889
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect)
Perform the local action of the BilinearFormIntegrator.
Definition: bilininteg.cpp:56
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect)
Perform the local action of the BilinearFormIntegrator.
Definition: bilininteg.cpp:517
MixedWeakGradDotIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:731
ScalarVectorProductInterpolator(Coefficient &sc)
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:717
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:558
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:22
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:297
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
VectorDivergenceIntegrator(Coefficient *_q)
MixedVectorIntegrator(Coefficient &q)
Definition: bilininteg.hpp:264
DiffusionIntegrator(MatrixCoefficient &q)
Construct a diffusion integrator with a matrix coefficient q.
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:660
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
MixedGradGradIntegrator(MatrixCoefficient &mq)
Definition: bilininteg.hpp:798
virtual void AssembleElementMatrix(const FiniteElement &, ElementTransformation &, DenseMatrix &)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:855
VectorFEMassIntegrator(MatrixCoefficient &mq)
MixedVectorMassIntegrator(VectorCoefficient &dq)
Definition: bilininteg.hpp:685
MixedScalarWeakCurlCrossIntegrator(VectorCoefficient &vq)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:610
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Given a particular Finite Element computes the element matrix elmat.
Definition: bilininteg.cpp:23
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:538
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual const char * FiniteElementTypeFailureMessage() const
void CalcPhysShape(ElementTransformation &Trans, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in physical space at the point ...
Definition: fe.cpp:170
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
MixedScalarWeakDivergenceIntegrator(VectorCoefficient &vq)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:437
DGDiffusionIntegrator(const double s, const double k)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:781
LumpedIntegrator(BilinearFormIntegrator *_bfi, int _own_bfi=1)
Definition: bilininteg.hpp:123
MixedCrossGradCurlIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:997
MixedVectorDivergenceIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:534
MassIntegrator(Coefficient &q, const IntegrationRule *ir=NULL)
Construct a mass integrator with coefficient q.
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:206
MixedCrossGradGradIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:844
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:517
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:926
MixedVectorIntegrator(VectorCoefficient &dq, bool diag=true)
Definition: bilininteg.hpp:266
MixedScalarCrossCurlIntegrator(VectorCoefficient &vq)
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Assemble an element matrix.
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:39
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:583
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:743
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:774
VectorMassIntegrator(Coefficient &q, int qo=0)
VectorFEMassIntegrator(Coefficient *_q)
DGElasticityIntegrator(double alpha_, double kappa_)
MixedScalarDerivativeIntegrator(Coefficient &q)
Definition: bilininteg.hpp:425
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
VectorFEMassIntegrator(MatrixCoefficient *_mq)
MixedDotProductIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:706
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, DenseMatrix &shape)
MixedVectorCurlIntegrator(Coefficient &q)
MixedWeakCurlCrossIntegrator(VectorCoefficient &vq)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:473
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:767
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void AssembleElementMatrix2(const FiniteElement &nd_fe, const FiniteElement &rt_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:188
VectorMassIntegrator(MatrixCoefficient &q, int qo=0)
Construct an integrator with matrix coefficient q.
MixedScalarWeakCrossProductIntegrator(VectorCoefficient &vq)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void AssembleElementMatrix(const FiniteElement &, ElementTransformation &, DenseMatrix &)
Given a particular Finite Element computes the element matrix elmat.
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:286
GroupConvectionIntegrator(VectorCoefficient &q, double a=1.0)
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:373
virtual void CalcVShape(const FiniteElement &vector_fe, ElementTransformation &Trans, DenseMatrix &shape)
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:869
MixedVectorGradientIntegrator(VectorCoefficient &dq)
static void AssembleBlock(const int dim, const int row_ndofs, const int col_ndofs, const int row_offset, const int col_offset, const double jmatcoef, const Vector &col_nL, const Vector &col_nM, const Vector &row_shape, const Vector &col_shape, const Vector &col_dshape_dnM, const DenseMatrix &col_dshape, DenseMatrix &elmat, DenseMatrix &jmat)
virtual void CalcShape(const FiniteElement &scalar_fe, ElementTransformation &Trans, Vector &shape)
virtual ~SumIntegrator()
Definition: bilininteg.cpp:118
virtual const char * FiniteElementTypeFailureMessage() const
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect)
Perform the local action of the BilinearFormIntegrator.
MixedVectorWeakDivergenceIntegrator(MatrixCoefficient &mq)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual void AssembleElementMatrix2(const FiniteElement &rt_fe, const FiniteElement &l2_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:438
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:900
virtual const char * FiniteElementTypeFailureMessage() const
DGDiffusionIntegrator(MatrixCoefficient &q, const double s, const double k)
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:631
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun)
Compute element energy: (1/2) (curl u, curl u)_E.
virtual void AssembleFaceGrad(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr, const Vector &elfun, DenseMatrix &elmat)
Assemble the local action of the gradient of the NonlinearFormIntegrator resulting from a face integr...
Definition: bilininteg.hpp:66
MixedCrossCurlIntegrator(VectorCoefficient &vq)
Vector data type.
Definition: vector.hpp:41
int GetDerivType() const
Definition: fe.hpp:136
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:546
MixedVectorProductIntegrator(VectorCoefficient &vq)
Definition: bilininteg.hpp:415
VectorScalarProductInterpolator(VectorCoefficient &vc)
virtual const char * FiniteElementTypeFailureMessage() const
ConvectionIntegrator(VectorCoefficient &q, double a=1.0)
virtual void AssembleFaceMatrix(const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Trans, DenseMatrix &elmat)
MixedScalarCurlIntegrator(Coefficient &q)
Definition: bilininteg.hpp:606
virtual void AssembleElementMatrix2(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:31
virtual void CalcTestShape(const FiniteElement &test_fe, ElementTransformation &Trans, DenseMatrix &shape)
Definition: bilininteg.hpp:949
DerivativeIntegrator(Coefficient &q, int i)
DGTraceIntegrator(Coefficient &_rho, VectorCoefficient &_u, double a, double b)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:709
CurlCurlIntegrator(Coefficient &q)
Construct a bilinear form integrator for Nedelec elements.
virtual const char * FiniteElementTypeFailureMessage() const
MixedGradGradIntegrator(Coefficient &q)
Definition: bilininteg.hpp:794
MixedVectorGradientIntegrator(MatrixCoefficient &mq)
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:818
virtual void AssembleFaceMatrix(const FiniteElement &trial_face_fe, const FiniteElement &test_fe1, const FiniteElement &test_fe2, FaceElementTransformations &Trans, DenseMatrix &elmat)
int GetRangeType() const
Definition: fe.hpp:130
virtual void CalcTrialShape(const FiniteElement &trial_fe, ElementTransformation &Trans, Vector &shape)
Definition: bilininteg.hpp:522
virtual void AssembleElementMatrix2(const FiniteElement &dom_fe, const FiniteElement &ran_fe, ElementTransformation &Trans, DenseMatrix &elmat)
virtual void AssembleElementMatrix(const FiniteElement &el, ElementTransformation &Trans, DenseMatrix &elmat)
Definition: bilininteg.cpp:365
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:974
virtual const char * FiniteElementTypeFailureMessage() const
Definition: bilininteg.hpp:937
MassIntegrator(const IntegrationRule *ir=NULL)
virtual void ComputeElementFlux(const FiniteElement &el, ElementTransformation &Trans, Vector &u, const FiniteElement &fluxelem, Vector &flux, int with_coef)
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:212
virtual double ComputeFluxEnergy(const FiniteElement &fluxelem, ElementTransformation &Trans, Vector &flux, Vector *d_energy=NULL)
Definition: bilininteg.cpp:653
Integrator for (Q u, v) for VectorFiniteElements.
MixedVectorWeakCurlIntegrator(Coefficient &q)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:847
TransposeIntegrator(BilinearFormIntegrator *_bfi, int _own_bfi=1)
Definition: bilininteg.hpp:95
MixedDivGradIntegrator(VectorCoefficient &vq)
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:429
virtual int GetIntegrationOrder(const FiniteElement &trial_fe, const FiniteElement &test_fe, ElementTransformation &Trans)
Definition: bilininteg.hpp:626
MixedScalarIntegrator(Coefficient &q)
Definition: bilininteg.hpp:196
DGElasticityIntegrator(Coefficient &lambda_, Coefficient &mu_, double alpha_, double kappa_)
Class for integrating (Q D_i(u), v); u and v are scalars.
InverseIntegrator(BilinearFormIntegrator *integ, int own_integ=1)
Definition: bilininteg.hpp:141
virtual bool VerifyFiniteElementTypes(const FiniteElement &trial_fe, const FiniteElement &test_fe) const
Definition: bilininteg.hpp:338
alpha (q . grad u, v)