MFEM  v4.6.0 Finite element discretization library
coefficient.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
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;
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
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.
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 /// A common base class for returning individual components of the domain's
262 /// Cartesian coordinates.
264 {
265 protected:
266  int comp;
267  mutable Vector transip;
268
269  /// @a comp_ index of the desired component (0 -> x, 1 -> y, 2 -> z)
270  CartesianCoefficient(int comp_) : comp(comp_), transip(3) {}
271
272 public:
273  /// Evaluate the coefficient at @a ip.
274  virtual double Eval(ElementTransformation &T,
275  const IntegrationPoint &ip);
276 };
277
278 /// Scalar coefficient which returns the x-component of the evaluation point
280 {
281 public:
283 };
284
285 /// Scalar coefficient which returns the y-component of the evaluation point
287 {
288 public:
290 };
291
292 /// Scalar coefficient which returns the z-component of the evaluation point
294 {
295 public:
297 };
298
299 /// Scalar coefficient which returns the radial distance from the axis of
300 /// the evaluation point in the cylindrical coordinate system
302 {
303 private:
304  mutable Vector transip;
305
306 public:
308
309  /// Evaluate the coefficient at @a ip.
310  virtual double Eval(ElementTransformation &T,
311  const IntegrationPoint &ip);
312 };
313
314 /// Scalar coefficient which returns the angular position or azimuth (often
315 /// denoted by theta) of the evaluation point in the cylindrical coordinate
316 /// system
318 {
319 private:
320  mutable Vector transip;
321
322 public:
324
325  /// Evaluate the coefficient at @a ip.
326  virtual double Eval(ElementTransformation &T,
327  const IntegrationPoint &ip);
328 };
329
330 /// Scalar coefficient which returns the height or altitude of
331 /// the evaluation point in the cylindrical coordinate system
333
334 /// Scalar coefficient which returns the radial distance from the origin of
335 /// the evaluation point in the spherical coordinate system
337 {
338 private:
339  mutable Vector transip;
340
341 public:
343
344  /// Evaluate the coefficient at @a ip.
345  virtual double Eval(ElementTransformation &T,
346  const IntegrationPoint &ip);
347 };
348
349 /// Scalar coefficient which returns the azimuthal angle (often denoted by phi)
350 /// of the evaluation point in the spherical coordinate system
352 {
353 private:
354  mutable Vector transip;
355
356 public:
358
359  /// Evaluate the coefficient at @a ip.
360  virtual double Eval(ElementTransformation &T,
361  const IntegrationPoint &ip);
362 };
363
364 /// Scalar coefficient which returns the polar angle (often denoted by theta)
365 /// of the evaluation point in the spherical coordinate system
367 {
368 private:
369  mutable Vector transip;
370
371 public:
372  SphericalPolarCoefficient() : transip(3) {}
373
374  /// Evaluate the coefficient at @a ip.
375  virtual double Eval(ElementTransformation &T,
376  const IntegrationPoint &ip);
377 };
378
379 class GridFunction;
380
381 /// Coefficient defined by a GridFunction. This coefficient is mesh dependent.
383 {
384 private:
385  const GridFunction *GridF;
386  int Component;
387
388 public:
389  GridFunctionCoefficient() : GridF(NULL), Component(1) { }
390  /** Construct GridFunctionCoefficient from a given GridFunction, and
391  optionally specify a component to use if it is a vector GridFunction. */
392  GridFunctionCoefficient (const GridFunction *gf, int comp = 1)
393  { GridF = gf; Component = comp; }
394
395  /// Set the internal GridFunction
396  void SetGridFunction(const GridFunction *gf) { GridF = gf; }
397
398  /// Get the internal GridFunction
399  const GridFunction * GetGridFunction() const { return GridF; }
400
401  /// Evaluate the coefficient at @a ip.
402  virtual double Eval(ElementTransformation &T,
403  const IntegrationPoint &ip);
404
405  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
407  ///
408  /// This function uses the efficient QuadratureFunction::ProjectGridFunction
409  /// to fill the QuadratureFunction.
411 };
412
413
414 /** @brief A coefficient that depends on 1 or 2 parent coefficients and a
415  transformation rule represented by a C-function.
416
417  \f$C(x,t) = T(Q1(x,t)) \f$ or \f$C(x,t) = T(Q1(x,t), Q2(x,t)) \f$
418
419  where T is the transformation rule, and Q1/Q2 are the parent coefficients.*/
421 {
422 private:
423  Coefficient * Q1;
424  Coefficient * Q2;
425  double (*Transform1)(double);
426  double (*Transform2)(double,double);
427
428 public:
429  TransformedCoefficient (Coefficient * q,double (*F)(double))
430  : Q1(q), Transform1(F) { Q2 = 0; Transform2 = 0; }
432  double (*F)(double,double))
433  : Q1(q1), Q2(q2), Transform2(F) { Transform1 = 0; }
434
435  /// Set the time for internally stored coefficients
436  void SetTime(double t);
437
438  /// Evaluate the coefficient at @a ip.
439  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
440 };
441
442 /** @brief Delta function coefficient optionally multiplied by a weight
443  coefficient and a scaled time dependent C-function.
444
445  \f$F(x,t) = w(x,t) s T(t) d(x - xc) \f$
446
447  where w is the optional weight coefficient, @a s is a scale factor
448  T is an optional time-dependent function and d is a delta function.
449
450  WARNING this cannot be used as a normal coefficient. The usual Eval
451  method is disabled. */
453 {
454 protected:
455  double center[3], scale, tol;
457  int sdim;
458  double (*tdf)(double);
459
460 public:
461
462  /// Construct a unit delta function centered at (0.0,0.0,0.0)
464  {
465  center[0] = center[1] = center[2] = 0.; scale = 1.; tol = 1e-12;
466  weight = NULL; sdim = 0; tdf = NULL;
467  }
468
469  /// Construct a delta function scaled by @a s and centered at (x,0.0,0.0)
470  DeltaCoefficient(double x, double s)
471  {
472  center[0] = x; center[1] = 0.; center[2] = 0.; scale = s; tol = 1e-12;
473  weight = NULL; sdim = 1; tdf = NULL;
474  }
475
476  /// Construct a delta function scaled by @a s and centered at (x,y,0.0)
477  DeltaCoefficient(double x, double y, double s)
478  {
479  center[0] = x; center[1] = y; center[2] = 0.; scale = s; tol = 1e-12;
480  weight = NULL; sdim = 2; tdf = NULL;
481  }
482
483  /// Construct a delta function scaled by @a s and centered at (x,y,z)
484  DeltaCoefficient(double x, double y, double z, double s)
485  {
486  center[0] = x; center[1] = y; center[2] = z; scale = s; tol = 1e-12;
487  weight = NULL; sdim = 3; tdf = NULL;
488  }
489
490  /// Set the time for internally stored coefficients
491  void SetTime(double t);
492
493  /// Set the center location of the delta function.
494  void SetDeltaCenter(const Vector& center);
495
496  /// Set the scale value multiplying the delta function.
497  void SetScale(double s_) { scale = s_; }
498
499  /// Set a time-dependent function that multiplies the Scale().
500  void SetFunction(double (*f)(double)) { tdf = f; }
501
502  /** @brief Set the tolerance used during projection onto GridFunction to
503  identify the Mesh vertex where the Center() of the delta function
504  lies. (default 1e-12)*/
505  void SetTol(double tol_) { tol = tol_; }
506
507  /// Set a weight Coefficient that multiplies the DeltaCoefficient.
508  /** The weight Coefficient multiplies the value returned by EvalDelta() but
509  not the value returned by Scale().
510  The weight Coefficient is also used as the L2-weight function when
511  projecting the DeltaCoefficient onto a GridFunction, so that the weighted
512  integral of the projection is exactly equal to the Scale(). */
513  void SetWeight(Coefficient *w) { weight = w; }
514
515  /// Return a pointer to a c-array representing the center of the delta
516  /// function.
517  const double *Center() { return center; }
518
519  /** @brief Return the scale factor times the optional time dependent
520  function. Returns \f$s T(t) \f$ with \f$T(t) = 1 \f$ when
521  not set by the user. */
522  double Scale() { return tdf ? (*tdf)(GetTime())*scale : scale; }
523
524  /// Return the tolerance used to identify the mesh vertices
526
527  /// See SetWeight() for description of the weight Coefficient.
528  Coefficient *Weight() { return weight; }
529
530  /// Write the center of the delta function into @a center.
532
533  /// The value of the function assuming we are evaluating at the delta center.
534  virtual double EvalDelta(ElementTransformation &T, const IntegrationPoint &ip);
535  /** @brief A DeltaFunction cannot be evaluated. Calling this method will
536  cause an MFEM error, terminating the application. */
537  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
538  { mfem_error("DeltaCoefficient::Eval"); return 0.; }
539  virtual ~DeltaCoefficient() { delete weight; }
540 };
541
542 /** @brief Derived coefficient that takes the value of the parent coefficient
543  for the active attributes and is zero otherwise. */
545 {
546 private:
547  Coefficient *c;
548  Array<int> active_attr;
549
550 public:
551  /** @brief Construct with a parent coefficient and an array with
552  ones marking the attributes on which this coefficient should be
553  active. */
555  { c = &c_; attr.Copy(active_attr); }
556
557  /// Set the time for internally stored coefficients
558  void SetTime(double t);
559
560  /// Evaluate the coefficient at @a ip.
561  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
562  { return active_attr[T.Attribute-1] ? c->Eval(T, ip, GetTime()) : 0.0; }
563 };
564
565 /// Base class for vector Coefficients that optionally depend on time and space.
567 {
568 protected:
569  int vdim;
570  double time;
571
572 public:
573  /// Initialize the VectorCoefficient with vector dimension @a vd.
574  VectorCoefficient(int vd) { vdim = vd; time = 0.; }
575
576  /// Set the time for time dependent coefficients
577  virtual void SetTime(double t) { time = t; }
578
579  /// Get the time for time dependent coefficients
580  double GetTime() { return time; }
581
582  /// Returns dimension of the vector.
583  int GetVDim() { return vdim; }
584
585  /** @brief Evaluate the vector coefficient in the element described by @a T
586  at the point @a ip, storing the result in @a V. */
587  /** @note When this method is called, the caller must make sure that the
588  IntegrationPoint associated with @a T is the same as @a ip. This can be
589  achieved by calling T.SetIntPoint(&ip). */
590  virtual void Eval(Vector &V, ElementTransformation &T,
591  const IntegrationPoint &ip) = 0;
592
593  /** @brief Evaluate the vector coefficient in the element described by @a T
594  at all points of @a ir, storing the result in @a M. */
595  /** The dimensions of @a M are GetVDim() by ir.GetNPoints() and they must be
596  set by the implementation of this method.
597
598  The general implementation provided by the base class (using the Eval
599  method for one IntegrationPoint at a time) can be overloaded for more
600  efficient implementation.
601
602  @note The IntegrationPoint associated with @a T is not used, and this
603  method will generally modify this IntegrationPoint associated with @a T.
604  */
605  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
606  const IntegrationRule &ir);
607
608  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
610  ///
611  /// The @a vdim of the VectorCoefficient should be equal to the @a vdim of
614
615  virtual ~VectorCoefficient() { }
616 };
617
618
619 /// Vector coefficient that is constant in space and time.
621 {
622 private:
623  Vector vec;
624 public:
625  /// Construct the coefficient with constant vector @a v.
627  : VectorCoefficient(v.Size()), vec(v) { }
629
630  /// Evaluate the vector coefficient at @a ip.
631  virtual void Eval(Vector &V, ElementTransformation &T,
632  const IntegrationPoint &ip) { V = vec; }
633
634  /// Return a reference to the constant vector in this class.
635  const Vector& GetVec() const { return vec; }
636 };
637
638 /** @brief A piecewise vector-valued coefficient with the pieces keyed off the
639  element attribute numbers.
640
641  A value of zero will be returned for any missing attribute numbers.
642
643  This object will not assume ownership of any VectorCoefficient
644  objects passed to it. Consequently, the caller must ensure that
645  the individual VectorCoefficient objects are not deleted while
646  this PWVectorCoefficient is still in use.
647
648  \note The keys may either be domain attribute numbers or boundary
649  attribute numbers. If the PWVectorCoefficient is used with a
650  domain integrator the keys are assumed to be domain attribute
651  numbers. Similarly, if the PWVectorCoefficient is used with a
652  boundary integrator the keys are assumed to be boundary attribute
653  numbers.
654 */
656 {
657 private:
658  /** Internal data structure to store pointers to the appropriate
659  coefficients for different regions of the mesh. The keys used
660  in the map are the mesh attribute numbers (either element
661  attribute or boundary element attribute depending upon
662  context). The values returned for any missing attributes will
663  be zero. The coefficient pointers may be NULL in which case a
664  value of zero is returned.
665
666  The VectorCoefficient objects contained in this map are NOT
667  owned by this PWVectorCoefficient object. This means that they
668  will not be deleted when this object is deleted also the caller
669  must ensure that the various VectorCoefficient objects are not
670  deleted while this PWVectorCoefficient is still needed.
671  */
672  std::map<int, VectorCoefficient*> pieces;
673
674  /** Convenience function to check for compatible array lengths,
675  loop over the arrays, and add their attribute/VectorCoefficient
676  pairs to the internal data structure.
677  */
678  void InitMap(const Array<int> & attr,
679  const Array<VectorCoefficient*> & coefs);
680
681 public:
682
683  /// Constructs a piecewise vector coefficient of dimension vd
684  explicit PWVectorCoefficient(int vd): VectorCoefficient(vd) {}
685
686  /// Construct the coefficient using arrays describing the pieces
687  /** \param vd - dimension of the vector-valued result
688  \param attr - an array of attribute numbers for each piece
689  \param coefs - the corresponding array of VectorCoefficient pointers
690  Any missing attributes or NULL coefficient pointers will result in a
691  zero vector being returned for that attribute.
692
693  \note Ownership of the VectorCoefficient objects will NOT be
694  transferred to this object.
695  */
696  PWVectorCoefficient(int vd, const Array<int> & attr,
697  const Array<VectorCoefficient*> & coefs)
698  : VectorCoefficient(vd) { InitMap(attr, coefs); }
699
700  /// Set the time for time dependent coefficients
701  virtual void SetTime(double t);
702
703  /// Replace a set of coefficients
704  void UpdateCoefficients(const Array<int> & attr,
705  const Array<VectorCoefficient*> & coefs)
706  { InitMap(attr, coefs); }
707
708  /// Replace a single Coefficient for a particular attribute
709  void UpdateCoefficient(int attr, VectorCoefficient & coef);
710
711  /// Remove a single VectorCoefficient for a particular attribute
712  void ZeroCoefficient(int attr)
713  { pieces.erase(attr); }
714
715  /// Evaluate the coefficient.
716  virtual void Eval(Vector &V, ElementTransformation &T,
717  const IntegrationPoint &ip);
719 };
720
721 /// A vector coefficient which returns the physical location of the
722 /// evaluation point in the Cartesian coordinate system.
724 {
725 public:
726
728
730  /// Evaluate the vector coefficient at @a ip.
731  virtual void Eval(Vector &V, ElementTransformation &T,
732  const IntegrationPoint &ip);
733
735 };
736
737 /// A general vector function coefficient
739 {
740 private:
741  std::function<void(const Vector &, Vector &)> Function;
742  std::function<void(const Vector &, double, Vector &)> TDFunction;
743  Coefficient *Q;
744
745 public:
746  /// Define a time-independent vector coefficient from a std function
747  /** \param dim - the size of the vector
748  \param F - time-independent function
749  \param q - optional scalar Coefficient to scale the vector coefficient */
751  std::function<void(const Vector &, Vector &)> F,
752  Coefficient *q = nullptr)
753  : VectorCoefficient(dim), Function(std::move(F)), Q(q)
754  { }
755
756  /// Define a time-dependent vector coefficient from a std function
757  /** \param dim - the size of the vector
758  \param TDF - time-dependent function
759  \param q - optional scalar Coefficient to scale the vector coefficient */
761  std::function<void(const Vector &, double, Vector &)> TDF,
762  Coefficient *q = nullptr)
763  : VectorCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
764  { }
765
767  /// Evaluate the vector coefficient at @a ip.
768  virtual void Eval(Vector &V, ElementTransformation &T,
769  const IntegrationPoint &ip);
770
772 };
773
774 /** @brief Vector coefficient defined by an array of scalar coefficients.
775  Coefficients that are not set will evaluate to zero in the vector. This
776  object takes ownership of the array of coefficients inside it and deletes
777  them at object destruction. */
779 {
780 private:
781  Array<Coefficient*> Coeff;
782  Array<bool> ownCoeff;
783
784 public:
785  /** @brief Construct vector of dim coefficients. The actual coefficients
786  still need to be added with Set(). */
787  explicit VectorArrayCoefficient(int dim);
788
789  /// Set the time for internally stored coefficients
790  void SetTime(double t);
791
792  /// Returns i'th coefficient.
793  Coefficient* GetCoeff(int i) { return Coeff[i]; }
794
795  /// Returns the entire array of coefficients.
796  Coefficient **GetCoeffs() { return Coeff; }
797
798  /// Sets coefficient in the vector.
799  void Set(int i, Coefficient *c, bool own=true);
800
801  /// Evaluates i'th component of the vector of coefficients and returns the
802  /// value.
803  double Eval(int i, ElementTransformation &T, const IntegrationPoint &ip)
804  { return Coeff[i] ? Coeff[i]->Eval(T, ip, GetTime()) : 0.0; }
805
807  /** @brief Evaluate the coefficient. Each element of vector V comes from the
808  associated array of scalar coefficients. */
809  virtual void Eval(Vector &V, ElementTransformation &T,
810  const IntegrationPoint &ip);
811
812  /// Destroys vector coefficient.
813  virtual ~VectorArrayCoefficient();
814 };
815
816 /// Vector coefficient defined by a vector GridFunction
818 {
819 protected:
821
822 public:
823  /** @brief Construct an empty coefficient. Calling Eval() before the grid
824  function is set will cause a segfault. */
826
827  /** @brief Construct the coefficient with grid function @a gf. The
828  grid function is not owned by the coefficient. */
830
831  /** @brief Set the grid function for this coefficient. Also sets the Vector
832  dimension to match that of the @a gf. */
833  void SetGridFunction(const GridFunction *gf);
834
835  /// Returns a pointer to the grid function in this Coefficient
836  const GridFunction * GetGridFunction() const { return GridFunc; }
837
838  /// Evaluate the vector coefficient at @a ip.
839  virtual void Eval(Vector &V, ElementTransformation &T,
840  const IntegrationPoint &ip);
841
842  /** @brief Evaluate the vector coefficients at all of the locations in the
843  integration rule and write the vectors into the columns of matrix @a
844  M. */
845  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
846  const IntegrationRule &ir);
847
848  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
850  ///
851  /// This function uses the efficient QuadratureFunction::ProjectGridFunction
852  /// to fill the QuadratureFunction.
854
856 };
857
858 /// Vector coefficient defined as the Gradient of a scalar GridFunction
860 {
861 protected:
863
864 public:
865
866  /** @brief Construct the coefficient with a scalar grid function @a gf. The
867  grid function is not owned by the coefficient. */
869
870  ///Set the scalar grid function.
871  void SetGridFunction(const GridFunction *gf);
872
873  ///Get the scalar grid function.
874  const GridFunction * GetGridFunction() const { return GridFunc; }
875
876  /// Evaluate the gradient vector coefficient at @a ip.
877  virtual void Eval(Vector &V, ElementTransformation &T,
878  const IntegrationPoint &ip);
879
880  /** @brief Evaluate the gradient vector coefficient at all of the locations
881  in the integration rule and write the vectors into columns of matrix @a
882  M. */
883  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
884  const IntegrationRule &ir);
885
887 };
888
889 /// Vector coefficient defined as the Curl of a vector GridFunction
891 {
892 protected:
894
895 public:
896  /** @brief Construct the coefficient with a vector grid function @a gf. The
897  grid function is not owned by the coefficient. */
899
900  /// Set the vector grid function.
901  void SetGridFunction(const GridFunction *gf);
902
903  /// Get the vector grid function.
904  const GridFunction * GetGridFunction() const { return GridFunc; }
905
907  /// Evaluate the vector curl coefficient at @a ip.
908  virtual void Eval(Vector &V, ElementTransformation &T,
909  const IntegrationPoint &ip);
910
912 };
913
914 /// Scalar coefficient defined as the Divergence of a vector GridFunction
916 {
917 protected:
919
920 public:
921  /** @brief Construct the coefficient with a vector grid function @a gf. The
922  grid function is not owned by the coefficient. */
924
925  /// Set the vector grid function.
926  void SetGridFunction(const GridFunction *gf) { GridFunc = gf; }
927
928  /// Get the vector grid function.
929  const GridFunction * GetGridFunction() const { return GridFunc; }
930
931  /// Evaluate the scalar divergence coefficient at @a ip.
932  virtual double Eval(ElementTransformation &T,
933  const IntegrationPoint &ip);
934
936 };
937
938 /** @brief Vector coefficient defined by a scalar DeltaCoefficient and a
939  constant vector direction.
940
941  WARNING this cannot be used as a normal coefficient. The usual Eval method
942  is disabled. */
944 {
945 protected:
948
949 public:
950  /// Construct with a vector of dimension @a vdim_.
952  : VectorCoefficient(vdim_), dir(vdim_), d() { }
953
954  /** @brief Construct with a Vector object representing the direction and a
955  unit delta function centered at (0.0,0.0,0.0) */
957  : VectorCoefficient(dir_.Size()), dir(dir_), d() { }
958
959  /** @brief Construct with a Vector object representing the direction and a
960  delta function scaled by @a s and centered at (x,0.0,0.0) */
961  VectorDeltaCoefficient(const Vector& dir_, double x, double s)
962  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,s) { }
963
964  /** @brief Construct with a Vector object representing the direction and a
965  delta function scaled by @a s and centered at (x,y,0.0) */
966  VectorDeltaCoefficient(const Vector& dir_, double x, double y, double s)
967  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,y,s) { }
968
969  /** @brief Construct with a Vector object representing the direction and a
970  delta function scaled by @a s and centered at (x,y,z) */
971  VectorDeltaCoefficient(const Vector& dir_, double x, double y, double z,
972  double s)
973  : VectorCoefficient(dir_.Size()), dir(dir_), d(x,y,z,s) { }
974
975  /// Set the time for internally stored coefficients
976  void SetTime(double t);
977
978  /// Replace the associated DeltaCoefficient with a new DeltaCoefficient.
979  /** The new DeltaCoefficient cannot have a specified weight Coefficient, i.e.
980  DeltaCoefficient::Weight() should return NULL. */
981  void SetDeltaCoefficient(const DeltaCoefficient& d_) { d = d_; }
982
983  /// Return the associated scalar DeltaCoefficient.
985
986  void SetScale(double s) { d.SetScale(s); }
987  void SetDirection(const Vector& d_);
988
989  void SetDeltaCenter(const Vector& center) { d.SetDeltaCenter(center); }
990  void GetDeltaCenter(Vector& center) { d.GetDeltaCenter(center); }
991
992  /** @brief Return the specified direction vector multiplied by the value
993  returned by DeltaCoefficient::EvalDelta() of the associated scalar
994  DeltaCoefficient. */
995  virtual void EvalDelta(Vector &V, ElementTransformation &T,
996  const IntegrationPoint &ip);
997
999  /** @brief A VectorDeltaFunction cannot be evaluated. Calling this method
1000  will cause an MFEM error, terminating the application. */
1001  virtual void Eval(Vector &V, ElementTransformation &T,
1002  const IntegrationPoint &ip)
1003  { mfem_error("VectorDeltaCoefficient::Eval"); }
1005 };
1006
1007 /** @brief Derived vector coefficient that has the value of the parent vector
1008  where it is active and is zero otherwise. */
1010 {
1011 private:
1012  VectorCoefficient *c;
1013  Array<int> active_attr;
1014
1015 public:
1016  /** @brief Construct with a parent vector coefficient and an array of zeros
1017  and ones representing the attributes for which this coefficient should be
1018  active. */
1020  : VectorCoefficient(vc.GetVDim())
1021  { c = &vc; attr.Copy(active_attr); }
1022
1023  /// Set the time for internally stored coefficients
1024  void SetTime(double t);
1025
1026  /// Evaluate the vector coefficient at @a ip.
1027  virtual void Eval(Vector &V, ElementTransformation &T,
1028  const IntegrationPoint &ip);
1029
1030  /** @brief Evaluate the vector coefficient at all of the locations in the
1031  integration rule and write the vectors into the columns of matrix @a
1032  M. */
1033  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1034  const IntegrationRule &ir);
1035 };
1036
1038
1039 /// Base class for Matrix Coefficients that optionally depend on time and space.
1041 {
1042 protected:
1044  double time;
1045  bool symmetric; // deprecated
1046
1047 public:
1048  /// Construct a dim x dim matrix coefficient.
1049  explicit MatrixCoefficient(int dim, bool symm=false)
1050  { height = width = dim; time = 0.; symmetric = symm; }
1051
1052  /// Construct a h x w matrix coefficient.
1053  MatrixCoefficient(int h, int w, bool symm=false) :
1054  height(h), width(w), time(0.), symmetric(symm) { }
1055
1056  /// Set the time for time dependent coefficients
1057  virtual void SetTime(double t) { time = t; }
1058
1059  /// Get the time for time dependent coefficients
1060  double GetTime() { return time; }
1061
1062  /// Get the height of the matrix.
1063  int GetHeight() const { return height; }
1064
1065  /// Get the width of the matrix.
1066  int GetWidth() const { return width; }
1067
1068  /// For backward compatibility get the width of the matrix.
1069  int GetVDim() const { return width; }
1070
1071  /** @deprecated Use SymmetricMatrixCoefficient instead */
1072  bool IsSymmetric() const { return symmetric; }
1073
1074  /** @brief Evaluate the matrix coefficient in the element described by @a T
1075  at the point @a ip, storing the result in @a K. */
1076  /** @note When this method is called, the caller must make sure that the
1077  IntegrationPoint associated with @a T is the same as @a ip. This can be
1078  achieved by calling T.SetIntPoint(&ip). */
1079  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1080  const IntegrationPoint &ip) = 0;
1081
1082  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
1083  /// the quadrature points. The matrix will be transposed or not according to
1084  /// the boolean argument @a transpose.
1085  ///
1086  /// The @a vdim of the QuadratureFunction should be equal to the height times
1087  /// the width of the matrix.
1088  virtual void Project(QuadratureFunction &qf, bool transpose=false);
1089
1090  /// (DEPRECATED) Evaluate a symmetric matrix coefficient.
1091  /** @brief Evaluate the upper triangular entries of the matrix coefficient
1092  in the symmetric case, similarly to Eval. Matrix entry (i,j) is stored
1093  in K[j - i + os_i] for 0 <= i <= j < width, os_0 = 0,
1094  os_{i+1} = os_i + width - i. That is, K = {M(0,0), ..., M(0,w-1),
1095  M(1,1), ..., M(1,w-1), ..., M(w-1,w-1) with w = width.
1096  @deprecated Use Eval() instead. */
1098  const IntegrationPoint &ip)
1099  { mfem_error("MatrixCoefficient::EvalSymmetric"); }
1100
1101  virtual ~MatrixCoefficient() { }
1102 };
1103
1104
1105 /// A matrix coefficient that is constant in space and time.
1107 {
1108 private:
1109  DenseMatrix mat;
1110 public:
1111  ///Construct using matrix @a m for the constant.
1113  : MatrixCoefficient(m.Height(), m.Width()), mat(m) { }
1115  /// Evaluate the matrix coefficient at @a ip.
1117  const IntegrationPoint &ip) { M = mat; }
1118  /// Return a reference to the constant matrix.
1119  const DenseMatrix& GetMatrix() { return mat; }
1120 };
1121
1122
1123 /** @brief A piecewise matrix-valued coefficient with the pieces keyed off the
1124  element attribute numbers.
1125
1126  A value of zero will be returned for any missing attribute numbers.
1127
1128  This object will not assume ownership of any MatrixCoefficient
1129  objects passed to it. Consequently, the caller must ensure that
1130  the individual MatrixCoefficient objects are not deleted while
1131  this PWMatrixCoefficient is still in use.
1132
1133  \note The keys may either be domain attribute numbers or boundary
1134  attribute numbers. If the PWMatrixCoefficient is used with a
1135  domain integrator the keys are assumed to be domain attribute
1136  numbers. Similarly, if the PWMatrixCoefficient is used with a
1137  boundary integrator the keys are assumed to be boundary attribute
1138  numbers.
1139 */
1141 {
1142 private:
1143  /** Internal data structure to store pointers to the appropriate
1144  coefficients for different regions of the mesh. The keys used
1145  in the map are the mesh attribute numbers (either element
1146  attribute or boundary element attribute depending upon
1147  context). The values returned for any missing attributes will
1148  be zero. The coefficient pointers may be NULL in which case a
1149  value of zero is returned.
1150
1151  The MatrixCoefficient objects contained in this map are NOT
1152  owned by this PWMatrixCoefficient object. This means that they
1153  will not be deleted when this object is deleted also the caller
1154  must ensure that the various MatrixCoefficient objects are not
1155  deleted while this PWMatrixCoefficient is still needed.
1156  */
1157  std::map<int, MatrixCoefficient*> pieces;
1158
1159  /** Convenience function to check for compatible array lengths,
1160  loop over the arrays, and add their attribute/MatrixCoefficient
1161  pairs to the internal data structure.
1162  */
1163  void InitMap(const Array<int> & attr,
1164  const Array<MatrixCoefficient*> & coefs);
1165
1166 public:
1167
1168  /// Constructs a piecewise matrix coefficient of dimension dim by dim
1169  explicit PWMatrixCoefficient(int dim, bool symm = false)
1170  : MatrixCoefficient(dim, symm) {}
1171
1172  /// Constructs a piecewise matrix coefficient of dimension h by w
1173  explicit PWMatrixCoefficient(int h, int w, bool symm = false)
1174  : MatrixCoefficient(h, w, symm) {}
1175
1176  /// Construct the coefficient using arrays describing the pieces
1177  /** \param dim - size of the square matrix-valued result
1178  \param attr - an array of attribute numbers for each piece
1179  \param coefs - the corresponding array of MatrixCoefficient pointers
1180  \param symm - true if the result will be symmetric, false otherwise
1181  Any missing attributes or NULL coefficient pointers will result in a
1182  zero matrix being returned.
1183
1184  \note Ownership of the MatrixCoefficient objects will NOT be
1185  transferred to this object.
1186  */
1188  const Array<MatrixCoefficient*> & coefs,
1189  bool symm=false)
1190  : MatrixCoefficient(dim, symm) { InitMap(attr, coefs); }
1191
1192  /// Construct the coefficient using arrays describing the pieces
1193  /** \param h - height of the matrix-valued result
1194  \param w - width of the matrix-valued result
1195  \param attr - an array of attribute numbers for each piece
1196  \param coefs - the corresponding array of MatrixCoefficient pointers
1197  \param symm - true if the result will be symmetric, false otherwise
1198  Any missing attributes or NULL coefficient pointers will result in a
1199  zero matrix being returned for that attribute.
1200
1201  \note Ownership of the MatrixCoefficient objects will NOT be
1202  transferred to this object.
1203  */
1204  PWMatrixCoefficient(int h, int w, const Array<int> & attr,
1205  const Array<MatrixCoefficient*> & coefs,
1206  bool symm=false)
1207  : MatrixCoefficient(h, w, symm) { InitMap(attr, coefs); }
1208
1209  /// Set the time for time dependent coefficients
1210  virtual void SetTime(double t);
1211
1212  /// Replace a set of coefficients
1213  void UpdateCoefficients(const Array<int> & attr,
1214  const Array<MatrixCoefficient*> & coefs)
1215  { InitMap(attr, coefs); }
1216
1217  /// Replace a single coefficient for a particular attribute
1218  void UpdateCoefficient(int attr, MatrixCoefficient & coef);
1219
1220  /// Remove a single MatrixCoefficient for a particular attribute
1221  void ZeroCoefficient(int attr)
1222  { pieces.erase(attr); }
1223
1224  /// Evaluate the coefficient.
1225  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1226  const IntegrationPoint &ip);
1227 };
1228
1229 /** @brief A matrix coefficient with an optional scalar coefficient multiplier
1230  \a q. The matrix function can either be represented by a std function or
1231  a constant matrix provided when constructing this object. */
1233 {
1234 private:
1235  std::function<void(const Vector &, DenseMatrix &)> Function;
1236  std::function<void(const Vector &, Vector &)> SymmFunction; // deprecated
1237  std::function<void(const Vector &, double, DenseMatrix &)> TDFunction;
1238
1239  Coefficient *Q;
1240  DenseMatrix mat;
1241
1242 public:
1243  /// Define a time-independent square matrix coefficient from a std function
1244  /** \param dim - the size of the matrix
1245  \param F - time-independent function
1246  \param q - optional scalar Coefficient to scale the matrix coefficient */
1248  std::function<void(const Vector &, DenseMatrix &)> F,
1249  Coefficient *q = nullptr)
1250  : MatrixCoefficient(dim), Function(std::move(F)), Q(q), mat(0)
1251  { }
1252
1253  /// Define a constant matrix coefficient times a scalar Coefficient
1254  /** \param m - constant matrix
1255  \param q - optional scalar Coefficient to scale the matrix coefficient */
1257  : MatrixCoefficient(m.Height(), m.Width()), Q(&q), mat(m)
1258  { }
1259
1260  /** @brief Define a time-independent symmetric square matrix coefficient from
1261  a std function */
1262  /** \param dim - the size of the matrix
1263  \param SymmF - function used in EvalSymmetric
1264  \param q - optional scalar Coefficient to scale the matrix coefficient
1265  @deprecated Use another constructor without setting SymmFunction. */
1267  std::function<void(const Vector &, Vector &)> SymmF,
1268  Coefficient *q = NULL)
1269  : MatrixCoefficient(dim, true), SymmFunction(std::move(SymmF)), Q(q), mat(0)
1270  { }
1271
1272  /// Define a time-dependent square matrix coefficient from a std function
1273  /** \param dim - the size of the matrix
1274  \param TDF - time-dependent function
1275  \param q - optional scalar Coefficient to scale the matrix coefficient */
1277  std::function<void(const Vector &, double, DenseMatrix &)> TDF,
1278  Coefficient *q = nullptr)
1279  : MatrixCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
1280  { }
1281
1282  /// Set the time for internally stored coefficients
1283  void SetTime(double t);
1284
1285  /// Evaluate the matrix coefficient at @a ip.
1286  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1287  const IntegrationPoint &ip);
1288
1289  /// (DEPRECATED) Evaluate the symmetric matrix coefficient at @a ip.
1290  /** @deprecated Use Eval() instead. */
1291  virtual void EvalSymmetric(Vector &K, ElementTransformation &T,
1292  const IntegrationPoint &ip);
1293
1295 };
1296
1297
1298 /** @brief Matrix coefficient defined by a matrix of scalar coefficients.
1299  Coefficients that are not set will evaluate to zero in the vector. The
1300  coefficient is stored as a flat Array with indexing (i,j) -> i*width+j. */
1302 {
1303 private:
1304  Array<Coefficient *> Coeff;
1305  Array<bool> ownCoeff;
1306
1307 public:
1308  /** @brief Construct a coefficient matrix of dimensions @a dim * @a dim. The
1309  actual coefficients still need to be added with Set(). */
1310  explicit MatrixArrayCoefficient (int dim);
1311
1312  /// Set the time for internally stored coefficients
1313  void SetTime(double t);
1314
1315  /// Get the coefficient located at (i,j) in the matrix.
1316  Coefficient* GetCoeff (int i, int j) { return Coeff[i*width+j]; }
1317
1318  /** @brief Set the coefficient located at (i,j) in the matrix. By default by
1319  default this will take ownership of the Coefficient passed in, but this
1320  can be overridden with the @a own parameter. */
1321  void Set(int i, int j, Coefficient * c, bool own=true);
1322
1324
1325  /// Evaluate coefficient located at (i,j) in the matrix using integration
1326  /// point @a ip.
1327  double Eval(int i, int j, ElementTransformation &T, const IntegrationPoint &ip)
1328  { return Coeff[i*width+j] ? Coeff[i*width+j] -> Eval(T, ip, GetTime()) : 0.0; }
1329
1330  /// Evaluate the matrix coefficient @a ip.
1331  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1332  const IntegrationPoint &ip);
1333
1334  virtual ~MatrixArrayCoefficient();
1335 };
1336
1337
1338 /** @brief Derived matrix coefficient that has the value of the parent matrix
1339  coefficient where it is active and is zero otherwise. */
1341 {
1342 private:
1343  MatrixCoefficient *c;
1344  Array<int> active_attr;
1345
1346 public:
1347  /** @brief Construct with a parent matrix coefficient and an array of zeros
1348  and ones representing the attributes for which this coefficient should be
1349  active. */
1351  : MatrixCoefficient(mc.GetHeight(), mc.GetWidth())
1352  { c = &mc; attr.Copy(active_attr); }
1353
1354  /// Set the time for internally stored coefficients
1355  void SetTime(double t);
1356
1357  /// Evaluate the matrix coefficient at @a ip.
1358  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1359  const IntegrationPoint &ip);
1360 };
1361
1362 /// Coefficients based on sums, products, or other functions of coefficients.
1363 ///@{
1364 /** @brief Scalar coefficient defined as the linear combination of two scalar
1365  coefficients or a scalar and a scalar coefficient */
1367 {
1368 private:
1369  double aConst;
1370  Coefficient * a;
1371  Coefficient * b;
1372
1373  double alpha;
1374  double beta;
1375
1376 public:
1377  /// Constructor with one coefficient. Result is alpha_ * A + beta_ * B
1379  double alpha_ = 1.0, double beta_ = 1.0)
1380  : aConst(A), a(NULL), b(&B), alpha(alpha_), beta(beta_) { }
1381
1382  /// Constructor with two coefficients. Result is alpha_ * A + beta_ * B.
1384  double alpha_ = 1.0, double beta_ = 1.0)
1385  : aConst(0.0), a(&A), b(&B), alpha(alpha_), beta(beta_) { }
1386
1387  /// Set the time for internally stored coefficients
1388  void SetTime(double t);
1389
1390  /// Reset the first term in the linear combination as a constant
1391  void SetAConst(double A) { a = NULL; aConst = A; }
1392  /// Return the first term in the linear combination
1393  double GetAConst() const { return aConst; }
1394
1395  /// Reset the first term in the linear combination
1396  void SetACoef(Coefficient &A) { a = &A; }
1397  /// Return the first term in the linear combination
1398  Coefficient * GetACoef() const { return a; }
1399
1400  /// Reset the second term in the linear combination
1401  void SetBCoef(Coefficient &B) { b = &B; }
1402  /// Return the second term in the linear combination
1403  Coefficient * GetBCoef() const { return b; }
1404
1405  /// Reset the factor in front of the first term in the linear combination
1406  void SetAlpha(double alpha_) { alpha = alpha_; }
1407  /// Return the factor in front of the first term in the linear combination
1408  double GetAlpha() const { return alpha; }
1409
1410  /// Reset the factor in front of the second term in the linear combination
1411  void SetBeta(double beta_) { beta = beta_; }
1412  /// Return the factor in front of the second term in the linear combination
1413  double GetBeta() const { return beta; }
1414
1415  /// Evaluate the coefficient at @a ip.
1416  virtual double Eval(ElementTransformation &T,
1417  const IntegrationPoint &ip)
1418  {
1419  return alpha * ((a == NULL ) ? aConst : a->Eval(T, ip) )
1420  + beta * b->Eval(T, ip);
1421  }
1422 };
1423
1424
1425 /// Base class for symmetric matrix coefficients that optionally depend on time and space.
1427 {
1428 protected:
1429  /// Internal matrix used when evaluating this coefficient as a DenseMatrix.
1431 public:
1432  /// Construct a dim x dim matrix coefficient.
1434  : MatrixCoefficient(dimension, true) { }
1435
1436  /// Get the size of the matrix.
1437  int GetSize() const { return height; }
1438
1439  /// @brief Fill the QuadratureFunction @a qf by evaluating the coefficient at
1441  ///
1442  /// @note As opposed to MatrixCoefficient::Project, this function stores only
1443  /// the @a symmetric part of the matrix at each quadrature point.
1444  ///
1445  /// The @a vdim of the coefficient should be equal to height*(height+1)/2.
1447
1448  /** @brief Evaluate the matrix coefficient in the element described by @a T
1449  at the point @a ip, storing the result as a symmetric matrix @a K. */
1450  /** @note When this method is called, the caller must make sure that the
1451  IntegrationPoint associated with @a T is the same as @a ip. This can be
1452  achieved by calling T.SetIntPoint(&ip). */
1453  virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T,
1454  const IntegrationPoint &ip) = 0;
1455
1456  /** @brief Evaluate the matrix coefficient in the element described by @a T
1457  at the point @a ip, storing the result as a dense matrix @a K. */
1458  /** This function allows the use of SymmetricMatrixCoefficient in situations
1459  where the symmetry is not taken advantage of.
1460
1461  @note When this method is called, the caller must make sure that the
1462  IntegrationPoint associated with @a T is the same as @a ip. This can be
1463  achieved by calling T.SetIntPoint(&ip). */
1464  virtual void Eval(DenseMatrix &K, ElementTransformation &T,
1465  const IntegrationPoint &ip);
1466
1467  /// Return a reference to the constant matrix.
1468  const DenseSymmetricMatrix& GetMatrix() { return mat; }
1469
1471 };
1472
1473
1474 /// A matrix coefficient that is constant in space and time.
1476 {
1477 private:
1479
1480 public:
1481  ///Construct using matrix @a m for the constant.
1483  : SymmetricMatrixCoefficient(m.Height()), mat(m) { }
1485  /// Evaluate the matrix coefficient at @a ip.
1487  const IntegrationPoint &ip) { M = mat; }
1488 };
1489
1490
1491 /** @brief A matrix coefficient with an optional scalar coefficient multiplier
1492  \a q. The matrix function can either be represented by a std function or
1493  a constant matrix provided when constructing this object. */
1495 {
1496 private:
1497  std::function<void(const Vector &, DenseSymmetricMatrix &)> Function;
1498  std::function<void(const Vector &, double, DenseSymmetricMatrix &)> TDFunction;
1499
1500  Coefficient *Q;
1502
1503 public:
1504  /// Define a time-independent symmetric matrix coefficient from a std function
1505  /** \param dim - the size of the matrix
1506  \param F - time-independent function
1507  \param q - optional scalar Coefficient to scale the matrix coefficient */
1509  std::function<void(const Vector &, DenseSymmetricMatrix &)> F,
1510  Coefficient *q = nullptr)
1511  : SymmetricMatrixCoefficient(dim), Function(std::move(F)), Q(q), mat(0)
1512  { }
1513
1514  /// Define a constant matrix coefficient times a scalar Coefficient
1515  /** \param m - constant matrix
1516  \param q - optional scalar Coefficient to scale the matrix coefficient */
1518  Coefficient &q)
1519  : SymmetricMatrixCoefficient(m.Height()), Q(&q), mat(m)
1520  { }
1521
1522  /// Define a time-dependent square matrix coefficient from a std function
1523  /** \param dim - the size of the matrix
1524  \param TDF - time-dependent function
1525  \param q - optional scalar Coefficient to scale the matrix coefficient */
1527  std::function<void(const Vector &, double, DenseSymmetricMatrix &)> TDF,
1528  Coefficient *q = nullptr)
1529  : SymmetricMatrixCoefficient(dim), TDFunction(std::move(TDF)), Q(q)
1530  { }
1531
1532  /// Set the time for internally stored coefficients
1533  void SetTime(double t);
1534
1536  /// Evaluate the matrix coefficient at @a ip.
1537  virtual void Eval(DenseSymmetricMatrix &K, ElementTransformation &T,
1538  const IntegrationPoint &ip);
1539
1541 };
1542
1543
1544 /** @brief Scalar coefficient defined as the product of two scalar coefficients
1545  or a scalar and a scalar coefficient. */
1547 {
1548 private:
1549  double aConst;
1550  Coefficient * a;
1551  Coefficient * b;
1552
1553 public:
1554  /// Constructor with one coefficient. Result is A * B.
1556  : aConst(A), a(NULL), b(&B) { }
1557
1558  /// Constructor with two coefficients. Result is A * B.
1560  : aConst(0.0), a(&A), b(&B) { }
1561
1562  /// Set the time for internally stored coefficients
1563  void SetTime(double t);
1564
1565  /// Reset the first term in the product as a constant
1566  void SetAConst(double A) { a = NULL; aConst = A; }
1567  /// Return the first term in the product
1568  double GetAConst() const { return aConst; }
1569
1570  /// Reset the first term in the product
1571  void SetACoef(Coefficient &A) { a = &A; }
1572  /// Return the first term in the product
1573  Coefficient * GetACoef() const { return a; }
1574
1575  /// Reset the second term in the product
1576  void SetBCoef(Coefficient &B) { b = &B; }
1577  /// Return the second term in the product
1578  Coefficient * GetBCoef() const { return b; }
1579
1580  /// Evaluate the coefficient at @a ip.
1581  virtual double Eval(ElementTransformation &T,
1582  const IntegrationPoint &ip)
1583  { return ((a == NULL ) ? aConst : a->Eval(T, ip) ) * b->Eval(T, ip); }
1584 };
1585
1586 /** @brief Scalar coefficient defined as the ratio of two scalars where one or
1587  both scalars are scalar coefficients. */
1589 {
1590 private:
1591  double aConst;
1592  double bConst;
1593  Coefficient * a;
1594  Coefficient * b;
1595
1596 public:
1597  /** Initialize a coefficient which returns A / B where @a A is a
1598  constant and @a B is a scalar coefficient */
1600  : aConst(A), bConst(1.0), a(NULL), b(&B) { }
1601  /** Initialize a coefficient which returns A / B where @a A and @a B are both
1602  scalar coefficients */
1604  : aConst(0.0), bConst(1.0), a(&A), b(&B) { }
1605  /** Initialize a coefficient which returns A / B where @a A is a
1606  scalar coefficient and @a B is a constant */
1608  : aConst(0.0), bConst(B), a(&A), b(NULL) { }
1609
1610  /// Set the time for internally stored coefficients
1611  void SetTime(double t);
1612
1613  /// Reset the numerator in the ratio as a constant
1614  void SetAConst(double A) { a = NULL; aConst = A; }
1615  /// Return the numerator of the ratio
1616  double GetAConst() const { return aConst; }
1617
1618  /// Reset the denominator in the ratio as a constant
1619  void SetBConst(double B) { b = NULL; bConst = B; }
1620  /// Return the denominator of the ratio
1621  double GetBConst() const { return bConst; }
1622
1623  /// Reset the numerator in the ratio
1624  void SetACoef(Coefficient &A) { a = &A; }
1625  /// Return the numerator of the ratio
1626  Coefficient * GetACoef() const { return a; }
1627
1628  /// Reset the denominator in the ratio
1629  void SetBCoef(Coefficient &B) { b = &B; }
1630  /// Return the denominator of the ratio
1631  Coefficient * GetBCoef() const { return b; }
1632
1633  /// Evaluate the coefficient
1634  virtual double Eval(ElementTransformation &T,
1635  const IntegrationPoint &ip)
1636  {
1637  double den = (b == NULL ) ? bConst : b->Eval(T, ip);
1638  MFEM_ASSERT(den != 0.0, "Division by zero in RatioCoefficient");
1639  return ((a == NULL ) ? aConst : a->Eval(T, ip) ) / den;
1640  }
1641 };
1642
1643 /// Scalar coefficient defined as a scalar raised to a power
1645 {
1646 private:
1647  Coefficient * a;
1648
1649  double p;
1650
1651 public:
1652  /// Construct with a coefficient and a constant power @a p_. Result is A^p.
1654  : a(&A), p(p_) { }
1655
1656  /// Set the time for internally stored coefficients
1657  void SetTime(double t);
1658
1659  /// Reset the base coefficient
1660  void SetACoef(Coefficient &A) { a = &A; }
1661  /// Return the base coefficient
1662  Coefficient * GetACoef() const { return a; }
1663
1664  /// Reset the exponent
1665  void SetExponent(double p_) { p = p_; }
1666  /// Return the exponent
1667  double GetExponent() const { return p; }
1668
1669  /// Evaluate the coefficient at @a ip.
1670  virtual double Eval(ElementTransformation &T,
1671  const IntegrationPoint &ip)
1672  { return pow(a->Eval(T, ip), p); }
1673 };
1674
1675
1676 /// Scalar coefficient defined as the inner product of two vector coefficients
1678 {
1679 private:
1680  VectorCoefficient * a;
1681  VectorCoefficient * b;
1682
1683  mutable Vector va;
1684  mutable Vector vb;
1685 public:
1686  /// Construct with the two vector coefficients. Result is \f$A \cdot B \f$.
1688
1689  /// Set the time for internally stored coefficients
1690  void SetTime(double t);
1691
1692  /// Reset the first vector in the inner product
1693  void SetACoef(VectorCoefficient &A) { a = &A; }
1694  /// Return the first vector coefficient in the inner product
1695  VectorCoefficient * GetACoef() const { return a; }
1696
1697  /// Reset the second vector in the inner product
1698  void SetBCoef(VectorCoefficient &B) { b = &B; }
1699  /// Return the second vector coefficient in the inner product
1700  VectorCoefficient * GetBCoef() const { return b; }
1701
1702  /// Evaluate the coefficient at @a ip.
1703  virtual double Eval(ElementTransformation &T,
1704  const IntegrationPoint &ip);
1705 };
1706
1707 /// Scalar coefficient defined as a cross product of two vectors in the xy-plane.
1709 {
1710 private:
1711  VectorCoefficient * a;
1712  VectorCoefficient * b;
1713
1714  mutable Vector va;
1715  mutable Vector vb;
1716
1717 public:
1718  /// Constructor with two vector coefficients. Result is \f$A_x B_y - A_y * B_x; \f$.
1720
1721  /// Set the time for internally stored coefficients
1722  void SetTime(double t);
1723
1724  /// Reset the first vector in the product
1725  void SetACoef(VectorCoefficient &A) { a = &A; }
1726  /// Return the first vector of the product
1727  VectorCoefficient * GetACoef() const { return a; }
1728
1729  /// Reset the second vector in the product
1730  void SetBCoef(VectorCoefficient &B) { b = &B; }
1731  /// Return the second vector of the product
1732  VectorCoefficient * GetBCoef() const { return b; }
1733
1734  /// Evaluate the coefficient at @a ip.
1735  virtual double Eval(ElementTransformation &T,
1736  const IntegrationPoint &ip);
1737 };
1738
1739 /// Scalar coefficient defined as the determinant of a matrix coefficient
1741 {
1742 private:
1743  MatrixCoefficient * a;
1744
1745  mutable DenseMatrix ma;
1746
1747 public:
1748  /// Construct with the matrix.
1750
1751  /// Set the time for internally stored coefficients
1752  void SetTime(double t);
1753
1754  /// Reset the matrix coefficient
1755  void SetACoef(MatrixCoefficient &A) { a = &A; }
1756  /// Return the matrix coefficient
1757  MatrixCoefficient * GetACoef() const { return a; }
1758
1759  /// Evaluate the determinant coefficient at @a ip.
1760  virtual double Eval(ElementTransformation &T,
1761  const IntegrationPoint &ip);
1762 };
1763
1764 /// Vector coefficient defined as the linear combination of two vectors
1766 {
1767 private:
1768  VectorCoefficient * ACoef;
1769  VectorCoefficient * BCoef;
1770
1771  Vector A;
1772  Vector B;
1773
1774  Coefficient * alphaCoef;
1775  Coefficient * betaCoef;
1776
1777  double alpha;
1778  double beta;
1779
1780  mutable Vector va;
1781
1782 public:
1783  /** Constructor with no coefficients.
1784  To be used with the various "Set" methods */
1786
1787  /** Constructor with two vector coefficients.
1788  Result is alpha_ * A + beta_ * B */
1790  double alpha_ = 1.0, double beta_ = 1.0);
1791
1792  /** Constructor with scalar coefficients.
1793  Result is alpha_ * A_ + beta_ * B_ */
1795  Coefficient &alpha_, Coefficient &beta_);
1796
1797  /// Set the time for internally stored coefficients
1798  void SetTime(double t);
1799
1800  /// Reset the first vector coefficient
1801  void SetACoef(VectorCoefficient &A_) { ACoef = &A_; }
1802  /// Return the first vector coefficient
1803  VectorCoefficient * GetACoef() const { return ACoef; }
1804
1805  /// Reset the second vector coefficient
1806  void SetBCoef(VectorCoefficient &B_) { BCoef = &B_; }
1807  /// Return the second vector coefficient
1808  VectorCoefficient * GetBCoef() const { return BCoef; }
1809
1810  /// Reset the factor in front of the first vector coefficient
1811  void SetAlphaCoef(Coefficient &A_) { alphaCoef = &A_; }
1812  /// Return the factor in front of the first vector coefficient
1813  Coefficient * GetAlphaCoef() const { return alphaCoef; }
1814
1815  /// Reset the factor in front of the second vector coefficient
1816  void SetBetaCoef(Coefficient &B_) { betaCoef = &B_; }
1817  /// Return the factor in front of the second vector coefficient
1818  Coefficient * GetBetaCoef() const { return betaCoef; }
1819
1820  /// Reset the first vector as a constant
1821  void SetA(const Vector &A_) { A = A_; ACoef = NULL; }
1822  /// Return the first vector constant
1823  const Vector & GetA() const { return A; }
1824
1825  /// Reset the second vector as a constant
1826  void SetB(const Vector &B_) { B = B_; BCoef = NULL; }
1827  /// Return the second vector constant
1828  const Vector & GetB() const { return B; }
1829
1830  /// Reset the factor in front of the first vector coefficient as a constant
1831  void SetAlpha(double alpha_) { alpha = alpha_; alphaCoef = NULL; }
1832  /// Return the factor in front of the first vector coefficient
1833  double GetAlpha() const { return alpha; }
1834
1835  /// Reset the factor in front of the second vector coefficient as a constant
1836  void SetBeta(double beta_) { beta = beta_; betaCoef = NULL; }
1837  /// Return the factor in front of the second vector coefficient
1838  double GetBeta() const { return beta; }
1839
1840  /// Evaluate the coefficient at @a ip.
1841  virtual void Eval(Vector &V, ElementTransformation &T,
1842  const IntegrationPoint &ip);
1844 };
1845
1846 /// Vector coefficient defined as a product of scalar and vector coefficients.
1848 {
1849 private:
1850  double aConst;
1851  Coefficient * a;
1852  VectorCoefficient * b;
1853
1854 public:
1855  /// Constructor with constant and vector coefficient. Result is A * B.
1857
1858  /// Constructor with two coefficients. Result is A * B.
1860
1861  /// Set the time for internally stored coefficients
1862  void SetTime(double t);
1863
1864  /// Reset the scalar factor as a constant
1865  void SetAConst(double A) { a = NULL; aConst = A; }
1866  /// Return the scalar factor
1867  double GetAConst() const { return aConst; }
1868
1869  /// Reset the scalar factor
1870  void SetACoef(Coefficient &A) { a = &A; }
1871  /// Return the scalar factor
1872  Coefficient * GetACoef() const { return a; }
1873
1874  /// Reset the vector factor
1875  void SetBCoef(VectorCoefficient &B) { b = &B; }
1876  /// Return the vector factor
1877  VectorCoefficient * GetBCoef() const { return b; }
1878
1879  /// Evaluate the coefficient at @a ip.
1880  virtual void Eval(Vector &V, ElementTransformation &T,
1881  const IntegrationPoint &ip);
1883 };
1884
1885 /// Vector coefficient defined as a normalized vector field (returns v/|v|)
1887 {
1888 private:
1889  VectorCoefficient * a;
1890
1891  double tol;
1892
1893 public:
1894  /** @brief Return a vector normalized to a length of one
1895
1896  This class evaluates the vector coefficient @a A and, if |A| > @a tol,
1897  returns the normalized vector A / |A|. If |A| <= @a tol, the zero
1898  vector is returned.
1899  */
1900  NormalizedVectorCoefficient(VectorCoefficient &A, double tol = 1e-6);
1901
1902  /// Set the time for internally stored coefficients
1903  void SetTime(double t);
1904
1905  /// Reset the vector coefficient
1906  void SetACoef(VectorCoefficient &A) { a = &A; }
1907  /// Return the vector coefficient
1908  VectorCoefficient * GetACoef() const { return a; }
1909
1910  /// Evaluate the coefficient at @a ip.
1911  virtual void Eval(Vector &V, ElementTransformation &T,
1912  const IntegrationPoint &ip);
1914 };
1915
1916 /// Vector coefficient defined as a cross product of two vectors
1918 {
1919 private:
1920  VectorCoefficient * a;
1921  VectorCoefficient * b;
1922
1923  mutable Vector va;
1924  mutable Vector vb;
1925
1926 public:
1927  /// Construct with the two coefficients. Result is A x B.
1929
1930  /// Set the time for internally stored coefficients
1931  void SetTime(double t);
1932
1933  /// Reset the first term in the product
1934  void SetACoef(VectorCoefficient &A) { a = &A; }
1935  /// Return the first term in the product
1936  VectorCoefficient * GetACoef() const { return a; }
1937
1938  /// Reset the second term in the product
1939  void SetBCoef(VectorCoefficient &B) { b = &B; }
1940  /// Return the second term in the product
1941  VectorCoefficient * GetBCoef() const { return b; }
1942
1943  /// Evaluate the coefficient at @a ip.
1944  virtual void Eval(Vector &V, ElementTransformation &T,
1945  const IntegrationPoint &ip);
1947 };
1948
1949 /** @brief Vector coefficient defined as a product of a matrix coefficient and
1950  a vector coefficient. */
1952 {
1953 private:
1954  MatrixCoefficient * a;
1955  VectorCoefficient * b;
1956
1957  mutable DenseMatrix ma;
1958  mutable Vector vb;
1959
1960 public:
1961  /// Constructor with two coefficients. Result is A*B.
1963
1964  /// Set the time for internally stored coefficients
1965  void SetTime(double t);
1966
1967  /// Reset the matrix coefficient
1968  void SetACoef(MatrixCoefficient &A) { a = &A; }
1969  /// Return the matrix coefficient
1970  MatrixCoefficient * GetACoef() const { return a; }
1971
1972  /// Reset the vector coefficient
1973  void SetBCoef(VectorCoefficient &B) { b = &B; }
1974  /// Return the vector coefficient
1975  VectorCoefficient * GetBCoef() const { return b; }
1976
1977  /// Evaluate the vector coefficient at @a ip.
1978  virtual void Eval(Vector &V, ElementTransformation &T,
1979  const IntegrationPoint &ip);
1981 };
1982
1983 /// Convenient alias for the MatrixVectorProductCoefficient
1985
1986 /// Constant matrix coefficient defined as the identity of dimension d
1988 {
1989 private:
1990  int dim;
1991
1992 public:
1993  /// Construct with the dimension of the square identity matrix.
1995  : MatrixCoefficient(d, d), dim(d) { }
1996
1997  /// Evaluate the matrix coefficient at @a ip.
1998  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
1999  const IntegrationPoint &ip);
2000 };
2001
2002 /// Matrix coefficient defined as the linear combination of two matrices
2004 {
2005 private:
2006  MatrixCoefficient * a;
2007  MatrixCoefficient * b;
2008
2009  double alpha;
2010  double beta;
2011
2012  mutable DenseMatrix ma;
2013
2014 public:
2015  /// Construct with the two coefficients. Result is alpha_ * A + beta_ * B.
2017  double alpha_ = 1.0, double beta_ = 1.0);
2018
2019  /// Set the time for internally stored coefficients
2020  void SetTime(double t);
2021
2022  /// Reset the first matrix coefficient
2023  void SetACoef(MatrixCoefficient &A) { a = &A; }
2024  /// Return the first matrix coefficient
2025  MatrixCoefficient * GetACoef() const { return a; }
2026
2027  /// Reset the second matrix coefficient
2028  void SetBCoef(MatrixCoefficient &B) { b = &B; }
2029  /// Return the second matrix coefficient
2030  MatrixCoefficient * GetBCoef() const { return b; }
2031
2032  /// Reset the factor in front of the first matrix coefficient
2033  void SetAlpha(double alpha_) { alpha = alpha_; }
2034  /// Return the factor in front of the first matrix coefficient
2035  double GetAlpha() const { return alpha; }
2036
2037  /// Reset the factor in front of the second matrix coefficient
2038  void SetBeta(double beta_) { beta = beta_; }
2039  /// Return the factor in front of the second matrix coefficient
2040  double GetBeta() const { return beta; }
2041
2042  /// Evaluate the matrix coefficient at @a ip.
2043  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2044  const IntegrationPoint &ip);
2045 };
2046
2047 /// Matrix coefficient defined as the product of two matrices
2049 {
2050 private:
2051  MatrixCoefficient * a;
2052  MatrixCoefficient * b;
2053
2054  mutable DenseMatrix ma;
2055  mutable DenseMatrix mb;
2056
2057 public:
2058  /// Construct with the two coefficients. Result is A * B.
2060
2061  /// Reset the first matrix coefficient
2062  void SetACoef(MatrixCoefficient &A) { a = &A; }
2063  /// Return the first matrix coefficient
2064  MatrixCoefficient * GetACoef() const { return a; }
2065
2066  /// Reset the second matrix coefficient
2067  void SetBCoef(MatrixCoefficient &B) { b = &B; }
2068  /// Return the second matrix coefficient
2069  MatrixCoefficient * GetBCoef() const { return b; }
2070
2071  /// Evaluate the matrix coefficient at @a ip.
2072  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2073  const IntegrationPoint &ip);
2074 };
2075
2076 /** @brief Matrix coefficient defined as a product of a scalar coefficient and a
2077  matrix coefficient.*/
2079 {
2080 private:
2081  double aConst;
2082  Coefficient * a;
2083  MatrixCoefficient * b;
2084
2085 public:
2086  /// Constructor with one coefficient. Result is A*B.
2088
2089  /// Constructor with two coefficients. Result is A*B.
2091
2092  /// Set the time for internally stored coefficients
2093  void SetTime(double t);
2094
2095  /// Reset the scalar factor as a constant
2096  void SetAConst(double A) { a = NULL; aConst = A; }
2097  /// Return the scalar factor
2098  double GetAConst() const { return aConst; }
2099
2100  /// Reset the scalar factor
2101  void SetACoef(Coefficient &A) { a = &A; }
2102  /// Return the scalar factor
2103  Coefficient * GetACoef() const { return a; }
2104
2105  /// Reset the matrix factor
2106  void SetBCoef(MatrixCoefficient &B) { b = &B; }
2107  /// Return the matrix factor
2108  MatrixCoefficient * GetBCoef() const { return b; }
2109
2110  /// Evaluate the matrix coefficient at @a ip.
2111  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2112  const IntegrationPoint &ip);
2113 };
2114
2115 /// Matrix coefficient defined as the transpose a matrix coefficient
2117 {
2118 private:
2119  MatrixCoefficient * a;
2120
2121 public:
2122  /// Construct with the matrix coefficient. Result is \f$A^T \f$.
2124
2125  /// Set the time for internally stored coefficients
2126  void SetTime(double t);
2127
2128  /// Reset the matrix coefficient
2129  void SetACoef(MatrixCoefficient &A) { a = &A; }
2130  /// Return the matrix coefficient
2131  MatrixCoefficient * GetACoef() const { return a; }
2132
2133  /// Evaluate the matrix coefficient at @a ip.
2134  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2135  const IntegrationPoint &ip);
2136 };
2137
2138 /// Matrix coefficient defined as the inverse a matrix coefficient.
2140 {
2141 private:
2142  MatrixCoefficient * a;
2143
2144 public:
2145  /// Construct with the matrix coefficient. Result is \f$A^{-1} \f$.
2147
2148  /// Set the time for internally stored coefficients
2149  void SetTime(double t);
2150
2151  /// Reset the matrix coefficient
2152  void SetACoef(MatrixCoefficient &A) { a = &A; }
2153  /// Return the matrix coefficient
2154  MatrixCoefficient * GetACoef() const { return a; }
2155
2156  /// Evaluate the matrix coefficient at @a ip.
2157  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2158  const IntegrationPoint &ip);
2159 };
2160
2161 /// Matrix coefficient defined as the outer product of two vector coefficients.
2163 {
2164 private:
2165  VectorCoefficient * a;
2166  VectorCoefficient * b;
2167
2168  mutable Vector va;
2169  mutable Vector vb;
2170
2171 public:
2172  /// Construct with two vector coefficients. Result is \f$A B^T \f$.
2174
2175  /// Set the time for internally stored coefficients
2176  void SetTime(double t);
2177
2178  /// Reset the first vector in the outer product
2179  void SetACoef(VectorCoefficient &A) { a = &A; }
2180  /// Return the first vector coefficient in the outer product
2181  VectorCoefficient * GetACoef() const { return a; }
2182
2183  /// Reset the second vector in the outer product
2184  void SetBCoef(VectorCoefficient &B) { b = &B; }
2185  /// Return the second vector coefficient in the outer product
2186  VectorCoefficient * GetBCoef() const { return b; }
2187
2188  /// Evaluate the matrix coefficient at @a ip.
2189  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2190  const IntegrationPoint &ip);
2191 };
2192
2193 /** @brief Matrix coefficient defined as -a k x k x, for a vector k and scalar a
2194
2195  This coefficient returns \f$a * (|k|^2 I - k \otimes k)\f$, where I is
2196  the identity matrix and \f$\otimes\f$ indicates the outer product. This
2197  can be evaluated for vectors of any dimension but in three
2198  dimensions it corresponds to computing the cross product with k twice.
2199 */
2201 {
2202 private:
2203  double aConst;
2204  Coefficient * a;
2205  VectorCoefficient * k;
2206
2207  mutable Vector vk;
2208
2209 public:
2212
2213  /// Set the time for internally stored coefficients
2214  void SetTime(double t);
2215
2216  /// Reset the scalar factor as a constant
2217  void SetAConst(double A) { a = NULL; aConst = A; }
2218  /// Return the scalar factor
2219  double GetAConst() const { return aConst; }
2220
2221  /// Reset the scalar factor
2222  void SetACoef(Coefficient &A) { a = &A; }
2223  /// Return the scalar factor
2224  Coefficient * GetACoef() const { return a; }
2225
2226  /// Reset the vector factor
2227  void SetKCoef(VectorCoefficient &K) { k = &K; }
2228  /// Return the vector factor
2229  VectorCoefficient * GetKCoef() const { return k; }
2230
2231  /// Evaluate the matrix coefficient at @a ip.
2232  virtual void Eval(DenseMatrix &M, ElementTransformation &T,
2233  const IntegrationPoint &ip);
2234 };
2235 ///@}
2236
2237 /** @brief Vector quadrature function coefficient which requires that the
2238  quadrature rules used for this vector coefficient be the same as those that
2239  live within the supplied QuadratureFunction. */
2241 {
2242 private:
2244  int index;
2245
2246 public:
2247  /// Constructor with a quadrature function as input
2249
2250  /** Set the starting index within the QuadFunc that'll be used to project
2251  outwards as well as the corresponding length. The projected length should
2252  have the bounds of 1 <= length <= (length QuadFunc - index). */
2253  void SetComponent(int index_, int length_);
2254
2256
2258  virtual void Eval(Vector &V, ElementTransformation &T,
2259  const IntegrationPoint &ip);
2260
2262
2264 };
2265
2267  rules used for this coefficient be the same as those that live within the
2270 {
2271 private:
2273
2274 public:
2275  /// Constructor with a quadrature function as input
2277
2279
2280  virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip);
2281
2283
2285 };
2286
2287 /// Flags that determine what storage optimizations to use in CoefficientVector
2288 enum class CoefficientStorage : int
2289 {
2290  FULL = 0, ///< Store the coefficient as a full QuadratureFunction.
2291  CONSTANTS = 1 << 0, ///< Store constants using only @a vdim entries.
2292  SYMMETRIC = 1 << 1, ///< Store the triangular part of symmetric matrices.
2293  COMPRESSED = CONSTANTS | SYMMETRIC ///< Enable all above compressions.
2294 };
2295
2297 {
2298  return CoefficientStorage(int(a) | int(b));
2299 }
2300
2302 {
2303  return int(a) & int(b);
2304 }
2305
2306
2307 /// @brief Class to represent a coefficient evaluated at quadrature points.
2308 ///
2309 /// In the general case, a CoefficientVector is the same as a QuadratureFunction
2310 /// with a coefficient projected onto it.
2311 ///
2312 /// This class allows for some "compression" of the coefficient data, according
2313 /// to the storage flags given by CoefficientStorage. For example, constant
2314 /// coefficients can be stored using only @a vdim values, and symmetric matrices
2315 /// can be stored using e.g. the upper triangular part of the matrix.
2317 {
2318 protected:
2319  CoefficientStorage storage; ///< Storage optimizations (see CoefficientStorage).
2320  int vdim; ///< Number of values per quadrature point.
2323 public:
2324  /// Create an empty CoefficientVector.
2327
2328  /// @brief Create a CoefficientVector from the given Coefficient and
2330  ///
2331  /// If @a coeff is NULL, it will be interpreted as a constant with value one.
2332  /// @sa CoefficientStorage for a description of @a storage_.
2335
2336  /// @brief Create a CoefficientVector from the given Coefficient and
2338  ///
2339  /// @sa CoefficientStorage for a description of @a storage_.
2342
2343  /// @brief Create a CoefficientVector from the given VectorCoefficient and
2345  ///
2346  /// @sa CoefficientStorage for a description of @a storage_.
2349
2350  /// @brief Create a CoefficientVector from the given MatrixCoefficient and
2352  ///
2353  /// @sa CoefficientStorage for a description of @a storage_.
2356
2357  /// @brief Evaluate the given Coefficient at the quadrature points defined by
2358  /// @ref qs.
2359  void Project(Coefficient &coeff);
2360
2361  /// @brief Evaluate the given VectorCoefficient at the quadrature points
2362  /// defined by @ref qs.
2363  ///
2364  /// @sa CoefficientVector for a description of the @a compress argument.
2365  void Project(VectorCoefficient &coeff);
2366
2367  /// @brief Evaluate the given MatrixCoefficient at the quadrature points
2368  /// defined by @ref qs.
2369  ///
2370  /// @sa CoefficientVector for a description of the @a compress argument.
2371  void Project(MatrixCoefficient &coeff, bool transpose=false);
2372
2373  /// @brief Project the transpose of @a coeff.
2374  ///
2375  /// @sa Project(MatrixCoefficient&, QuadratureSpace&, bool, bool)
2376  void ProjectTranspose(MatrixCoefficient &coeff);
2377
2378  /// Make this vector a reference to the given QuadratureFunction.
2380
2381  /// Set this vector to the given constant.
2382  void SetConstant(double constant);
2383
2384  /// Set this vector to the given constant vector.
2385  void SetConstant(const Vector &constant);
2386
2387  /// Set this vector to the given constant matrix.
2388  void SetConstant(const DenseMatrix &constant);
2389
2390  /// Set this vector to the given constant symmetric matrix.
2391  void SetConstant(const DenseSymmetricMatrix &constant);
2392
2393  /// Return the number of values per quadrature point.
2394  int GetVDim() const;
2395
2397 };
2398
2399 /** @brief Compute the Lp norm of a function f.
2400  \f$\| f \|_{Lp} = ( \int_\Omega | f |^p d\Omega)^{1/p} \f$ */
2401 double ComputeLpNorm(double p, Coefficient &coeff, Mesh &mesh,
2402  const IntegrationRule *irs[]);
2403
2404 /** @brief Compute the Lp norm of a vector function f = {f_i}_i=1...N.
2405  \f$\| f \|_{Lp} = ( \sum_i \| f_i \|_{Lp}^p )^{1/p} \f$ */
2406 double ComputeLpNorm(double p, VectorCoefficient &coeff, Mesh &mesh,
2407  const IntegrationRule *irs[]);
2408
2409 #ifdef MFEM_USE_MPI
2410 /** @brief Compute the global Lp norm of a function f.
2411  \f$\| f \|_{Lp} = ( \int_\Omega | f |^p d\Omega)^{1/p} \f$ */
2412 double ComputeGlobalLpNorm(double p, Coefficient &coeff, ParMesh &pmesh,
2413  const IntegrationRule *irs[]);
2414
2415 /** @brief Compute the global Lp norm of a vector function f = {f_i}_i=1...N.
2416  \f$\| f \|_{Lp} = ( \sum_i \| f_i \|_{Lp}^p )^{1/p} \f$ */
2417 double ComputeGlobalLpNorm(double p, VectorCoefficient &coeff, ParMesh &pmesh,
2418  const IntegrationRule *irs[]);
2419 #endif
2420
2421 }
2422
2423 #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:96
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.
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:517
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:197
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.
Scalar coefficient which returns the y-component of the evaluation point.
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.
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.
Vector beta
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)
CartesianZCoefficient CylindricalZCoefficient
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.
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.
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...
void ZeroCoefficient(int attr)
Remove a single VectorCoefficient for a particular attribute.
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)
std::function< double(const Vector &)> f(double mass_coeff)
Definition: lor_mms.hpp:30
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.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
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.
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...
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.
Make this vector a reference to the given QuadratureFunction.
MatrixCoefficient * GetBCoef() const
Return the second matrix coefficient.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
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.
Definition: qspace.hpp:25
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.
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 double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
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
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.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
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 double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient at ip.
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...
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.
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:31
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.
Scalar coefficient which returns the x-component of the evaluation point.
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.
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:864
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.
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.
CartesianCoefficient(int comp_)
comp_ index of the desired component (0 -> x, 1 -> y, 2 -> z)
MatrixCoefficient * GetACoef() const
Return the matrix coefficient.
Scalar coefficient defined as the determinant of a matrix coefficient.
Scalar coefficient which returns the z-component of the evaluation point.
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:58
TransformedCoefficient(Coefficient *q1, Coefficient *q2, double(*F)(double, double))
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.
Fill the QuadratureFunction qf with the constant value.
Definition: coefficient.cpp:71
const Vector & GetB() const
Return the second vector constant.
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.
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.
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector 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.