MFEM  v4.5.2
Finite element discretization library
coefficient.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_COEFFICIENT
13 #define MFEM_COEFFICIENT
14 
15 #include <functional>
16 
17 #include "../config/config.hpp"
18 #include "../linalg/linalg.hpp"
19 #include "intrules.hpp"
20 #include "eltrans.hpp"
21 
22 namespace mfem
23 {
24 
25 class Mesh;
26 class QuadratureSpaceBase;
27 class QuadratureFunction;
28 
29 #ifdef MFEM_USE_MPI
30 class ParMesh;
31 #endif
32 
33 
34 /** @brief Base class Coefficients that optionally depend on space and time.
35  These are used by the BilinearFormIntegrator, LinearFormIntegrator, and
36  NonlinearFormIntegrator classes to represent the physical coefficients in
37  the PDEs that are being discretized. This class can also be used in a more
38  general way to represent functions that don't necessarily belong to a FE
39  space, e.g., to project onto GridFunctions to use as initial conditions,
40  exact solutions, etc. See, e.g., ex4 or ex22 for these uses. */
42 {
43 protected:
44  double time;
45 
46 public:
47  Coefficient() { time = 0.; }
48 
49  /// Set the time for time dependent coefficients
50  virtual void SetTime(double t) { time = t; }
51 
52  /// Get the time for time dependent coefficients
53  double GetTime() { return time; }
54 
55  /** @brief Evaluate the coefficient in the element described by @a T at the
56  point @a ip. */
57  /** @note When this method is called, the caller must make sure that the
58  IntegrationPoint associated with @a T is the same as @a ip. This can be
59  achieved by calling T.SetIntPoint(&ip). */
60  virtual double Eval(ElementTransformation &T,
61  const IntegrationPoint &ip) = 0;
62 
63  /** @brief Evaluate the coefficient in the element described by @a T at the
64  point @a ip at time @a t. */
65  /** @note When this method is called, the caller must make sure that the
66  IntegrationPoint associated with @a T is the same as @a ip. This can be
67  achieved by calling T.SetIntPoint(&ip). */
69  const IntegrationPoint &ip, double t)
70  {
71  SetTime(t);
72  return Eval(T, ip);
73  }
74 
75  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
76  /// the quadrature points.
77  virtual void Project(QuadratureFunction &qf);
78 
79  virtual ~Coefficient() { }
80 };
81 
82 
83 /// A coefficient that is constant across space and time
85 {
86 public:
87  double constant;
88 
89  /// c is value of constant function
90  explicit ConstantCoefficient(double c = 1.0) { constant=c; }
91 
92  /// Evaluate the coefficient at @a ip.
93  virtual double Eval(ElementTransformation &T,
94  const IntegrationPoint &ip)
95  { return (constant); }
96 
97  /// Fill the QuadratureFunction @a qf with the constant value.
98  void Project(QuadratureFunction &qf);
99 };
100 
101 /** @brief A piecewise constant coefficient with the constants keyed
102  off the element attribute numbers. */
104 {
105 private:
106  Vector constants;
107 
108 public:
109 
110  /// Constructs a piecewise constant coefficient in NumOfSubD subdomains
111  explicit PWConstCoefficient(int NumOfSubD = 0) : constants(NumOfSubD)
112  { constants = 0.0; }
113 
114  /// Construct the constant coefficient using a vector of constants.
115  /** @a c should be a vector defined by attributes, so for region with
116  attribute @a i @a c[i-1] is the coefficient in that region */
118  { constants.SetSize(c.Size()); constants=c; }
119 
120  /// Update the constants with vector @a c.
121  void UpdateConstants(Vector &c) { constants.SetSize(c.Size()); constants=c; }
122 
123  /// Return a reference to the i-th constant
124  double &operator()(int i) { return constants(i-1); }
125 
126  /// Set the constants for all attributes to constant @a c.
127  void operator=(double c) { constants = c; }
128 
129  /// Returns the number of constants representing different attributes.
130  int GetNConst() { return constants.Size(); }
131 
132  /// Evaluate the coefficient.
133  virtual double Eval(ElementTransformation &T,
134  const IntegrationPoint &ip);
135 };
136 
137 /** @brief A piecewise coefficient with the pieces keyed off the element
138  attribute numbers.
139 
140  A value of zero will be returned for any missing attribute numbers.
141 
142  This object will not assume ownership of any Coefficient objects
143  passed to it. Consequently, the caller must ensure that the
144  individual Coefficient objects are not deleted while this
145  PWCoefficient is still in use.
146 
147  \note The keys may either be domain attribute numbers or boundary
148  attribute numbers. If the PWCoefficient is used with a domain
149  integrator the keys are assumed to be domain attribute
150  numbers. Similarly, if the PWCoefficient is used with a boundary
151  integrator the keys are assumed to be boundary attribute numbers.
152 */
154 {
155 private:
156  /** Internal data structure to store pointers to the appropriate
157  coefficients for different regions of the mesh. The keys used
158  in the map are the mesh attribute numbers (either element
159  attribute or boundary element attribute depending upon
160  context). The values returned for any missing attributes will
161  be zero. The coefficient pointers may be NULL in which case a
162  value of zero is returned.
163 
164  The Coefficient objects contained in this map are NOT owned by
165  this PWCoefficient object. This means that they will not be
166  deleted when this object is deleted also the caller must ensure
167  that the various Coefficient objects are not deleted while this
168  PWCoefficient is still needed.
169  */
170  std::map<int, Coefficient*> pieces;
171 
172  /** Convenience function to check for compatible array lengths,
173  loop over the arrays, and add their attribute/Coefficient pairs
174  to the internal data structure.
175  */
176  void InitMap(const Array<int> & attr,
177  const Array<Coefficient*> & coefs);
178 
179 public:
180 
181  /// Constructs a piecewise coefficient
182  explicit PWCoefficient() {}
183 
184  /// Construct the coefficient using arrays describing the pieces
185  /** \param attr - an array of attribute numbers for each piece
186  \param coefs - the corresponding array of Coefficient pointers
187  Any missing attributes or NULL coefficient pointers will result in a
188  value of zero being returned for that attribute.
189 
190  \note Ownership of the Coefficient objects will NOT be
191  transferred to this object.
192  */
193  PWCoefficient(const Array<int> & attr,
194  const Array<Coefficient*> & coefs)
195  { InitMap(attr, coefs); }
196 
197  /// Set the time for time dependent coefficients
198  virtual void SetTime(double t);
199 
200  /// Replace a set of coefficients
201  void UpdateCoefficients(const Array<int> & attr,
202  const Array<Coefficient*> & coefs)
203  { InitMap(attr, coefs); }
204 
205  /// Replace a single Coefficient for a particular attribute
206  void UpdateCoefficient(int attr, Coefficient & coef)
207  { pieces[attr] = &coef; }
208 
209  /// Remove a single Coefficient for a particular attribute
210  void ZeroCoefficient(int attr)
211  { pieces.erase(attr); }
212 
213  /// Evaluate the coefficient.
214  virtual double Eval(ElementTransformation &T,
215  const IntegrationPoint &ip);
216 };
217 
218 /// A general function coefficient
220 {
221 protected:
222  std::function<double(const Vector &)> Function;
223  std::function<double(const Vector &, double)> TDFunction;
224 
225 public:
226  /// Define a time-independent coefficient from a std function
227  /** \param F time-independent std::function */
228  FunctionCoefficient(std::function<double(const Vector &)> F)
229  : Function(std::move(F))
230  { }
231 
232  /// Define a time-dependent coefficient from a std function
233  /** \param TDF time-dependent function */
234  FunctionCoefficient(std::function<double(const Vector &, double)> TDF)
235  : TDFunction(std::move(TDF))
236  { }
237 
238  /// (DEPRECATED) Define a time-independent coefficient from a C-function
239  /** @deprecated Use the method where the C-function, @a f, uses a const
240  Vector argument instead of Vector. */
241  MFEM_DEPRECATED FunctionCoefficient(double (*f)(Vector &))
242  {
243  Function = reinterpret_cast<double(*)(const Vector&)>(f);
244  TDFunction = NULL;
245  }
246 
247  /// (DEPRECATED) Define a time-dependent coefficient from a C-function
248  /** @deprecated Use the method where the C-function, @a tdf, uses a const
249  Vector argument instead of Vector. */
250  MFEM_DEPRECATED FunctionCoefficient(double (*tdf)(Vector &, double))
251  {
252  Function = NULL;
253  TDFunction = reinterpret_cast<double(*)(const Vector&,double)>(tdf);
254  }
255 
256  /// Evaluate the coefficient at @a ip.
257  virtual double Eval(ElementTransformation &T,
258  const IntegrationPoint &ip);
259 };
260 
261 class GridFunction;
262 
263 /// Coefficient defined by a GridFunction. This coefficient is mesh dependent.
265 {
266 private:
267  const GridFunction *GridF;
268  int Component;
269 
270 public:
271  GridFunctionCoefficient() : GridF(NULL), Component(1) { }
272  /** Construct GridFunctionCoefficient from a given GridFunction, and
273  optionally specify a component to use if it is a vector GridFunction. */
274  GridFunctionCoefficient (const GridFunction *gf, int comp = 1)
275  { GridF = gf; Component = comp; }
276 
277  /// Set the internal GridFunction
278  void SetGridFunction(const GridFunction *gf) { GridF = gf; }
279 
280  /// Get the internal GridFunction
281  const GridFunction * GetGridFunction() const { return GridF; }
282 
283  /// Evaluate the coefficient at @a ip.
284  virtual double Eval(ElementTransformation &T,
285  const IntegrationPoint &ip);
286 
287  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
288  /// the quadrature points.
289  ///
290  /// This function uses the efficient QuadratureFunction::ProjectGridFunction
291  /// to fill the QuadratureFunction.
292  virtual void Project(QuadratureFunction &qf);
293 };
294 
295 
296 /** @brief A coefficient that depends on 1 or 2 parent coefficients and a
297  transformation rule represented by a C-function.
298 
299  \f$ C(x,t) = T(Q1(x,t)) \f$ or \f$ C(x,t) = T(Q1(x,t), Q2(x,t)) \f$
300 
301  where T is the transformation rule, and Q1/Q2 are the parent coefficients.*/
303 {
304 private:
305  Coefficient * Q1;
306  Coefficient * Q2;
307  double (*Transform1)(double);
308  double (*Transform2)(double,double);
309 
310 public:
311  TransformedCoefficient (Coefficient * q,double (*F)(double))
312  : Q1(q), Transform1(F) { Q2 = 0; Transform2 = 0; }
314  double (*F)(double,double))
315  : Q1(q1), Q2(q2), Transform2(F) { Transform1 = 0; }
316 
317  /// Set the time for internally stored coefficients
318  void SetTime(double t);
319 
320  /// Evaluate the coefficient at @a ip.
321  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
322 };
323 
324 /** @brief Delta function coefficient optionally multiplied by a weight
325  coefficient and a scaled time dependent C-function.
326 
327  \f$ F(x,t) = w(x,t) s T(t) d(x - xc) \f$
328 
329  where w is the optional weight coefficient, @a s is a scale factor
330  T is an optional time-dependent function and d is a delta function.
331 
332  WARNING this cannot be used as a normal coefficient. The usual Eval
333  method is disabled. */
335 {
336 protected:
337  double center[3], scale, tol;
339  int sdim;
340  double (*tdf)(double);
341 
342 public:
343 
344  /// Construct a unit delta function centered at (0.0,0.0,0.0)
346  {
347  center[0] = center[1] = center[2] = 0.; scale = 1.; tol = 1e-12;
348  weight = NULL; sdim = 0; tdf = NULL;
349  }
350 
351  /// Construct a delta function scaled by @a s and centered at (x,0.0,0.0)
352  DeltaCoefficient(double x, double s)
353  {
354  center[0] = x; center[1] = 0.; center[2] = 0.; scale = s; tol = 1e-12;
355  weight = NULL; sdim = 1; tdf = NULL;
356  }
357 
358  /// Construct a delta function scaled by @a s and centered at (x,y,0.0)
359  DeltaCoefficient(double x, double y, double s)
360  {
361  center[0] = x; center[1] = y; center[2] = 0.; scale = s; tol = 1e-12;
362  weight = NULL; sdim = 2; tdf = NULL;
363  }
364 
365  /// Construct a delta function scaled by @a s and centered at (x,y,z)
366  DeltaCoefficient(double x, double y, double z, double s)
367  {
368  center[0] = x; center[1] = y; center[2] = z; scale = s; tol = 1e-12;
369  weight = NULL; sdim = 3; tdf = NULL;
370  }
371 
372  /// Set the time for internally stored coefficients
373  void SetTime(double t);
374 
375  /// Set the center location of the delta function.
376  void SetDeltaCenter(const Vector& center);
377 
378  /// Set the scale value multiplying the delta function.
379  void SetScale(double s_) { scale = s_; }
380 
381  /// Set a time-dependent function that multiplies the Scale().
382  void SetFunction(double (*f)(double)) { tdf = f; }
383 
384  /** @brief Set the tolerance used during projection onto GridFunction to
385  identify the Mesh vertex where the Center() of the delta function
386  lies. (default 1e-12)*/
387  void SetTol(double tol_) { tol = tol_; }
388 
389  /// Set a weight Coefficient that multiplies the DeltaCoefficient.
390  /** The weight Coefficient multiplies the value returned by EvalDelta() but
391  not the value returned by Scale().
392  The weight Coefficient is also used as the L2-weight function when
393  projecting the DeltaCoefficient onto a GridFunction, so that the weighted
394  integral of the projection is exactly equal to the Scale(). */
395  void SetWeight(Coefficient *w) { weight = w; }
396 
397  /// Return a pointer to a c-array representing the center of the delta
398  /// function.
399  const double *Center() { return center; }
400 
401  /** @brief Return the scale factor times the optional time dependent
402  function. Returns \f$ s T(t) \f$ with \f$ T(t) = 1 \f$ when
403  not set by the user. */
404  double Scale() { return tdf ? (*tdf)(GetTime())*scale : scale; }
405 
406  /// Return the tolerance used to identify the mesh vertices
407  double Tol() { return tol; }
408 
409  /// See SetWeight() for description of the weight Coefficient.
410  Coefficient *Weight() { return weight; }
411 
412  /// Write the center of the delta function into @a center.
414 
415  /// The value of the function assuming we are evaluating at the delta center.
416  virtual double EvalDelta(ElementTransformation &T, const IntegrationPoint &ip);
417  /** @brief A DeltaFunction cannot be evaluated. Calling this method will
418  cause an MFEM error, terminating the application. */
419  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
420  { mfem_error("DeltaCoefficient::Eval"); return 0.; }
421  virtual ~DeltaCoefficient() { delete weight; }
422 };
423 
424 /** @brief Derived coefficient that takes the value of the parent coefficient
425  for the active attributes and is zero otherwise. */
427 {
428 private:
429  Coefficient *c;
430  Array<int> active_attr;
431 
432 public:
433  /** @brief Construct with a parent coefficient and an array with
434  ones marking the attributes on which this coefficient should be
435  active. */
437  { c = &c_; attr.Copy(active_attr); }
438 
439  /// Set the time for internally stored coefficients
440  void SetTime(double t);
441 
442  /// Evaluate the coefficient at @a ip.
443  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
444  { return active_attr[T.Attribute-1] ? c->Eval(T, ip, GetTime()) : 0.0; }
445 };
446 
447 /// Base class for vector Coefficients that optionally depend on time and space.
449 {
450 protected:
451  int vdim;
452  double time;
453 
454 public:
455  /// Initialize the VectorCoefficient with vector dimension @a vd.
456  VectorCoefficient(int vd) { vdim = vd; time = 0.; }
457 
458  /// Set the time for time dependent coefficients
459  virtual void SetTime(double t) { time = t; }
460 
461  /// Get the time for time dependent coefficients
462  double GetTime() { return time; }
463 
464  /// Returns dimension of the vector.
465  int GetVDim() { return vdim; }
466 
467  /** @brief Evaluate the vector coefficient in the element described by @a T
468  at the point @a ip, storing the result in @a V. */
469  /** @note When this method is called, the caller must make sure that the
470  IntegrationPoint associated with @a T is the same as @a ip. This can be
471  achieved by calling T.SetIntPoint(&ip). */
472  virtual void Eval(Vector &V, ElementTransformation &T,
473  const IntegrationPoint &ip) = 0;
474 
475  /** @brief Evaluate the vector coefficient in the element described by @a T
476  at all points of @a ir, storing the result in @a M. */
477  /** The dimensions of @a M are GetVDim() by ir.GetNPoints() and they must be
478  set by the implementation of this method.
479 
480  The general implementation provided by the base class (using the Eval
481  method for one IntegrationPoint at a time) can be overloaded for more
482  efficient implementation.
483 
484  @note The IntegrationPoint associated with @a T is not used, and this
485  method will generally modify this IntegrationPoint associated with @a T.
486  */
487  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
488  const IntegrationRule &ir);
489 
490  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
491  /// the quadrature points.
492  ///
493  /// The @a vdim of the VectorCoefficient should be equal to the @a vdim of
494  /// the QuadratureFunction.
495  virtual void Project(QuadratureFunction &qf);
496 
497  virtual ~VectorCoefficient() { }
498 };
499 
500 
501 /// Vector coefficient that is constant in space and time.
503 {
504 private:
505  Vector vec;
506 public:
507  /// Construct the coefficient with constant vector @a v.
509  : VectorCoefficient(v.Size()), vec(v) { }
511 
512  /// Evaluate the vector coefficient at @a ip.
513  virtual void Eval(Vector &V, ElementTransformation &T,
514  const IntegrationPoint &ip) { V = vec; }
515 
516  /// Return a reference to the constant vector in this class.
517  const Vector& GetVec() const { return vec; }
518 };
519 
520 /** @brief A piecewise vector-valued coefficient with the pieces keyed off the
521  element attribute numbers.
522 
523  A value of zero will be returned for any missing attribute numbers.
524 
525  This object will not assume ownership of any VectorCoefficient
526  objects passed to it. Consequently, the caller must ensure that
527  the individual VectorCoefficient objects are not deleted while
528  this PWVectorCoefficient is still in use.
529 
530  \note The keys may either be domain attribute numbers or boundary
531  attribute numbers. If the PWVectorCoefficient is used with a
532  domain integrator the keys are assumed to be domain attribute
533  numbers. Similarly, if the PWVectorCoefficient is used with a
534  boundary integrator the keys are assumed to be boundary attribute
535  numbers.
536 */
538 {
539 private:
540  /** Internal data structure to store pointers to the appropriate
541  coefficients for different regions of the mesh. The keys used
542  in the map are the mesh attribute numbers (either element
543  attribute or boundary element attribute depending upon
544  context). The values returned for any missing attributes will
545  be zero. The coefficient pointers may be NULL in which case a
546  value of zero is returned.
547 
548  The VectorCoefficient objects contained in this map are NOT
549  owned by this PWVectorCoefficient object. This means that they
550  will not be deleted when this object is deleted also the caller
551  must ensure that the various VectorCoefficient objects are not
552  deleted while this PWVectorCoefficient is still needed.
553  */
554  std::map<int, VectorCoefficient*> pieces;
555 
556  /** Convenience function to check for compatible array lengths,
557  loop over the arrays, and add their attribute/VectorCoefficient
558  pairs to the internal data structure.
559  */
560  void InitMap(const Array<int> & attr,
561  const Array<VectorCoefficient*> & coefs);
562 
563 public:
564 
565  /// Constructs a piecewise vector coefficient of dimension vd
566  explicit PWVectorCoefficient(int vd): VectorCoefficient(vd) {}
567 
568  /// Construct the coefficient using arrays describing the pieces
569  /** \param vd - dimension of the vector-valued result
570  \param attr - an array of attribute numbers for each piece
571  \param coefs - the corresponding array of VectorCoefficient pointers
572  Any missing attributes or NULL coefficient pointers will result in a
573  zero vector being returned for that attribute.
574 
575  \note Ownership of the VectorCoefficient objects will NOT be
576  transferred to this object.
577  */
578  PWVectorCoefficient(int vd, const Array<int> & attr,
579  const Array<VectorCoefficient*> & coefs)
580  : VectorCoefficient(vd) { InitMap(attr, coefs); }
581 
582  /// Set the time for time dependent coefficients
583  virtual void SetTime(double t);
584 
585  /// Replace a set of coefficients
586  void UpdateCoefficients(const Array<int> & attr,
587  const Array<VectorCoefficient*> & coefs)
588  { InitMap(attr, coefs); }
589 
590  /// Replace a single Coefficient for a particular attribute
591  void UpdateCoefficient(int attr, VectorCoefficient & coef);
592 
593  /// Remove a single VectorCoefficient for a particular attribute
594  void ZeroCoefficient(int attr)
595  { pieces.erase(attr); }
596 
597  /// Evaluate the coefficient.
598  virtual void Eval(Vector &V, ElementTransformation &T,
599  const IntegrationPoint &ip);
601 };
602 
603 /// A general vector function coefficient
605 {
606 private:
607  std::function<void(const Vector &, Vector &)> Function;
608  std::function<void(const Vector &, double, Vector &)> TDFunction;
609  Coefficient *Q;
610 
611 public:
612  /// Define a time-independent vector coefficient from a std function
613  /** \param dim - the size of the vector
614  \param F - time-independent function
615  \param q - optional scalar Coefficient to scale the vector coefficient */
617  std::function<void(const Vector &, Vector &)> F,
618  Coefficient *q = nullptr)
619  : VectorCoefficient(dim), Function(std::move(F)), Q(q)
620  { }
621 
622  /// Define a time-dependent vector coefficient from a std function
623  /** \param dim - the size of the vector
624  \param TDF - time-dependent function
625  \param q - optional scalar Coefficient to scale the vector coefficient */
627  std::function<void(const Vector &, double, Vector &)> TDF,
628  Coefficient *q = nullptr)
629  : VectorCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
630  { }
631 
633  /// Evaluate the vector coefficient at @a ip.
634  virtual void Eval(Vector &V, ElementTransformation &T,
635  const IntegrationPoint &ip);
636 
638 };
639 
640 /** @brief Vector coefficient defined by an array of scalar coefficients.
641  Coefficients that are not set will evaluate to zero in the vector. This
642  object takes ownership of the array of coefficients inside it and deletes
643  them at object destruction. */
645 {
646 private:
647  Array<Coefficient*> Coeff;
648  Array<bool> ownCoeff;
649 
650 public:
651  /** @brief Construct vector of dim coefficients. The actual coefficients
652  still need to be added with Set(). */
653  explicit VectorArrayCoefficient(int dim);
654 
655  /// Set the time for internally stored coefficients
656  void SetTime(double t);
657 
658  /// Returns i'th coefficient.
659  Coefficient* GetCoeff(int i) { return Coeff[i]; }
660 
661  /// Returns the entire array of coefficients.
662  Coefficient **GetCoeffs() { return Coeff; }
663 
664  /// Sets coefficient in the vector.
665  void Set(int i, Coefficient *c, bool own=true);
666 
667  /// Evaluates i'th component of the vector of coefficients and returns the
668  /// value.
669  double Eval(int i, ElementTransformation &T, const IntegrationPoint &ip)
670  { return Coeff[i] ? Coeff[i]->Eval(T, ip, GetTime()) : 0.0; }
671 
673  /** @brief Evaluate the coefficient. Each element of vector V comes from the
674  associated array of scalar coefficients. */
675  virtual void Eval(Vector &V, ElementTransformation &T,
676  const IntegrationPoint &ip);
677 
678  /// Destroys vector coefficient.
679  virtual ~VectorArrayCoefficient();
680 };
681 
682 /// Vector coefficient defined by a vector GridFunction
684 {
685 protected:
687 
688 public:
689  /** @brief Construct an empty coefficient. Calling Eval() before the grid
690  function is set will cause a segfault. */
692 
693  /** @brief Construct the coefficient with grid function @a gf. The
694  grid function is not owned by the coefficient. */
696 
697  /** @brief Set the grid function for this coefficient. Also sets the Vector
698  dimension to match that of the @a gf. */
699  void SetGridFunction(const GridFunction *gf);
700 
701  /// Returns a pointer to the grid function in this Coefficient
702  const GridFunction * GetGridFunction() const { return GridFunc; }
703 
704  /// Evaluate the vector coefficient at @a ip.
705  virtual void Eval(Vector &V, ElementTransformation &T,
706  const IntegrationPoint &ip);
707 
708  /** @brief Evaluate the vector coefficients at all of the locations in the
709  integration rule and write the vectors into the columns of matrix @a
710  M. */
711  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
712  const IntegrationRule &ir);
713 
714  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
715  /// the quadrature points.
716  ///
717  /// This function uses the efficient QuadratureFunction::ProjectGridFunction
718  /// to fill the QuadratureFunction.
719  virtual void Project(QuadratureFunction &qf);
720 
722 };
723 
724 /// Vector coefficient defined as the Gradient of a scalar GridFunction
726 {
727 protected:
729 
730 public:
731 
732  /** @brief Construct the coefficient with a scalar grid function @a gf. The
733  grid function is not owned by the coefficient. */
735 
736  ///Set the scalar grid function.
737  void SetGridFunction(const GridFunction *gf);
738 
739  ///Get the scalar grid function.
740  const GridFunction * GetGridFunction() const { return GridFunc; }
741 
742  /// Evaluate the gradient vector coefficient at @a ip.
743  virtual void Eval(Vector &V, ElementTransformation &T,
744  const IntegrationPoint &ip);
745 
746  /** @brief Evaluate the gradient vector coefficient at all of the locations
747  in the integration rule and write the vectors into columns of matrix @a
748  M. */
749  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
750  const IntegrationRule &ir);
751 
753 };
754 
755 /// Vector coefficient defined as the Curl of a vector GridFunction
757 {
758 protected:
760 
761 public:
762  /** @brief Construct the coefficient with a vector grid function @a gf. The
763  grid function is not owned by the coefficient. */
765 
766  /// Set the vector grid function.
767  void SetGridFunction(const GridFunction *gf);
768 
769  /// Get the vector grid function.
770  const GridFunction * GetGridFunction() const { return GridFunc; }
771 
773  /// Evaluate the vector curl coefficient at @a ip.
774  virtual void Eval(Vector &V, ElementTransformation &T,
775  const IntegrationPoint &ip);
776 
778 };
779 
780 /// Scalar coefficient defined as the Divergence of a vector GridFunction
782 {
783 protected:
785 
786 public:
787  /** @brief Construct the coefficient with a vector grid function @a gf. The
788  grid function is not owned by the coefficient. */
790 
791  /// Set the vector grid function.
792  void SetGridFunction(const GridFunction *gf) { GridFunc = gf; }
793 
794  /// Get the vector grid function.
795  const GridFunction * GetGridFunction() const { return GridFunc; }
796 
797  /// Evaluate the scalar divergence coefficient at @a ip.
798  virtual double Eval(ElementTransformation &T,
799  const IntegrationPoint &ip);
800 
802 };
803 
804 /** @brief Vector coefficient defined by a scalar DeltaCoefficient and a
805  constant vector direction.
806 
807  WARNING this cannot be used as a normal coefficient. The usual Eval method
808  is disabled. */
810 {
811 protected:
814 
815 public:
816  /// Construct with a vector of dimension @a vdim_.
818  : VectorCoefficient(vdim_), dir(vdim_), d() { }
819 
820  /** @brief Construct with a Vector object representing the direction and a
821  unit delta function centered at (0.0,0.0,0.0) */
823  : VectorCoefficient(dir_.Size()), dir(dir_), d() { }
824 
825  /** @brief Construct with a Vector object representing the direction and a
826  delta function scaled by @a s and centered at (x,0.0,0.0) */
827  VectorDeltaCoefficient(const Vector& dir_, double x, double s)
828  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,s) { }
829 
830  /** @brief Construct with a Vector object representing the direction and a
831  delta function scaled by @a s and centered at (x,y,0.0) */
832  VectorDeltaCoefficient(const Vector& dir_, double x, double y, double s)
833  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,y,s) { }
834 
835  /** @brief Construct with a Vector object representing the direction and a
836  delta function scaled by @a s and centered at (x,y,z) */
837  VectorDeltaCoefficient(const Vector& dir_, double x, double y, double z,
838  double s)
839  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,y,z,s) { }
840 
841  /// Set the time for internally stored coefficients
842  void SetTime(double t);
843 
844  /// Replace the associated DeltaCoefficient with a new DeltaCoefficient.
845  /** The new DeltaCoefficient cannot have a specified weight Coefficient, i.e.
846  DeltaCoefficient::Weight() should return NULL. */
847  void SetDeltaCoefficient(const DeltaCoefficient& d_) { d = d_; }
848 
849  /// Return the associated scalar DeltaCoefficient.
851 
852  void SetScale(double s) { d.SetScale(s); }
853  void SetDirection(const Vector& d_);
854 
855  void SetDeltaCenter(const Vector& center) { d.SetDeltaCenter(center); }
856  void GetDeltaCenter(Vector& center) { d.GetDeltaCenter(center); }
857 
858  /** @brief Return the specified direction vector multiplied by the value
859  returned by DeltaCoefficient::EvalDelta() of the associated scalar
860  DeltaCoefficient. */
861  virtual void EvalDelta(Vector &V, ElementTransformation &T,
862  const IntegrationPoint &ip);
863 
865  /** @brief A VectorDeltaFunction cannot be evaluated. Calling this method
866  will cause an MFEM error, terminating the application. */
867  virtual void Eval(Vector &V, ElementTransformation &T,
868  const IntegrationPoint &ip)
869  { mfem_error("VectorDeltaCoefficient::Eval"); }
871 };
872 
873 /** @brief Derived vector coefficient that has the value of the parent vector
874  where it is active and is zero otherwise. */
876 {
877 private:
879  Array<int> active_attr;
880 
881 public:
882  /** @brief Construct with a parent vector coefficient and an array of zeros
883  and ones representing the attributes for which this coefficient should be
884  active. */
886  : VectorCoefficient(vc.GetVDim())
887  { c = &vc; attr.Copy(active_attr); }
888 
889  /// Set the time for internally stored coefficients
890  void SetTime(double t);
891 
892  /// Evaluate the vector coefficient at @a ip.
893  virtual void Eval(Vector &V, ElementTransformation &T,
894  const IntegrationPoint &ip);
895 
896  /** @brief Evaluate the vector coefficient at all of the locations in the
897  integration rule and write the vectors into the columns of matrix @a
898  M. */
899  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
900  const IntegrationRule &ir);
901 };
902 
904 
905 /// Base class for Matrix Coefficients that optionally depend on time and space.
907 {
908 protected:
909  int height, width;
910  double time;
911  bool symmetric; // deprecated
912 
913 public:
914  /// Construct a dim x dim matrix coefficient.
915  explicit MatrixCoefficient(int dim, bool symm=false)
916  { height = width = dim; time = 0.; symmetric = symm; }
917 
918  /// Construct a h x w matrix coefficient.
919  MatrixCoefficient(int h, int w, bool symm=false) :
920  height(h), width(w), time(0.), symmetric(symm) { }
921 
922  /// Set the time for time dependent coefficients
923  virtual void SetTime(double t) { time = t; }
924 
925  /// Get the time for time dependent coefficients
926  double GetTime() { return time; }
927 
928  /// Get the height of the matrix.
929  int GetHeight() const { return height; }
930 
931  /// Get the width of the matrix.
932  int GetWidth() const { return width; }
933 
934  /// For backward compatibility get the width of the matrix.
935  int GetVDim() const { return width; }
936 
937  /** @deprecated Use SymmetricMatrixCoefficient instead */
938  bool IsSymmetric() const { return symmetric; }
939 
940  /** @brief Evaluate the matrix coefficient in the element described by @a T
941  at the point @a ip, storing the result in @a K. */
942  /** @note When this method is called, the caller must make sure that the
943  IntegrationPoint associated with @a T is the same as @a ip. This can be
944  achieved by calling T.SetIntPoint(&ip). */
945  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
946  const IntegrationPoint &ip) = 0;
947 
948  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
949  /// the quadrature points. The matrix will be transposed or not according to
950  /// the boolean argument @a transpose.
951  ///
952  /// The @a vdim of the QuadratureFunction should be equal to the height times
953  /// the width of the matrix.
954  virtual void Project(QuadratureFunction &qf, bool transpose=false);
955 
956  /// (DEPRECATED) Evaluate a symmetric matrix coefficient.
957  /** @brief Evaluate the upper triangular entries of the matrix coefficient
958  in the symmetric case, similarly to Eval. Matrix entry (i,j) is stored
959  in K[j - i + os_i] for 0 <= i <= j < width, os_0 = 0,
960  os_{i+1} = os_i + width - i. That is, K = {M(0,0), ..., M(0,w-1),
961  M(1,1), ..., M(1,w-1), ..., M(w-1,w-1) with w = width.
962  @deprecated Use Eval() instead. */
964  const IntegrationPoint &ip)
965  { mfem_error("MatrixCoefficient::EvalSymmetric"); }
966 
967  virtual ~MatrixCoefficient() { }
968 };
969 
970 
971 /// A matrix coefficient that is constant in space and time.
973 {
974 private:
975  DenseMatrix mat;
976 public:
977  ///Construct using matrix @a m for the constant.
979  : MatrixCoefficient(m.Height(), m.Width()), mat(m) { }
981  /// Evaluate the matrix coefficient at @a ip.
983  const IntegrationPoint &ip) { M = mat; }
984  /// Return a reference to the constant matrix.
985  const DenseMatrix& GetMatrix() { return mat; }
986 };
987 
988 
989 /** @brief A piecewise matrix-valued coefficient with the pieces keyed off the
990  element attribute numbers.
991 
992  A value of zero will be returned for any missing attribute numbers.
993 
994  This object will not assume ownership of any MatrixCoefficient
995  objects passed to it. Consequently, the caller must ensure that
996  the individual MatrixCoefficient objects are not deleted while
997  this PWMatrixCoefficient is still in use.
998 
999  \note The keys may either be domain attribute numbers or boundary
1000  attribute numbers. If the PWMatrixCoefficient is used with a
1001  domain integrator the keys are assumed to be domain attribute
1002  numbers. Similarly, if the PWMatrixCoefficient is used with a
1003  boundary integrator the keys are assumed to be boundary attribute
1004  numbers.
1005 */
1007 {
1008 private:
1009  /** Internal data structure to store pointers to the appropriate
1010  coefficients for different regions of the mesh. The keys used
1011  in the map are the mesh attribute numbers (either element
1012  attribute or boundary element attribute depending upon
1013  context). The values returned for any missing attributes will
1014  be zero. The coefficient pointers may be NULL in which case a
1015  value of zero is returned.
1016 
1017  The MatrixCoefficient objects contained in this map are NOT
1018  owned by this PWMatrixCoefficient object. This means that they
1019  will not be deleted when this object is deleted also the caller
1020  must ensure that the various MatrixCoefficient objects are not
1021  deleted while this PWMatrixCoefficient is still needed.
1022  */
1023  std::map<int, MatrixCoefficient*> pieces;
1024 
1025  /** Convenience function to check for compatible array lengths,
1026  loop over the arrays, and add their attribute/MatrixCoefficient
1027  pairs to the internal data structure.
1028  */
1029  void InitMap(const Array<int> & attr,
1030  const Array<MatrixCoefficient*> & coefs);
1031 
1032 public:
1033 
1034  /// Constructs a piecewise matrix coefficient of dimension dim by dim
1035  explicit PWMatrixCoefficient(int dim, bool symm = false)
1036  : MatrixCoefficient(dim, symm) {}
1037 
1038  /// Constructs a piecewise matrix coefficient of dimension h by w
1039  explicit PWMatrixCoefficient(int h, int w, bool symm = false)
1040  : MatrixCoefficient(h, w, symm) {}
1041 
1042  /// Construct the coefficient using arrays describing the pieces
1043  /** \param dim - size of the square matrix-valued result
1044  \param attr - an array of attribute numbers for each piece
1045  \param coefs - the corresponding array of MatrixCoefficient pointers
1046  \param symm - true if the result will be symmetric, false otherwise
1047  Any missing attributes or NULL coefficient pointers will result in a
1048  zero matrix being returned.
1049 
1050  \note Ownership of the MatrixCoefficient objects will NOT be
1051  transferred to this object.
1052  */
1054  const Array<MatrixCoefficient*> & coefs,
1055  bool symm=false)
1056  : MatrixCoefficient(dim, symm) { InitMap(attr, coefs); }
1057 
1058  /// Construct the coefficient using arrays describing the pieces
1059  /** \param h - height of the matrix-valued result
1060  \param w - width of the matrix-valued result
1061  \param attr - an array of attribute numbers for each piece
1062  \param coefs - the corresponding array of MatrixCoefficient pointers
1063  \param symm - true if the result will be symmetric, false otherwise
1064  Any missing attributes or NULL coefficient pointers will result in a
1065  zero matrix being returned for that attribute.
1066 
1067  \note Ownership of the MatrixCoefficient objects will NOT be
1068  transferred to this object.
1069  */
1070  PWMatrixCoefficient(int h, int w, const Array<int> & attr,
1071  const Array<MatrixCoefficient*> & coefs,
1072  bool symm=false)
1073  : MatrixCoefficient(h, w, symm) { InitMap(attr, coefs); }
1074 
1075  /// Set the time for time dependent coefficients
1076  virtual void SetTime(double t);
1077 
1078  /// Replace a set of coefficients
1079  void UpdateCoefficients(const Array<int> & attr,
1080  const Array<MatrixCoefficient*> & coefs)
1081  { InitMap(attr, coefs); }
1082 
1083  /// Replace a single coefficient for a particular attribute
1084  void UpdateCoefficient(int attr, MatrixCoefficient & coef);
1085 
1086  /// Remove a single MatrixCoefficient for a particular attribute
1087  void ZeroCoefficient(int attr)
1088  { pieces.erase(attr); }
1089 
1090  /// Evaluate the coefficient.
1091  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1092  const IntegrationPoint &ip);
1093 };
1094 
1095 /** @brief A matrix coefficient with an optional scalar coefficient multiplier
1096  \a q. The matrix function can either be represented by a std function or
1097  a constant matrix provided when constructing this object. */
1099 {
1100 private:
1101  std::function<void(const Vector &, DenseMatrix &)> Function;
1102  std::function<void(const Vector &, Vector &)> SymmFunction; // deprecated
1103  std::function<void(const Vector &, double, DenseMatrix &)> TDFunction;
1104 
1105  Coefficient *Q;
1106  DenseMatrix mat;
1107 
1108 public:
1109  /// Define a time-independent square matrix coefficient from a std function
1110  /** \param dim - the size of the matrix
1111  \param F - time-independent function
1112  \param q - optional scalar Coefficient to scale the matrix coefficient */
1114  std::function<void(const Vector &, DenseMatrix &)> F,
1115  Coefficient *q = nullptr)
1116  : MatrixCoefficient(dim), Function(std::move(F)), Q(q), mat(0)
1117  { }
1118 
1119  /// Define a constant matrix coefficient times a scalar Coefficient
1120  /** \param m - constant matrix
1121  \param q - optional scalar Coefficient to scale the matrix coefficient */
1123  : MatrixCoefficient(m.Height(), m.Width()), Q(&q), mat(m)
1124  { }
1125 
1126  /** @brief Define a time-independent symmetric square matrix coefficient from
1127  a std function */
1128  /** \param dim - the size of the matrix
1129  \param SymmF - function used in EvalSymmetric
1130  \param q - optional scalar Coefficient to scale the matrix coefficient
1131  @deprecated Use another constructor without setting SymmFunction. */
1133  std::function<void(const Vector &, Vector &)> SymmF,
1134  Coefficient *q = NULL)
1135  : MatrixCoefficient(dim, true), SymmFunction(std::move(SymmF)), Q(q), mat(0)
1136  { }
1137 
1138  /// Define a time-dependent square matrix coefficient from a std function
1139  /** \param dim - the size of the matrix
1140  \param TDF - time-dependent function
1141  \param q - optional scalar Coefficient to scale the matrix coefficient */
1143  std::function<void(const Vector &, double, DenseMatrix &)> TDF,
1144  Coefficient *q = nullptr)
1145  : MatrixCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
1146  { }
1147 
1148  /// Set the time for internally stored coefficients
1149  void SetTime(double t);
1150 
1151  /// Evaluate the matrix coefficient at @a ip.
1152  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1153  const IntegrationPoint &ip);
1154 
1155  /// (DEPRECATED) Evaluate the symmetric matrix coefficient at @a ip.
1156  /** @deprecated Use Eval() instead. */
1157  virtual void EvalSymmetric(Vector &K, ElementTransformation &T,
1158  const IntegrationPoint &ip);
1159 
1161 };
1162 
1163 
1164 /** @brief Matrix coefficient defined by a matrix of scalar coefficients.
1165  Coefficients that are not set will evaluate to zero in the vector. The
1166  coefficient is stored as a flat Array with indexing (i,j) -> i*width+j. */
1168 {
1169 private:
1170  Array<Coefficient *> Coeff;
1171  Array<bool> ownCoeff;
1172 
1173 public:
1174  /** @brief Construct a coefficient matrix of dimensions @a dim * @a dim. The
1175  actual coefficients still need to be added with Set(). */
1176  explicit MatrixArrayCoefficient (int dim);
1177 
1178  /// Set the time for internally stored coefficients
1179  void SetTime(double t);
1180 
1181  /// Get the coefficient located at (i,j) in the matrix.
1182  Coefficient* GetCoeff (int i, int j) { return Coeff[i*width+j]; }
1183 
1184  /** @brief Set the coefficient located at (i,j) in the matrix. By default by
1185  default this will take ownership of the Coefficient passed in, but this
1186  can be overridden with the @a own parameter. */
1187  void Set(int i, int j, Coefficient * c, bool own=true);
1188 
1190 
1191  /// Evaluate coefficient located at (i,j) in the matrix using integration
1192  /// point @a ip.
1193  double Eval(int i, int j, ElementTransformation &T, const IntegrationPoint &ip)
1194  { return Coeff[i*width+j] ? Coeff[i*width+j] -> Eval(T, ip, GetTime()) : 0.0; }
1195 
1196  /// Evaluate the matrix coefficient @a ip.
1197  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1198  const IntegrationPoint &ip);
1199 
1200  virtual ~MatrixArrayCoefficient();
1201 };
1202 
1203 
1204 /** @brief Derived matrix coefficient that has the value of the parent matrix
1205  coefficient where it is active and is zero otherwise. */
1207 {
1208 private:
1209  MatrixCoefficient *c;
1210  Array<int> active_attr;
1211 
1212 public:
1213  /** @brief Construct with a parent matrix coefficient and an array of zeros
1214  and ones representing the attributes for which this coefficient should be
1215  active. */
1217  : MatrixCoefficient(mc.GetHeight(), mc.GetWidth())
1218  { c = &mc; attr.Copy(active_attr); }
1219 
1220  /// Set the time for internally stored coefficients
1221  void SetTime(double t);
1222 
1223  /// Evaluate the matrix coefficient at @a ip.
1224  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1225  const IntegrationPoint &ip);
1226 };
1227 
1228 /// Coefficients based on sums, products, or other functions of coefficients.
1229 ///@{
1230 /** @brief Scalar coefficient defined as the linear combination of two scalar
1231  coefficients or a scalar and a scalar coefficient */
1233 {
1234 private:
1235  double aConst;
1236  Coefficient * a;
1237  Coefficient * b;
1238 
1239  double alpha;
1240  double beta;
1241 
1242 public:
1243  /// Constructor with one coefficient. Result is alpha_ * A + beta_ * B
1245  double alpha_ = 1.0, double beta_ = 1.0)
1246  : aConst(A), a(NULL), b(&B), alpha(alpha_), beta(beta_) { }
1247 
1248  /// Constructor with two coefficients. Result is alpha_ * A + beta_ * B.
1250  double alpha_ = 1.0, double beta_ = 1.0)
1251  : aConst(0.0), a(&A), b(&B), alpha(alpha_), beta(beta_) { }
1252 
1253  /// Set the time for internally stored coefficients
1254  void SetTime(double t);
1255 
1256  /// Reset the first term in the linear combination as a constant
1257  void SetAConst(double A) { a = NULL; aConst = A; }
1258  /// Return the first term in the linear combination
1259  double GetAConst() const { return aConst; }
1260 
1261  /// Reset the first term in the linear combination
1262  void SetACoef(Coefficient &A) { a = &A; }
1263  /// Return the first term in the linear combination
1264  Coefficient * GetACoef() const { return a; }
1265 
1266  /// Reset the second term in the linear combination
1267  void SetBCoef(Coefficient &B) { b = &B; }
1268  /// Return the second term in the linear combination
1269  Coefficient * GetBCoef() const { return b; }
1270 
1271  /// Reset the factor in front of the first term in the linear combination
1272  void SetAlpha(double alpha_) { alpha = alpha_; }
1273  /// Return the factor in front of the first term in the linear combination
1274  double GetAlpha() const { return alpha; }
1275 
1276  /// Reset the factor in front of the second term in the linear combination
1277  void SetBeta(double beta_) { beta = beta_; }
1278  /// Return the factor in front of the second term in the linear combination
1279  double GetBeta() const { return beta; }
1280 
1281  /// Evaluate the coefficient at @a ip.
1282  virtual double Eval(ElementTransformation &T,
1283  const IntegrationPoint &ip)
1284  {
1285  return alpha * ((a == NULL ) ? aConst : a->Eval(T, ip) )
1286  + beta * b->Eval(T, ip);
1287  }
1288 };
1289 
1290 
1291 /// Base class for symmetric matrix coefficients that optionally depend on time and space.
1293 {
1294 protected:
1295  /// Internal matrix used when evaluating this coefficient as a DenseMatrix.
1297 public:
1298  /// Construct a dim x dim matrix coefficient.
1300  : MatrixCoefficient(dimension, true) { }
1301 
1302  /// Get the size of the matrix.
1303  int GetSize() const { return height; }
1304 
1305  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
1306  /// the quadrature points.
1307  ///
1308  /// @note As opposed to MatrixCoefficient::Project, this function stores only
1309  /// the @a symmetric part of the matrix at each quadrature point.
1310  ///
1311  /// The @a vdim of the coefficient should be equal to height*(height+1)/2.
1312  virtual void ProjectSymmetric(QuadratureFunction &qf);
1313 
1314  /** @brief Evaluate the matrix coefficient in the element described by @a T
1315  at the point @a ip, storing the result as a symmetric matrix @a K. */
1316  /** @note When this method is called, the caller must make sure that the
1317  IntegrationPoint associated with @a T is the same as @a ip. This can be
1318  achieved by calling T.SetIntPoint(&ip). */
1319  virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T,
1320  const IntegrationPoint &ip) = 0;
1321 
1322  /** @brief Evaluate the matrix coefficient in the element described by @a T
1323  at the point @a ip, storing the result as a dense matrix @a K. */
1324  /** This function allows the use of SymmetricMatrixCoefficient in situations
1325  where the symmetry is not taken advantage of.
1326 
1327  @note When this method is called, the caller must make sure that the
1328  IntegrationPoint associated with @a T is the same as @a ip. This can be
1329  achieved by calling T.SetIntPoint(&ip). */
1330  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1331  const IntegrationPoint &ip);
1332 
1333  /// Return a reference to the constant matrix.
1334  const DenseSymmetricMatrix& GetMatrix() { return mat; }
1335 
1337 };
1338 
1339 
1340 /// A matrix coefficient that is constant in space and time.
1342 {
1343 private:
1345 
1346 public:
1347  ///Construct using matrix @a m for the constant.
1349  : SymmetricMatrixCoefficient(m.Height()), mat(m) { }
1351  /// Evaluate the matrix coefficient at @a ip.
1353  const IntegrationPoint &ip) { M = mat; }
1354 };
1355 
1356 
1357 /** @brief A matrix coefficient with an optional scalar coefficient multiplier
1358  \a q. The matrix function can either be represented by a std function or
1359  a constant matrix provided when constructing this object. */
1361 {
1362 private:
1363  std::function<void(const Vector &, DenseSymmetricMatrix &)> Function;
1364  std::function<void(const Vector &, double, DenseSymmetricMatrix &)> TDFunction;
1365 
1366  Coefficient *Q;
1368 
1369 public:
1370  /// Define a time-independent symmetric matrix coefficient from a std function
1371  /** \param dim - the size of the matrix
1372  \param F - time-independent function
1373  \param q - optional scalar Coefficient to scale the matrix coefficient */
1375  std::function<void(const Vector &, DenseSymmetricMatrix &)> F,
1376  Coefficient *q = nullptr)
1377  : SymmetricMatrixCoefficient(dim), Function(std::move(F)), Q(q), mat(0)
1378  { }
1379 
1380  /// Define a constant matrix coefficient times a scalar Coefficient
1381  /** \param m - constant matrix
1382  \param q - optional scalar Coefficient to scale the matrix coefficient */
1384  Coefficient &q)
1385  : SymmetricMatrixCoefficient(m.Height()), Q(&q), mat(m)
1386  { }
1387 
1388  /// Define a time-dependent square matrix coefficient from a std function
1389  /** \param dim - the size of the matrix
1390  \param TDF - time-dependent function
1391  \param q - optional scalar Coefficient to scale the matrix coefficient */
1393  std::function<void(const Vector &, double, DenseSymmetricMatrix &)> TDF,
1394  Coefficient *q = nullptr)
1395  : SymmetricMatrixCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
1396  { }
1397 
1398  /// Set the time for internally stored coefficients
1399  void SetTime(double t);
1400 
1402  /// Evaluate the matrix coefficient at @a ip.
1403  virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T,
1404  const IntegrationPoint &ip);
1405 
1407 };
1408 
1409 
1410 /** @brief Scalar coefficient defined as the product of two scalar coefficients
1411  or a scalar and a scalar coefficient. */
1413 {
1414 private:
1415  double aConst;
1416  Coefficient * a;
1417  Coefficient * b;
1418 
1419 public:
1420  /// Constructor with one coefficient. Result is A * B.
1422  : aConst(A), a(NULL), b(&B) { }
1423 
1424  /// Constructor with two coefficients. Result is A * B.
1426  : aConst(0.0), a(&A), b(&B) { }
1427 
1428  /// Set the time for internally stored coefficients
1429  void SetTime(double t);
1430 
1431  /// Reset the first term in the product as a constant
1432  void SetAConst(double A) { a = NULL; aConst = A; }
1433  /// Return the first term in the product
1434  double GetAConst() const { return aConst; }
1435 
1436  /// Reset the first term in the product
1437  void SetACoef(Coefficient &A) { a = &A; }
1438  /// Return the first term in the product
1439  Coefficient * GetACoef() const { return a; }
1440 
1441  /// Reset the second term in the product
1442  void SetBCoef(Coefficient &B) { b = &B; }
1443  /// Return the second term in the product
1444  Coefficient * GetBCoef() const { return b; }
1445 
1446  /// Evaluate the coefficient at @a ip.
1447  virtual double Eval(ElementTransformation &T,
1448  const IntegrationPoint &ip)
1449  { return ((a == NULL ) ? aConst : a->Eval(T, ip) ) * b->Eval(T, ip); }
1450 };
1451 
1452 /** @brief Scalar coefficient defined as the ratio of two scalars where one or
1453  both scalars are scalar coefficients. */
1455 {
1456 private:
1457  double aConst;
1458  double bConst;
1459  Coefficient * a;
1460  Coefficient * b;
1461 
1462 public:
1463  /** Initialize a coefficient which returns A / B where @a A is a
1464  constant and @a B is a scalar coefficient */
1466  : aConst(A), bConst(1.0), a(NULL), b(&B) { }
1467  /** Initialize a coefficient which returns A / B where @a A and @a B are both
1468  scalar coefficients */
1470  : aConst(0.0), bConst(1.0), a(&A), b(&B) { }
1471  /** Initialize a coefficient which returns A / B where @a A is a
1472  scalar coefficient and @a B is a constant */
1474  : aConst(0.0), bConst(B), a(&A), b(NULL) { }
1475 
1476  /// Set the time for internally stored coefficients
1477  void SetTime(double t);
1478 
1479  /// Reset the numerator in the ratio as a constant
1480  void SetAConst(double A) { a = NULL; aConst = A; }
1481  /// Return the numerator of the ratio
1482  double GetAConst() const { return aConst; }
1483 
1484  /// Reset the denominator in the ratio as a constant
1485  void SetBConst(double B) { b = NULL; bConst = B; }
1486  /// Return the denominator of the ratio
1487  double GetBConst() const { return bConst; }
1488 
1489  /// Reset the numerator in the ratio
1490  void SetACoef(Coefficient &A) { a = &A; }
1491  /// Return the numerator of the ratio
1492  Coefficient * GetACoef() const { return a; }
1493 
1494  /// Reset the denominator in the ratio
1495  void SetBCoef(Coefficient &B) { b = &B; }
1496  /// Return the denominator of the ratio
1497  Coefficient * GetBCoef() const { return b; }
1498 
1499  /// Evaluate the coefficient
1500  virtual double Eval(ElementTransformation &T,
1501  const IntegrationPoint &ip)
1502  {
1503  double den = (b == NULL ) ? bConst : b->Eval(T, ip);
1504  MFEM_ASSERT(den != 0.0, "Division by zero in RatioCoefficient");
1505  return ((a == NULL ) ? aConst : a->Eval(T, ip) ) / den;
1506  }
1507 };
1508 
1509 /// Scalar coefficient defined as a scalar raised to a power
1511 {
1512 private:
1513  Coefficient * a;
1514 
1515  double p;
1516 
1517 public:
1518  /// Construct with a coefficient and a constant power @a p_. Result is A^p.
1520  : a(&A), p(p_) { }
1521 
1522  /// Set the time for internally stored coefficients
1523  void SetTime(double t);
1524 
1525  /// Reset the base coefficient
1526  void SetACoef(Coefficient &A) { a = &A; }
1527  /// Return the base coefficient
1528  Coefficient * GetACoef() const { return a; }
1529 
1530  /// Reset the exponent
1531  void SetExponent(double p_) { p = p_; }
1532  /// Return the exponent
1533  double GetExponent() const { return p; }
1534 
1535  /// Evaluate the coefficient at @a ip.
1536  virtual double Eval(ElementTransformation &T,
1537  const IntegrationPoint &ip)
1538  { return pow(a->Eval(T, ip), p); }
1539 };
1540 
1541 
1542 /// Scalar coefficient defined as the inner product of two vector coefficients
1544 {
1545 private:
1546  VectorCoefficient * a;
1547  VectorCoefficient * b;
1548 
1549  mutable Vector va;
1550  mutable Vector vb;
1551 public:
1552  /// Construct with the two vector coefficients. Result is \f$ A \cdot B \f$.
1554 
1555  /// Set the time for internally stored coefficients
1556  void SetTime(double t);
1557 
1558  /// Reset the first vector in the inner product
1559  void SetACoef(VectorCoefficient &A) { a = &A; }
1560  /// Return the first vector coefficient in the inner product
1561  VectorCoefficient * GetACoef() const { return a; }
1562 
1563  /// Reset the second vector in the inner product
1564  void SetBCoef(VectorCoefficient &B) { b = &B; }
1565  /// Return the second vector coefficient in the inner product
1566  VectorCoefficient * GetBCoef() const { return b; }
1567 
1568  /// Evaluate the coefficient at @a ip.
1569  virtual double Eval(ElementTransformation &T,
1570  const IntegrationPoint &ip);
1571 };
1572 
1573 /// Scalar coefficient defined as a cross product of two vectors in the xy-plane.
1575 {
1576 private:
1577  VectorCoefficient * a;
1578  VectorCoefficient * b;
1579 
1580  mutable Vector va;
1581  mutable Vector vb;
1582 
1583 public:
1584  /// Constructor with two vector coefficients. Result is \f$ A_x B_y - A_y * B_x; \f$.
1586 
1587  /// Set the time for internally stored coefficients
1588  void SetTime(double t);
1589 
1590  /// Reset the first vector in the product
1591  void SetACoef(VectorCoefficient &A) { a = &A; }
1592  /// Return the first vector of the product
1593  VectorCoefficient * GetACoef() const { return a; }
1594 
1595  /// Reset the second vector in the product
1596  void SetBCoef(VectorCoefficient &B) { b = &B; }
1597  /// Return the second vector of the product
1598  VectorCoefficient * GetBCoef() const { return b; }
1599 
1600  /// Evaluate the coefficient at @a ip.
1601  virtual double Eval(ElementTransformation &T,
1602  const IntegrationPoint &ip);
1603 };
1604 
1605 /// Scalar coefficient defined as the determinant of a matrix coefficient
1607 {
1608 private:
1609  MatrixCoefficient * a;
1610 
1611  mutable DenseMatrix ma;
1612 
1613 public:
1614  /// Construct with the matrix.
1616 
1617  /// Set the time for internally stored coefficients
1618  void SetTime(double t);
1619 
1620  /// Reset the matrix coefficient
1621  void SetACoef(MatrixCoefficient &A) { a = &A; }
1622  /// Return the matrix coefficient
1623  MatrixCoefficient * GetACoef() const { return a; }
1624 
1625  /// Evaluate the determinant coefficient at @a ip.
1626  virtual double Eval(ElementTransformation &T,
1627  const IntegrationPoint &ip);
1628 };
1629 
1630 /// Vector coefficient defined as the linear combination of two vectors
1632 {
1633 private:
1634  VectorCoefficient * ACoef;
1635  VectorCoefficient * BCoef;
1636 
1637  Vector A;
1638  Vector B;
1639 
1640  Coefficient * alphaCoef;
1641  Coefficient * betaCoef;
1642 
1643  double alpha;
1644  double beta;
1645 
1646  mutable Vector va;
1647 
1648 public:
1649  /** Constructor with no coefficients.
1650  To be used with the various "Set" methods */
1652 
1653  /** Constructor with two vector coefficients.
1654  Result is alpha_ * A + beta_ * B */
1656  double alpha_ = 1.0, double beta_ = 1.0);
1657 
1658  /** Constructor with scalar coefficients.
1659  Result is alpha_ * A_ + beta_ * B_ */
1661  Coefficient &alpha_, Coefficient &beta_);
1662 
1663  /// Set the time for internally stored coefficients
1664  void SetTime(double t);
1665 
1666  /// Reset the first vector coefficient
1667  void SetACoef(VectorCoefficient &A_) { ACoef = &A_; }
1668  /// Return the first vector coefficient
1669  VectorCoefficient * GetACoef() const { return ACoef; }
1670 
1671  /// Reset the second vector coefficient
1672  void SetBCoef(VectorCoefficient &B_) { BCoef = &B_; }
1673  /// Return the second vector coefficient
1674  VectorCoefficient * GetBCoef() const { return BCoef; }
1675 
1676  /// Reset the factor in front of the first vector coefficient
1677  void SetAlphaCoef(Coefficient &A_) { alphaCoef = &A_; }
1678  /// Return the factor in front of the first vector coefficient
1679  Coefficient * GetAlphaCoef() const { return alphaCoef; }
1680 
1681  /// Reset the factor in front of the second vector coefficient
1682  void SetBetaCoef(Coefficient &B_) { betaCoef = &B_; }
1683  /// Return the factor in front of the second vector coefficient
1684  Coefficient * GetBetaCoef() const { return betaCoef; }
1685 
1686  /// Reset the first vector as a constant
1687  void SetA(const Vector &A_) { A = A_; ACoef = NULL; }
1688  /// Return the first vector constant
1689  const Vector & GetA() const { return A; }
1690 
1691  /// Reset the second vector as a constant
1692  void SetB(const Vector &B_) { B = B_; BCoef = NULL; }
1693  /// Return the second vector constant
1694  const Vector & GetB() const { return B; }
1695 
1696  /// Reset the factor in front of the first vector coefficient as a constant
1697  void SetAlpha(double alpha_) { alpha = alpha_; alphaCoef = NULL; }
1698  /// Return the factor in front of the first vector coefficient
1699  double GetAlpha() const { return alpha; }
1700 
1701  /// Reset the factor in front of the second vector coefficient as a constant
1702  void SetBeta(double beta_) { beta = beta_; betaCoef = NULL; }
1703  /// Return the factor in front of the second vector coefficient
1704  double GetBeta() const { return beta; }
1705 
1706  /// Evaluate the coefficient at @a ip.
1707  virtual void Eval(Vector &V, ElementTransformation &T,
1708  const IntegrationPoint &ip);
1710 };
1711 
1712 /// Vector coefficient defined as a product of scalar and vector coefficients.
1714 {
1715 private:
1716  double aConst;
1717  Coefficient * a;
1718  VectorCoefficient * b;
1719 
1720 public:
1721  /// Constructor with constant and vector coefficient. Result is A * B.
1723 
1724  /// Constructor with two coefficients. Result is A * B.
1726 
1727  /// Set the time for internally stored coefficients
1728  void SetTime(double t);
1729 
1730  /// Reset the scalar factor as a constant
1731  void SetAConst(double A) { a = NULL; aConst = A; }
1732  /// Return the scalar factor
1733  double GetAConst() const { return aConst; }
1734 
1735  /// Reset the scalar factor
1736  void SetACoef(Coefficient &A) { a = &A; }
1737  /// Return the scalar factor
1738  Coefficient * GetACoef() const { return a; }
1739 
1740  /// Reset the vector factor
1741  void SetBCoef(VectorCoefficient &B) { b = &B; }
1742  /// Return the vector factor
1743  VectorCoefficient * GetBCoef() const { return b; }
1744 
1745  /// Evaluate the coefficient at @a ip.
1746  virtual void Eval(Vector &V, ElementTransformation &T,
1747  const IntegrationPoint &ip);
1749 };
1750 
1751 /// Vector coefficient defined as a normalized vector field (returns v/|v|)
1753 {
1754 private:
1755  VectorCoefficient * a;
1756 
1757  double tol;
1758 
1759 public:
1760  /** @brief Return a vector normalized to a length of one
1761 
1762  This class evaluates the vector coefficient @a A and, if |A| > @a tol,
1763  returns the normalized vector A / |A|. If |A| <= @a tol, the zero
1764  vector is returned.
1765  */
1766  NormalizedVectorCoefficient(VectorCoefficient &A, double tol = 1e-6);
1767 
1768  /// Set the time for internally stored coefficients
1769  void SetTime(double t);
1770 
1771  /// Reset the vector coefficient
1772  void SetACoef(VectorCoefficient &A) { a = &A; }
1773  /// Return the vector coefficient
1774  VectorCoefficient * GetACoef() const { return a; }
1775 
1776  /// Evaluate the coefficient at @a ip.
1777  virtual void Eval(Vector &V, ElementTransformation &T,
1778  const IntegrationPoint &ip);
1780 };
1781 
1782 /// Vector coefficient defined as a cross product of two vectors
1784 {
1785 private:
1786  VectorCoefficient * a;
1787  VectorCoefficient * b;
1788 
1789  mutable Vector va;
1790  mutable Vector vb;
1791 
1792 public:
1793  /// Construct with the two coefficients. Result is A x B.
1795 
1796  /// Set the time for internally stored coefficients
1797  void SetTime(double t);
1798 
1799  /// Reset the first term in the product
1800  void SetACoef(VectorCoefficient &A) { a = &A; }
1801  /// Return the first term in the product
1802  VectorCoefficient * GetACoef() const { return a; }
1803 
1804  /// Reset the second term in the product
1805  void SetBCoef(VectorCoefficient &B) { b = &B; }
1806  /// Return the second term in the product
1807  VectorCoefficient * GetBCoef() const { return b; }
1808 
1809  /// Evaluate the coefficient at @a ip.
1810  virtual void Eval(Vector &V, ElementTransformation &T,
1811  const IntegrationPoint &ip);
1813 };
1814 
1815 /** @brief Vector coefficient defined as a product of a matrix coefficient and
1816  a vector coefficient. */
1818 {
1819 private:
1820  MatrixCoefficient * a;
1821  VectorCoefficient * b;
1822 
1823  mutable DenseMatrix ma;
1824  mutable Vector vb;
1825 
1826 public:
1827  /// Constructor with two coefficients. Result is A*B.
1829 
1830  /// Set the time for internally stored coefficients
1831  void SetTime(double t);
1832 
1833  /// Reset the matrix coefficient
1834  void SetACoef(MatrixCoefficient &A) { a = &A; }
1835  /// Return the matrix coefficient
1836  MatrixCoefficient * GetACoef() const { return a; }
1837 
1838  /// Reset the vector coefficient
1839  void SetBCoef(VectorCoefficient &B) { b = &B; }
1840  /// Return the vector coefficient
1841  VectorCoefficient * GetBCoef() const { return b; }
1842 
1843  /// Evaluate the vector coefficient at @a ip.
1844  virtual void Eval(Vector &V, ElementTransformation &T,
1845  const IntegrationPoint &ip);
1847 };
1848 
1849 /// Convenient alias for the MatrixVectorProductCoefficient
1851 
1852 /// Constant matrix coefficient defined as the identity of dimension d
1854 {
1855 private:
1856  int dim;
1857 
1858 public:
1859  /// Construct with the dimension of the square identity matrix.
1861  : MatrixCoefficient(d, d), dim(d) { }
1862 
1863  /// Evaluate the matrix coefficient at @a ip.
1864  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1865  const IntegrationPoint &ip);
1866 };
1867 
1868 /// Matrix coefficient defined as the linear combination of two matrices
1870 {
1871 private:
1872  MatrixCoefficient * a;
1873  MatrixCoefficient * b;
1874 
1875  double alpha;
1876  double beta;
1877 
1878  mutable DenseMatrix ma;
1879 
1880 public:
1881  /// Construct with the two coefficients. Result is alpha_ * A + beta_ * B.
1883  double alpha_ = 1.0, double beta_ = 1.0);
1884 
1885  /// Set the time for internally stored coefficients
1886  void SetTime(double t);
1887 
1888  /// Reset the first matrix coefficient
1889  void SetACoef(MatrixCoefficient &A) { a = &A; }
1890  /// Return the first matrix coefficient
1891  MatrixCoefficient * GetACoef() const { return a; }
1892 
1893  /// Reset the second matrix coefficient
1894  void SetBCoef(MatrixCoefficient &B) { b = &B; }
1895  /// Return the second matrix coefficient
1896  MatrixCoefficient * GetBCoef() const { return b; }
1897 
1898  /// Reset the factor in front of the first matrix coefficient
1899  void SetAlpha(double alpha_) { alpha = alpha_; }
1900  /// Return the factor in front of the first matrix coefficient
1901  double GetAlpha() const { return alpha; }
1902 
1903  /// Reset the factor in front of the second matrix coefficient
1904  void SetBeta(double beta_) { beta = beta_; }
1905  /// Return the factor in front of the second matrix coefficient
1906  double GetBeta() const { return beta; }
1907 
1908  /// Evaluate the matrix coefficient at @a ip.
1909  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1910  const IntegrationPoint &ip);
1911 };
1912 
1913 /// Matrix coefficient defined as the product of two matrices
1915 {
1916 private:
1917  MatrixCoefficient * a;
1918  MatrixCoefficient * b;
1919 
1920  mutable DenseMatrix ma;
1921  mutable DenseMatrix mb;
1922 
1923 public:
1924  /// Construct with the two coefficients. Result is A * B.
1926 
1927  /// Reset the first matrix coefficient
1928  void SetACoef(MatrixCoefficient &A) { a = &A; }
1929  /// Return the first matrix coefficient
1930  MatrixCoefficient * GetACoef() const { return a; }
1931 
1932  /// Reset the second matrix coefficient
1933  void SetBCoef(MatrixCoefficient &B) { b = &B; }
1934  /// Return the second matrix coefficient
1935  MatrixCoefficient * GetBCoef() const { return b; }
1936 
1937  /// Evaluate the matrix coefficient at @a ip.
1938  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1939  const IntegrationPoint &ip);
1940 };
1941 
1942 /** @brief Matrix coefficient defined as a product of a scalar coefficient and a
1943  matrix coefficient.*/
1945 {
1946 private:
1947  double aConst;
1948  Coefficient * a;
1949  MatrixCoefficient * b;
1950 
1951 public:
1952  /// Constructor with one coefficient. Result is A*B.
1954 
1955  /// Constructor with two coefficients. Result is A*B.
1957 
1958  /// Set the time for internally stored coefficients
1959  void SetTime(double t);
1960 
1961  /// Reset the scalar factor as a constant
1962  void SetAConst(double A) { a = NULL; aConst = A; }
1963  /// Return the scalar factor
1964  double GetAConst() const { return aConst; }
1965 
1966  /// Reset the scalar factor
1967  void SetACoef(Coefficient &A) { a = &A; }
1968  /// Return the scalar factor
1969  Coefficient * GetACoef() const { return a; }
1970 
1971  /// Reset the matrix factor
1972  void SetBCoef(MatrixCoefficient &B) { b = &B; }
1973  /// Return the matrix factor
1974  MatrixCoefficient * GetBCoef() const { return b; }
1975 
1976  /// Evaluate the matrix coefficient at @a ip.
1977  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1978  const IntegrationPoint &ip);
1979 };
1980 
1981 /// Matrix coefficient defined as the transpose a matrix coefficient
1983 {
1984 private:
1985  MatrixCoefficient * a;
1986 
1987 public:
1988  /// Construct with the matrix coefficient. Result is \f$ A^T \f$.
1990 
1991  /// Set the time for internally stored coefficients
1992  void SetTime(double t);
1993 
1994  /// Reset the matrix coefficient
1995  void SetACoef(MatrixCoefficient &A) { a = &A; }
1996  /// Return the matrix coefficient
1997  MatrixCoefficient * GetACoef() const { return a; }
1998 
1999  /// Evaluate the matrix coefficient at @a ip.
2000  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2001  const IntegrationPoint &ip);
2002 };
2003 
2004 /// Matrix coefficient defined as the inverse a matrix coefficient.
2006 {
2007 private:
2008  MatrixCoefficient * a;
2009 
2010 public:
2011  /// Construct with the matrix coefficient. Result is \f$ A^{-1} \f$.
2013 
2014  /// Set the time for internally stored coefficients
2015  void SetTime(double t);
2016 
2017  /// Reset the matrix coefficient
2018  void SetACoef(MatrixCoefficient &A) { a = &A; }
2019  /// Return the matrix coefficient
2020  MatrixCoefficient * GetACoef() const { return a; }
2021 
2022  /// Evaluate the matrix coefficient at @a ip.
2023  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2024  const IntegrationPoint &ip);
2025 };
2026 
2027 /// Matrix coefficient defined as the outer product of two vector coefficients.
2029 {
2030 private:
2031  VectorCoefficient * a;
2032  VectorCoefficient * b;
2033 
2034  mutable Vector va;
2035  mutable Vector vb;
2036 
2037 public:
2038  /// Construct with two vector coefficients. Result is \f$ A B^T \f$.
2040 
2041  /// Set the time for internally stored coefficients
2042  void SetTime(double t);
2043 
2044  /// Reset the first vector in the outer product
2045  void SetACoef(VectorCoefficient &A) { a = &A; }
2046  /// Return the first vector coefficient in the outer product
2047  VectorCoefficient * GetACoef() const { return a; }
2048 
2049  /// Reset the second vector in the outer product
2050  void SetBCoef(VectorCoefficient &B) { b = &B; }
2051  /// Return the second vector coefficient in the outer product
2052  VectorCoefficient * GetBCoef() const { return b; }
2053 
2054  /// Evaluate the matrix coefficient at @a ip.
2055  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2056  const IntegrationPoint &ip);
2057 };
2058 
2059 /** @brief Matrix coefficient defined as -a k x k x, for a vector k and scalar a
2060 
2061  This coefficient returns \f$a * (|k|^2 I - k \otimes k)\f$, where I is
2062  the identity matrix and \f$\otimes\f$ indicates the outer product. This
2063  can be evaluated for vectors of any dimension but in three
2064  dimensions it corresponds to computing the cross product with k twice.
2065 */
2067 {
2068 private:
2069  double aConst;
2070  Coefficient * a;
2071  VectorCoefficient * k;
2072 
2073  mutable Vector vk;
2074 
2075 public:
2078 
2079  /// Set the time for internally stored coefficients
2080  void SetTime(double t);
2081 
2082  /// Reset the scalar factor as a constant
2083  void SetAConst(double A) { a = NULL; aConst = A; }
2084  /// Return the scalar factor
2085  double GetAConst() const { return aConst; }
2086 
2087  /// Reset the scalar factor
2088  void SetACoef(Coefficient &A) { a = &A; }
2089  /// Return the scalar factor
2090  Coefficient * GetACoef() const { return a; }
2091 
2092  /// Reset the vector factor
2093  void SetKCoef(VectorCoefficient &K) { k = &K; }
2094  /// Return the vector factor
2095  VectorCoefficient * GetKCoef() const { return k; }
2096 
2097  /// Evaluate the matrix coefficient at @a ip.
2098  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2099  const IntegrationPoint &ip);
2100 };
2101 ///@}
2102 
2103 /** @brief Vector quadrature function coefficient which requires that the
2104  quadrature rules used for this vector coefficient be the same as those that
2105  live within the supplied QuadratureFunction. */
2107 {
2108 private:
2109  const QuadratureFunction &QuadF; //do not own
2110  int index;
2111 
2112 public:
2113  /// Constructor with a quadrature function as input
2115 
2116  /** Set the starting index within the QuadFunc that'll be used to project
2117  outwards as well as the corresponding length. The projected length should
2118  have the bounds of 1 <= length <= (length QuadFunc - index). */
2119  void SetComponent(int index_, int length_);
2120 
2121  const QuadratureFunction& GetQuadFunction() const { return QuadF; }
2122 
2124  virtual void Eval(Vector &V, ElementTransformation &T,
2125  const IntegrationPoint &ip);
2126 
2127  virtual void Project(QuadratureFunction &qf);
2128 
2130 };
2131 
2132 /** @brief Quadrature function coefficient which requires that the quadrature
2133  rules used for this coefficient be the same as those that live within the
2134  supplied QuadratureFunction. */
2136 {
2137 private:
2138  const QuadratureFunction &QuadF;
2139 
2140 public:
2141  /// Constructor with a quadrature function as input
2143 
2144  const QuadratureFunction& GetQuadFunction() const { return QuadF; }
2145 
2146  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
2147 
2148  virtual void Project(QuadratureFunction &qf);
2149 
2151 };
2152 
2153 /// Flags that determine what storage optimizations to use in CoefficientVector
2154 enum class CoefficientStorage : int
2155 {
2156  FULL = 0, ///< Store the coefficient as a full QuadratureFunction.
2157  CONSTANTS = 1 << 0, ///< Store constants using only @a vdim entries.
2158  SYMMETRIC = 1 << 1, ///< Store the triangular part of symmetric matrices.
2159  COMPRESSED = CONSTANTS | SYMMETRIC ///< Enable all above compressions.
2160 };
2161 
2163 {
2164  return CoefficientStorage(int(a) | int(b));
2165 }
2166 
2168 {
2169  return int(a) & int(b);
2170 }
2171 
2172 
2173 /// @brief Class to represent a coefficient evaluated at quadrature points.
2174 ///
2175 /// In the general case, a CoefficientVector is the same as a QuadratureFunction
2176 /// with a coefficient projected onto it.
2177 ///
2178 /// This class allows for some "compression" of the coefficient data, according
2179 /// to the storage flags given by CoefficientStorage. For example, constant
2180 /// coefficients can be stored using only @a vdim values, and symmetric matrices
2181 /// can be stored using e.g. the upper triangular part of the matrix.
2183 {
2184 protected:
2185  CoefficientStorage storage; ///< Storage optimizations (see CoefficientStorage).
2186  int vdim; ///< Number of values per quadrature point.
2187  QuadratureSpaceBase &qs; ///< Associated QuadratureSpaceBase.
2188  QuadratureFunction *qf; ///< Internal QuadratureFunction (owned, may be NULL).
2189 public:
2190  /// Create an empty CoefficientVector.
2193 
2194  /// @brief Create a CoefficientVector from the given Coefficient and
2195  /// QuadratureSpaceBase.
2196  ///
2197  /// If @a coeff is NULL, it will be interpreted as a constant with value one.
2198  /// @sa CoefficientStorage for a description of @a storage_.
2201 
2202  /// @brief Create a CoefficientVector from the given Coefficient and
2203  /// QuadratureSpaceBase.
2204  ///
2205  /// @sa CoefficientStorage for a description of @a storage_.
2208 
2209  /// @brief Create a CoefficientVector from the given VectorCoefficient and
2210  /// QuadratureSpaceBase.
2211  ///
2212  /// @sa CoefficientStorage for a description of @a storage_.
2215 
2216  /// @brief Create a CoefficientVector from the given MatrixCoefficient and
2217  /// QuadratureSpaceBase.
2218  ///
2219  /// @sa CoefficientStorage for a description of @a storage_.
2222 
2223  /// @brief Evaluate the given Coefficient at the quadrature points defined by
2224  /// @ref qs.
2225  void Project(Coefficient &coeff);
2226 
2227  /// @brief Evaluate the given VectorCoefficient at the quadrature points
2228  /// defined by @ref qs.
2229  ///
2230  /// @sa CoefficientVector for a description of the @a compress argument.
2231  void Project(VectorCoefficient &coeff);
2232 
2233  /// @brief Evaluate the given MatrixCoefficient at the quadrature points
2234  /// defined by @ref qs.
2235  ///
2236  /// @sa CoefficientVector for a description of the @a compress argument.
2237  void Project(MatrixCoefficient &coeff, bool transpose=false);
2238 
2239  /// @brief Project the transpose of @a coeff.
2240  ///
2241  /// @sa Project(MatrixCoefficient&, QuadratureSpace&, bool, bool)
2242  void ProjectTranspose(MatrixCoefficient &coeff);
2243 
2244  /// Make this vector a reference to the given QuadratureFunction.
2245  void MakeRef(const QuadratureFunction &qf_);
2246 
2247  /// Set this vector to the given constant.
2248  void SetConstant(double constant);
2249 
2250  /// Set this vector to the given constant vector.
2251  void SetConstant(const Vector &constant);
2252 
2253  /// Set this vector to the given constant matrix.
2254  void SetConstant(const DenseMatrix &constant);
2255 
2256  /// Set this vector to the given constant symmetric matrix.
2257  void SetConstant(const DenseSymmetricMatrix &constant);
2258 
2259  /// Return the number of values per quadrature point.
2260  int GetVDim() const;
2261 
2263 };
2264 
2265 /** @brief Compute the Lp norm of a function f.
2266  \f$ \| f \|_{Lp} = ( \int_\Omega | f |^p d\Omega)^{1/p} \f$ */
2267 double ComputeLpNorm(double p, Coefficient &coeff, Mesh &mesh,
2268  const IntegrationRule *irs[]);
2269 
2270 /** @brief Compute the Lp norm of a vector function f = {f_i}_i=1...N.
2271  \f$ \| f \|_{Lp} = ( \sum_i \| f_i \|_{Lp}^p )^{1/p} \f$ */
2272 double ComputeLpNorm(double p, VectorCoefficient &coeff, Mesh &mesh,
2273  const IntegrationRule *irs[]);
2274 
2275 #ifdef MFEM_USE_MPI
2276 /** @brief Compute the global Lp norm of a function f.
2277  \f$ \| f \|_{Lp} = ( \int_\Omega | f |^p d\Omega)^{1/p} \f$ */
2278 double ComputeGlobalLpNorm(double p, Coefficient &coeff, ParMesh &pmesh,
2279  const IntegrationRule *irs[]);
2280 
2281 /** @brief Compute the global Lp norm of a vector function f = {f_i}_i=1...N.
2282  \f$ \| f \|_{Lp} = ( \sum_i \| f_i \|_{Lp}^p )^{1/p} \f$ */
2283 double ComputeGlobalLpNorm(double p, VectorCoefficient &coeff, ParMesh &pmesh,
2284  const IntegrationRule *irs[]);
2285 #endif
2286 
2287 }
2288 
2289 #endif
FunctionCoefficient(std::function< double(const Vector &)> F)
Define a time-independent coefficient from a std function.
Vector coefficient defined as a cross product of two vectors.
int GetHeight() const
Get the height of the matrix.
A coefficient that depends on 1 or 2 parent coefficients and a transformation rule represented by a C...
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient at ip.
CoefficientStorage operator|(CoefficientStorage a, CoefficientStorage b)
RatioCoefficient(Coefficient &A, Coefficient &B)
A matrix coefficient that is constant in space and time.
Matrix coefficient defined as the product of two matrices.
void SetACoef(Coefficient &A)
Reset the scalar factor.
double Eval(ElementTransformation &T, const IntegrationPoint &ip, double t)
Evaluate the coefficient in the element described by T at the point ip at time t. ...
Definition: coefficient.hpp:68
CoefficientStorage
Flags that determine what storage optimizations to use in CoefficientVector.
void GetDeltaCenter(Vector &center)
MatrixCoefficient * GetACoef() const
Return the matrix coefficient.
void SetACoef(MatrixCoefficient &A)
Reset the matrix coefficient.
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
Vector coefficient defined by an array of scalar coefficients. Coefficients that are not set will eva...
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient at ip.
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
MatrixFunctionCoefficient(int dim, std::function< void(const Vector &, double, DenseMatrix &)> TDF, Coefficient *q=nullptr)
Define a time-dependent square matrix coefficient from a std function.
A piecewise coefficient with the pieces keyed off the element attribute numbers.
SymmetricMatrixFunctionCoefficient(int dim, std::function< void(const Vector &, DenseSymmetricMatrix &)> F, Coefficient *q=nullptr)
Define a time-independent symmetric matrix coefficient from a std function.
void Set(int i, int j, Coefficient *c, bool own=true)
Set the coefficient located at (i,j) in the matrix. By default by default this will take ownership of...
IdentityMatrixCoefficient(int d)
Construct with the dimension of the square identity matrix.
PWConstCoefficient(Vector &c)
Construct the constant coefficient using a vector of constants.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
const GridFunction * GetGridFunction() const
Get the scalar grid function.
void SetBCoef(VectorCoefficient &B_)
Reset the second vector coefficient.
const DenseMatrix & GetMatrix()
Return a reference to the constant matrix.
void SetWeight(Coefficient *w)
Set a weight Coefficient that multiplies the DeltaCoefficient.
Base class for vector Coefficients that optionally depend on time and space.
A coefficient that is constant across space and time.
Definition: coefficient.hpp:84
Coefficient ** GetCoeffs()
Returns the entire array of coefficients.
double GetAConst() const
Return the first term in the product.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector curl coefficient at ip.
void SetAConst(double A)
Reset the first term in the product as a constant.
const QuadratureFunction & GetQuadFunction() const
TransposeMatrixCoefficient(MatrixCoefficient &A)
Construct with the matrix coefficient. Result is .
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the vector coefficient in the element described by T at the point ip, storing the result in ...
VectorCrossProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
Construct with the two coefficients. Result is A x B.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
Scalar coefficient defined as the inner product of two vector coefficients.
NormalizedVectorCoefficient(VectorCoefficient &A, double tol=1e-6)
Return a vector normalized to a length of one.
int vdim
Number of values per quadrature point.
void SetBCoef(MatrixCoefficient &B)
Reset the matrix factor.
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:512
void SetKCoef(VectorCoefficient &K)
Reset the vector factor.
MatrixVectorProductCoefficient(MatrixCoefficient &A, VectorCoefficient &B)
Constructor with two coefficients. Result is A*B.
double Eval(int i, ElementTransformation &T, const IntegrationPoint &ip)
VectorDeltaCoefficient(int vdim_)
Construct with a vector of dimension vdim_.
void SetDeltaCenter(const Vector &center)
Set the center location of the delta function.
MatrixCoefficient(int h, int w, bool symm=false)
Construct a h x w matrix coefficient.
double Scale()
Return the scale factor times the optional time dependent function. Returns with when not set by th...
void SetTime(double t)
Set the time for internally stored coefficients.
int GetSize() const
Get the size of the matrix.
Vector quadrature function coefficient which requires that the quadrature rules used for this vector ...
void SetACoef(Coefficient &A)
Reset the numerator in the ratio.
Vector coefficient that is constant in space and time.
PWCoefficient(const Array< int > &attr, const Array< Coefficient *> &coefs)
Construct the coefficient using arrays describing the pieces.
PWMatrixCoefficient(int h, int w, const Array< int > &attr, const Array< MatrixCoefficient *> &coefs, bool symm=false)
Construct the coefficient using arrays describing the pieces.
void SetTime(double t)
Set the time for internally stored coefficients.
const Vector & GetVec() const
Return a reference to the constant vector in this class.
void SetTime(double t)
Set the time for internally stored coefficients.
Scalar coefficient defined as a cross product of two vectors in the xy-plane.
double GetTime()
Get the time for time dependent coefficients.
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
Scalar coefficient defined as the product of two scalar coefficients or a scalar and a scalar coeffic...
void SetTime(double t)
Set the time for internally stored coefficients.
void SetBCoef(VectorCoefficient &B)
Reset the second term in the product.
void SetTime(double t)
Set the time for internally stored coefficients.
MatrixCoefficient * GetACoef() const
Return the first matrix coefficient.
void SetBeta(double beta_)
Reset the factor in front of the second matrix coefficient.
int Size() const
Returns the size of the vector.
Definition: vector.hpp:199
MatrixCoefficient * GetBCoef() const
Return the second matrix coefficient.
Coefficient defined by a GridFunction. This coefficient is mesh dependent.
double GetAConst() const
Return the first term in the linear combination.
double & operator()(int i)
Return a reference to the i-th constant.
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
Delta function coefficient optionally multiplied by a weight coefficient and a scaled time dependent ...
void SetExponent(double p_)
Reset the exponent.
void SetAlphaCoef(Coefficient &A_)
Reset the factor in front of the first vector coefficient.
A matrix coefficient with an optional scalar coefficient multiplier q. The matrix function can either...
void SetFunction(double(*f)(double))
Set a time-dependent function that multiplies the Scale().
Coefficient * GetACoef() const
Return the first term in the linear combination.
Scalar coefficient defined as the ratio of two scalars where one or both scalars are scalar coefficie...
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient at ip.
CurlGridFunctionCoefficient(const GridFunction *gf)
Construct the coefficient with a vector grid function gf. The grid function is not owned by the coeff...
Matrix coefficient defined by a matrix of scalar coefficients. Coefficients that are not set will eva...
void UpdateCoefficient(int attr, VectorCoefficient &coef)
Replace a single Coefficient for a particular attribute.
GradientGridFunctionCoefficient(const GridFunction *gf)
Construct the coefficient with a scalar grid function gf. The grid function is not owned by the coeff...
CrossCrossCoefficient(double A, VectorCoefficient &K)
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void SetB(const Vector &B_)
Reset the second vector as a constant.
void SetDeltaCoefficient(const DeltaCoefficient &d_)
Replace the associated DeltaCoefficient with a new DeltaCoefficient.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
double Tol()
Return the tolerance used to identify the mesh vertices.
void SetBCoef(MatrixCoefficient &B)
Reset the second matrix coefficient.
SymmetricMatrixFunctionCoefficient(const DenseSymmetricMatrix &m, Coefficient &q)
Define a constant matrix coefficient times a scalar Coefficient.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void SetGridFunction(const GridFunction *gf)
Set the vector grid function.
double GetAConst() const
Return the scalar factor.
Vector coefficient defined as the linear combination of two vectors.
MFEM_DEPRECATED FunctionCoefficient(double(*f)(Vector &))
(DEPRECATED) Define a time-independent coefficient from a C-function
Matrix coefficient defined as a product of a scalar coefficient and a matrix coefficient.
STL namespace.
void SetConstant(double constant)
Set this vector to the given constant.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
int operator &(CoefficientStorage a, CoefficientStorage b)
virtual void Eval(DenseSymmetricMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
VectorCoefficient * GetBCoef() const
Return the second vector coefficient.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
A DeltaFunction cannot be evaluated. Calling this method will cause an MFEM error, terminating the application.
int GetVDim() const
For backward compatibility get the width of the matrix.
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
Coefficient * GetBCoef() const
Return the denominator of the ratio.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
Definition: coefficient.hpp:50
void SetBCoef(VectorCoefficient &B)
Reset the second vector in the inner product.
virtual void EvalSymmetric(Vector &K, ElementTransformation &T, const IntegrationPoint &ip)
(DEPRECATED) Evaluate a symmetric matrix coefficient.
Matrix coefficient defined as -a k x k x, for a vector k and scalar a.
DeltaCoefficient(double x, double y, double z, double s)
Construct a delta function scaled by s and centered at (x,y,z)
VectorCoefficient * GetKCoef() const
Return the vector factor.
void SetTime(double t)
Set the time for internally stored coefficients.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
MatrixFunctionCoefficient(int dim, std::function< void(const Vector &, DenseMatrix &)> F, Coefficient *q=nullptr)
Define a time-independent square matrix coefficient from a std function.
int GetWidth() const
Get the width of the matrix.
const Vector & GetA() const
Return the first vector constant.
RestrictedCoefficient(Coefficient &c_, Array< int > &attr)
Construct with a parent coefficient and an array with ones marking the attributes on which this coeff...
Matrix coefficient defined as the linear combination of two matrices.
Class to represent a coefficient evaluated at quadrature points.
DeltaCoefficient & GetDeltaCoefficient()
Return the associated scalar DeltaCoefficient.
void SetACoef(MatrixCoefficient &A)
Reset the matrix coefficient.
DeltaCoefficient(double x, double s)
Construct a delta function scaled by s and centered at (x,0.0,0.0)
DeterminantCoefficient(MatrixCoefficient &A)
Construct with the matrix.
SumCoefficient(double A, Coefficient &B, double alpha_=1.0, double beta_=1.0)
Constructor with one coefficient. Result is alpha_ * A + beta_ * B.
Matrix coefficient defined as the inverse a matrix coefficient.
void SetACoef(Coefficient &A)
Reset the first term in the product.
ScalarMatrixProductCoefficient(double A, MatrixCoefficient &B)
Constructor with one coefficient. Result is A*B.
FunctionCoefficient(std::function< double(const Vector &, double)> TDF)
Define a time-dependent coefficient from a std function.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
A VectorDeltaFunction cannot be evaluated. Calling this method will cause an MFEM error...
const QuadratureFunction & GetQuadFunction() const
void ZeroCoefficient(int attr)
Remove a single VectorCoefficient for a particular attribute.
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
RatioCoefficient(double A, Coefficient &B)
void Project(Coefficient &coeff)
Evaluate the given Coefficient at the quadrature points defined by qs.
Coefficient * GetCoeff(int i, int j)
Get the coefficient located at (i,j) in the matrix.
Matrix coefficient defined as the outer product of two vector coefficients.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
Coefficient * GetACoef() const
Return the first term in the product.
PWVectorCoefficient(int vd)
Constructs a piecewise vector coefficient of dimension vd.
void SetBCoef(Coefficient &B)
Reset the second term in the product.
MFEM_DEPRECATED FunctionCoefficient(double(*tdf)(Vector &, double))
(DEPRECATED) Define a time-dependent coefficient from a C-function
void SetACoef(MatrixCoefficient &A)
Reset the first matrix coefficient.
const GridFunction * GetGridFunction() const
Returns a pointer to the grid function in this Coefficient.
double Eval(int i, int j, ElementTransformation &T, const IntegrationPoint &ip)
void SetAConst(double A)
Reset the scalar factor as a constant.
double(* tdf)(double)
Derived vector coefficient that has the value of the parent vector where it is active and is zero oth...
InnerProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
Construct with the two vector coefficients. Result is .
void SetTime(double t)
Set the time for internally stored coefficients.
void SetTime(double t)
Set the time for internally stored coefficients.
MatrixCoefficient(int dim, bool symm=false)
Construct a dim x dim matrix coefficient.
void SetTime(double t)
Set the time for internally stored coefficients.
void SetTime(double t)
Set the time for internally stored coefficients.
Quadrature function coefficient which requires that the quadrature rules used for this coefficient be...
void SetAConst(double A)
Reset the scalar factor as a constant.
void SetTime(double t)
Set the time for internally stored coefficients.
void UpdateCoefficient(int attr, Coefficient &coef)
Replace a single Coefficient for a particular attribute.
double ComputeLpNorm(double p, Coefficient &coeff, Mesh &mesh, const IntegrationRule *irs[])
Compute the Lp norm of a function f. .
std::function< double(const Vector &, double)> TDFunction
void SetTime(double t)
Set the time for internally stored coefficients.
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
VectorConstantCoefficient(const Vector &v)
Construct the coefficient with constant vector v.
ConstantCoefficient(double c=1.0)
c is value of constant function
Definition: coefficient.hpp:90
void SetTime(double t)
Set the time for internally stored coefficients.
void SetBCoef(VectorCoefficient &B)
Reset the vector coefficient.
Store constants using only vdim entries.
void SetBCoef(VectorCoefficient &B)
Reset the second vector in the outer product.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
A matrix coefficient with an optional scalar coefficient multiplier q. The matrix function can either...
const double * Center()
double GetTime()
Get the time for time dependent coefficients.
void SetA(const Vector &A_)
Reset the first vector as a constant.
VectorCoefficient * GetACoef() const
Return the first vector coefficient.
Derived coefficient that takes the value of the parent coefficient for the active attributes and is z...
QuadratureFunctionCoefficient(QuadratureFunction &qf)
Constructor with a quadrature function as input.
void mfem_error(const char *msg)
Function called when an error is encountered. Used by the macros MFEM_ABORT, MFEM_ASSERT, MFEM_VERIFY.
Definition: error.cpp:154
void SetDeltaCenter(const Vector &center)
Coefficient * Weight()
See SetWeight() for description of the weight Coefficient.
double b
Definition: lissajous.cpp:42
VectorCoefficient * GetBCoef() const
Return the second vector of the product.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
RatioCoefficient(Coefficient &A, double B)
MatrixFunctionCoefficient(int dim, std::function< void(const Vector &, Vector &)> SymmF, Coefficient *q=NULL)
Define a time-independent symmetric square matrix coefficient from a std function.
MatrixRestrictedCoefficient(MatrixCoefficient &mc, Array< int > &attr)
Construct with a parent matrix coefficient and an array of zeros and ones representing the attributes...
void UpdateCoefficients(const Array< int > &attr, const Array< MatrixCoefficient *> &coefs)
Replace a set of coefficients.
std::function< double(const Vector &)> Function
MatrixCoefficient * GetACoef() const
Return the matrix coefficient.
Enable all above compressions.
Vector coefficient defined as the Gradient of a scalar GridFunction.
Vector coefficient defined as the Curl of a vector GridFunction.
OuterProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
Construct with two vector coefficients. Result is .
void SetGridFunction(const GridFunction *gf)
Set the vector grid function.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void Set(int i, Coefficient *c, bool own=true)
Sets coefficient in the vector.
void SetBeta(double beta_)
Reset the factor in front of the second term in the linear combination.
DeltaCoefficient()
Construct a unit delta function centered at (0.0,0.0,0.0)
void SetTol(double tol_)
Set the tolerance used during projection onto GridFunction to identify the Mesh vertex where the Cent...
void ZeroCoefficient(int attr)
Remove a single MatrixCoefficient for a particular attribute.
void SetACoef(VectorCoefficient &A)
Reset the first vector in the inner product.
PWMatrixCoefficient(int dim, const Array< int > &attr, const Array< MatrixCoefficient *> &coefs, bool symm=false)
Construct the coefficient using arrays describing the pieces.
double GetBeta() const
Return the factor in front of the second matrix coefficient.
void MakeRef(const QuadratureFunction &qf_)
Make this vector a reference to the given QuadratureFunction.
MatrixCoefficient * GetBCoef() const
Return the second matrix coefficient.
void SetAlpha(double alpha_)
Reset the factor in front of the first matrix coefficient.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void UpdateCoefficient(int attr, MatrixCoefficient &coef)
Replace a single coefficient for a particular attribute.
virtual ~VectorArrayCoefficient()
Destroys vector coefficient.
double GetAlpha() const
Return the factor in front of the first vector coefficient.
void SetTime(double t)
Set the time for internally stored coefficients.
void SetACoef(VectorCoefficient &A_)
Reset the first vector coefficient.
void SetTime(double t)
Set the time for internally stored coefficients.
void SetBeta(double beta_)
Reset the factor in front of the second vector coefficient as a constant.
void SetTime(double t)
Set the time for internally stored coefficients.
VectorFunctionCoefficient(int dim, std::function< void(const Vector &, Vector &)> F, Coefficient *q=nullptr)
Define a time-independent vector coefficient from a std function.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
MatrixCoefficient * GetACoef() const
Return the first matrix coefficient.
Store the triangular part of symmetric matrices.
int GetVDim()
Returns dimension of the vector.
VectorCoefficient * GetBCoef() const
Return the second vector coefficient in the inner product.
VectorCoefficient * GetBCoef() const
Return the second vector coefficient in the outer product.
VectorCoefficient * GetACoef() const
Return the first term in the product.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void SetACoef(Coefficient &A)
Reset the scalar factor.
PWMatrixCoefficient(int dim, bool symm=false)
Constructs a piecewise matrix coefficient of dimension dim by dim.
VectorCoefficient * GetBCoef() const
Return the vector coefficient.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void SetACoef(VectorCoefficient &A)
Reset the first term in the product.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient in the element described by T at the point ip, storing the result in ...
VectorCoefficient * GetACoef() const
Return the first vector of the product.
Coefficient * GetACoef() const
Return the numerator of the ratio.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the determinant coefficient at ip.
VectorDeltaCoefficient(const Vector &dir_)
Construct with a Vector object representing the direction and a unit delta function centered at (0...
void SetTime(double t)
Set the time for internally stored coefficients.
void SetComponent(int index_, int length_)
void SetBCoef(VectorCoefficient &B)
Reset the vector factor.
void operator=(double c)
Set the constants for all attributes to constant c.
VectorCoefficient DiagonalMatrixCoefficient
A general vector function coefficient.
Abstract base class for QuadratureSpace and FaceQuadratureSpace.
Definition: qspace.hpp:24
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
Definition: coefficient.hpp:93
void GetDeltaCenter(Vector &center)
Write the center of the delta function into center.
MatrixFunctionCoefficient(const DenseMatrix &m, Coefficient &q)
Define a constant matrix coefficient times a scalar Coefficient.
Scalar coefficient defined as the Divergence of a vector GridFunction.
int GetVDim() const
Return the number of values per quadrature point.
MatrixVectorProductCoefficient MatVecCoefficient
Convenient alias for the MatrixVectorProductCoefficient.
virtual void ProjectSymmetric(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
void SetACoef(Coefficient &A)
Reset the base coefficient.
void SetACoef(MatrixCoefficient &A)
Reset the matrix coefficient.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
VectorDeltaCoefficient(const Vector &dir_, double x, double y, double z, double s)
Construct with a Vector object representing the direction and a delta function scaled by s and center...
void SetGridFunction(const GridFunction *gf)
Set the scalar grid function.
Coefficient * GetACoef() const
Return the scalar factor.
Coefficient * GetBCoef() const
Return the second term in the linear combination.
VectorDeltaCoefficient(const Vector &dir_, double x, double s)
Construct with a Vector object representing the direction and a delta function scaled by s and center...
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:41
void SetBCoef(VectorCoefficient &B)
Reset the second vector in the product.
virtual ~Coefficient()
Definition: coefficient.hpp:79
VectorQuadratureFunctionCoefficient(QuadratureFunction &qf)
Constructor with a quadrature function as input.
Vector coefficient defined as a product of scalar and vector coefficients.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
void UpdateConstants(Vector &c)
Update the constants with vector c.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
void SetBetaCoef(Coefficient &B_)
Reset the factor in front of the second vector coefficient.
Base class for Matrix Coefficients that optionally depend on time and space.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
CoefficientStorage storage
Storage optimizations (see CoefficientStorage).
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
Definition: coefficient.cpp:76
Coefficient * GetCoeff(int i)
Returns i&#39;th coefficient.
virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result as ...
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient at ip.
virtual void EvalDelta(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Return the specified direction vector multiplied by the value returned by DeltaCoefficient::EvalDelta...
MatrixProductCoefficient(MatrixCoefficient &A, MatrixCoefficient &B)
Construct with the two coefficients. Result is A * B.
double GetExponent() const
Return the exponent.
double ComputeGlobalLpNorm(double p, Coefficient &coeff, ParMesh &pmesh, const IntegrationRule *irs[])
Compute the global Lp norm of a function f. .
void SetAConst(double A)
Reset the first term in the linear combination as a constant.
PWConstCoefficient(int NumOfSubD=0)
Constructs a piecewise constant coefficient in NumOfSubD subdomains.
Coefficient * GetACoef() const
Return the scalar factor.
PowerCoefficient(Coefficient &A, double p_)
Construct with a coefficient and a constant power p_. Result is A^p.
MatrixCoefficient * GetBCoef() const
Return the matrix factor.
PWMatrixCoefficient(int h, int w, bool symm=false)
Constructs a piecewise matrix coefficient of dimension h by w.
virtual double EvalDelta(ElementTransformation &T, const IntegrationPoint &ip)
The value of the function assuming we are evaluating at the delta center.
SymmetricMatrixFunctionCoefficient(int dim, std::function< void(const Vector &, double, DenseSymmetricMatrix &)> TDF, Coefficient *q=nullptr)
Define a time-dependent square matrix coefficient from a std function.
SymmetricMatrixConstantCoefficient(const DenseSymmetricMatrix &m)
Construct using matrix m for the constant.
Vector coefficient defined as a product of a matrix coefficient and a vector coefficient.
A piecewise matrix-valued coefficient with the pieces keyed off the element attribute numbers...
double GetAlpha() const
Return the factor in front of the first matrix coefficient.
Derived matrix coefficient that has the value of the parent matrix coefficient where it is active and...
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
Vector coefficient defined by a scalar DeltaCoefficient and a constant vector direction.
VectorCoefficient * GetBCoef() const
Return the vector factor.
PWCoefficient()
Constructs a piecewise coefficient.
PWVectorCoefficient(int vd, const Array< int > &attr, const Array< VectorCoefficient *> &coefs)
Construct the coefficient using arrays describing the pieces.
Store the coefficient as a full QuadratureFunction.
double a
Definition: lissajous.cpp:41
MatrixArrayCoefficient(int dim)
Construct a coefficient matrix of dimensions dim * dim. The actual coefficients still need to be adde...
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the gradient vector coefficient at ip.
A piecewise constant coefficient with the constants keyed off the element attribute numbers...
void SetACoef(VectorCoefficient &A)
Reset the vector coefficient.
void SetACoef(MatrixCoefficient &A)
Reset the matrix coefficient.
Coefficient * GetACoef() const
Return the base coefficient.
Class for integration point with weight.
Definition: intrules.hpp:25
void SetTime(double t)
Set the time for internally stored coefficients.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the scalar divergence coefficient at ip.
void SetAlpha(double alpha_)
Reset the factor in front of the first vector coefficient as a constant.
void SetAConst(double A)
Reset the numerator in the ratio as a constant.
void SetScale(double s_)
Set the scale value multiplying the delta function.
void ProjectTranspose(MatrixCoefficient &coeff)
Project the transpose of coeff.
Base class for symmetric matrix coefficients that optionally depend on time and space.
DivergenceGridFunctionCoefficient(const GridFunction *gf)
Construct the coefficient with a vector grid function gf. The grid function is not owned by the coeff...
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
SymmetricMatrixCoefficient(int dimension)
Construct a dim x dim matrix coefficient.
TransformedCoefficient(Coefficient *q, double(*F)(double))
MatrixCoefficient * GetACoef() const
Return the matrix coefficient.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
void SetTime(double t)
Set the time for internally stored coefficients.
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
Vector coefficient defined as a normalized vector field (returns v/|v|)
VectorCoefficient(int vd)
Initialize the VectorCoefficient with vector dimension vd.
int GetNConst()
Returns the number of constants representing different attributes.
int dim
Definition: ex24.cpp:53
void SetDirection(const Vector &d_)
void SetAConst(double A)
Reset the scalar factor as a constant.
VectorFunctionCoefficient(int dim, std::function< void(const Vector &, double, Vector &)> TDF, Coefficient *q=nullptr)
Define a time-dependent vector coefficient from a std function.
const GridFunction * GetGridFunction() const
Get the internal GridFunction.
void SetACoef(VectorCoefficient &A)
Reset the first vector in the product.
const DenseSymmetricMatrix & GetMatrix()
Return a reference to the constant matrix.
VectorCoefficient * GetACoef() const
Return the first vector coefficient in the outer product.
QuadratureFunction * qf
Internal QuadratureFunction (owned, may be NULL).
void UpdateCoefficients(const Array< int > &attr, const Array< Coefficient *> &coefs)
Replace a set of coefficients.
Coefficient * GetBetaCoef() const
Return the factor in front of the second vector coefficient.
double GetAConst() const
Return the scalar factor.
void SetGridFunction(const GridFunction *gf)
Set the grid function for this coefficient. Also sets the Vector dimension to match that of the gf...
Matrix coefficient defined as the transpose a matrix coefficient.
void Copy(Array &copy) const
Create a copy of the internal array to the provided copy.
Definition: array.hpp:861
void ZeroCoefficient(int attr)
Remove a single Coefficient for a particular attribute.
A piecewise vector-valued coefficient with the pieces keyed off the element attribute numbers...
constexpr int dimension
This example only works in 3D. Kernels for 2D are not implemented.
Definition: hooke.cpp:45
void SetBConst(double B)
Reset the denominator in the ratio as a constant.
void SetBCoef(Coefficient &B)
Reset the second term in the linear combination.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
A matrix coefficient that is constant in space and time.
DeltaCoefficient(double x, double y, double s)
Construct a delta function scaled by s and centered at (x,y,0.0)
void SetTime(double t)
Set the time for internally stored coefficients.
const GridFunction * GetGridFunction() const
Get the vector grid function.
Coefficient * GetAlphaCoef() const
Return the factor in front of the first vector coefficient.
MatrixSumCoefficient(MatrixCoefficient &A, MatrixCoefficient &B, double alpha_=1.0, double beta_=1.0)
Construct with the two coefficients. Result is alpha_ * A + beta_ * B.
Scalar coefficient defined as a scalar raised to a power.
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
Definition: coefficient.cpp:51
const GridFunction * GetGridFunction() const
Get the vector grid function.
double GetBeta() const
Return the factor in front of the second term in the linear combination.
RefCoord t[3]
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the coefficient in the element described by T at the point ip.
VectorArrayCoefficient(int dim)
Construct vector of dim coefficients. The actual coefficients still need to be added with Set()...
virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
double GetTime()
Get the time for time dependent coefficients.
Definition: coefficient.hpp:53
DenseSymmetricMatrix mat
Internal matrix used when evaluating this coefficient as a DenseMatrix.
MatrixCoefficient * GetACoef() const
Return the matrix coefficient.
Scalar coefficient defined as the determinant of a matrix coefficient.
const double alpha
Definition: ex15.cpp:369
A general function coefficient.
void SetACoef(VectorCoefficient &A)
Reset the first vector in the outer product.
ProductCoefficient(Coefficient &A, Coefficient &B)
Constructor with two coefficients. Result is A * B.
Vector data type.
Definition: vector.hpp:60
TransformedCoefficient(Coefficient *q1, Coefficient *q2, double(*F)(double, double))
virtual void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points.
virtual void SetTime(double t)
Set the time for time dependent coefficients.
double GetBConst() const
Return the denominator of the ratio.
GridFunctionCoefficient(const GridFunction *gf, int comp=1)
virtual void Project(QuadratureFunction &qf, bool transpose=false)
Fill the QuadratureFunction qf by evaluating the coefficient at the quadrature points. The matrix will be transposed or not according to the boolean argument transpose.
Vector coefficient defined by a vector GridFunction.
VectorCoefficient * GetACoef() const
Return the vector coefficient.
Coefficient * GetBCoef() const
Return the second term in the product.
void Project(QuadratureFunction &qf)
Fill the QuadratureFunction qf with the constant value.
Definition: coefficient.cpp:71
const Vector & GetB() const
Return the second vector constant.
CoefficientVector(QuadratureSpaceBase &qs_, CoefficientStorage storage_=CoefficientStorage::FULL)
Create an empty CoefficientVector.
VectorRestrictedCoefficient(VectorCoefficient &vc, Array< int > &attr)
Construct with a parent vector coefficient and an array of zeros and ones representing the attributes...
void SetTime(double t)
Set the time for internally stored coefficients.
VectorRotProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
Constructor with two vector coefficients. Result is .
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
void SetACoef(Coefficient &A)
Reset the scalar factor.
RefCoord s[3]
MatrixConstantCoefficient(const DenseMatrix &m)
Construct using matrix m for the constant.
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient at ip.
QuadratureSpaceBase & qs
Associated QuadratureSpaceBase.
void SetAlpha(double alpha_)
Reset the factor in front of the first term in the linear combination.
Constant matrix coefficient defined as the identity of dimension d.
virtual void EvalSymmetric(Vector &K, ElementTransformation &T, const IntegrationPoint &ip)
(DEPRECATED) Evaluate the symmetric matrix coefficient at ip.
void SetTime(double t)
Set the time for internally stored coefficients.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient at ip.
void SetGridFunction(const GridFunction *gf)
Set the internal GridFunction.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
void SetBCoef(Coefficient &B)
Reset the denominator in the ratio.
double GetBeta() const
Return the factor in front of the second vector coefficient.
void SetTime(double t)
Set the time for internally stored coefficients.
double GetAConst() const
Return the scalar factor.
void SetACoef(MatrixCoefficient &A)
Reset the first matrix coefficient.
Class for parallel meshes.
Definition: pmesh.hpp:32
InverseMatrixCoefficient(MatrixCoefficient &A)
Construct with the matrix coefficient. Result is .
void SetTime(double t)
Set the time for internally stored coefficients.
Scalar coefficient defined as the linear combination of two scalar coefficients or a scalar and a sca...
ProductCoefficient(double A, Coefficient &B)
Constructor with one coefficient. Result is A * B.
VectorCoefficient * GetACoef() const
Return the first vector coefficient in the inner product.
Represents values or vectors of values at quadrature points on a mesh.
Definition: qfunction.hpp:23
void SetTime(double t)
Set the time for internally stored coefficients.
void SetTime(double t)
Set the time for internally stored coefficients.
void UpdateCoefficients(const Array< int > &attr, const Array< VectorCoefficient *> &coefs)
Replace a set of coefficients.
VectorDeltaCoefficient(const Vector &dir_, double x, double y, double s)
Construct with a Vector object representing the direction and a delta function scaled by s and center...
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
double GetAConst() const
Return the numerator of the ratio.
void SetACoef(Coefficient &A)
Reset the first term in the linear combination.
double GetAlpha() const
Return the factor in front of the first term in the linear combination.
void SetBCoef(MatrixCoefficient &B)
Reset the second matrix coefficient.
SumCoefficient(Coefficient &A, Coefficient &B, double alpha_=1.0, double beta_=1.0)
Constructor with two coefficients. Result is alpha_ * A + beta_ * B.
VectorCoefficient * GetBCoef() const
Return the second term in the product.
Coefficient * GetACoef() const
Return the scalar factor.
void SetTime(double t)
Set the time for internally stored coefficients.
ScalarVectorProductCoefficient(double A, VectorCoefficient &B)
Constructor with constant and vector coefficient. Result is A * B.
VectorGridFunctionCoefficient()
Construct an empty coefficient. Calling Eval() before the grid function is set will cause a segfault...
void SetTime(double t)
Set the time for internally stored coefficients.