MFEM  v4.1.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
coefficient.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 // Implementation of Coefficient class
13 
14 #include "fem.hpp"
15 
16 #include <cmath>
17 #include <limits>
18 
19 namespace mfem
20 {
21 
22 using namespace std;
23 
25  const IntegrationPoint & ip)
26 {
27  int att = T.Attribute;
28  return (constants(att-1));
29 }
30 
32  const IntegrationPoint & ip)
33 {
34  double x[3];
35  Vector transip(x, 3);
36 
37  T.Transform(ip, transip);
38 
39  if (Function)
40  {
41  return ((*Function)(transip));
42  }
43  else
44  {
45  return (*TDFunction)(transip, GetTime());
46  }
47 }
48 
50  const IntegrationPoint &ip)
51 {
52  return GridF -> GetValue (T.ElementNo, ip, Component);
53 }
54 
56  const IntegrationPoint &ip)
57 {
58  if (Q2)
59  {
60  return (*Transform2)(Q1->Eval(T, ip, GetTime()),
61  Q2->Eval(T, ip, GetTime()));
62  }
63  else
64  {
65  return (*Transform1)(Q1->Eval(T, ip, GetTime()));
66  }
67 }
68 
70 {
71  MFEM_VERIFY(vcenter.Size() <= 3,
72  "SetDeltaCenter::Maximum number of dim supported is 3")
73  for (int i = 0; i < vcenter.Size(); i++) { center[i] = vcenter[i]; }
74  sdim = vcenter.Size();
75 }
76 
78 {
79  vcenter.SetSize(sdim);
80  vcenter = center;
81 }
82 
84  const IntegrationPoint &ip)
85 {
86  double w = Scale();
87  return weight ? weight->Eval(T, ip, GetTime())*w : w;
88 }
89 
91  const IntegrationRule &ir)
92 {
93  Vector Mi;
94  M.SetSize(vdim, ir.GetNPoints());
95  for (int i = 0; i < ir.GetNPoints(); i++)
96  {
97  M.GetColumnReference(i, Mi);
98  const IntegrationPoint &ip = ir.IntPoint(i);
99  T.SetIntPoint(&ip);
100  Eval(Mi, T, ip);
101  }
102 }
103 
105  const IntegrationPoint &ip)
106 {
107  double x[3];
108  Vector transip(x, 3);
109 
110  T.Transform(ip, transip);
111 
112  V.SetSize(vdim);
113  if (Function)
114  {
115  (*Function)(transip, V);
116  }
117  else
118  {
119  (*TDFunction)(transip, GetTime(), V);
120  }
121  if (Q)
122  {
123  V *= Q->Eval(T, ip, GetTime());
124  }
125 }
126 
128  : VectorCoefficient(dim), Coeff(dim), ownCoeff(dim)
129 {
130  for (int i = 0; i < dim; i++)
131  {
132  Coeff[i] = NULL;
133  ownCoeff[i] = true;
134  }
135 }
136 
137 void VectorArrayCoefficient::Set(int i, Coefficient *c, bool own)
138 {
139  if (ownCoeff[i]) { delete Coeff[i]; }
140  Coeff[i] = c;
141  ownCoeff[i] = own;
142 }
143 
145 {
146  for (int i = 0; i < vdim; i++)
147  {
148  if (ownCoeff[i]) { delete Coeff[i]; }
149  }
150 }
151 
153  const IntegrationPoint &ip)
154 {
155  V.SetSize(vdim);
156  for (int i = 0; i < vdim; i++)
157  {
158  V(i) = this->Eval(i, T, ip);
159  }
160 }
161 
163  GridFunction *gf)
164  : VectorCoefficient ((gf) ? gf -> VectorDim() : 0)
165 {
166  GridFunc = gf;
167 }
168 
170 {
171  GridFunc = gf; vdim = (gf) ? gf -> VectorDim() : 0;
172 }
173 
175  const IntegrationPoint &ip)
176 {
177  GridFunc->GetVectorValue(T.ElementNo, ip, V);
178 }
179 
182 {
183  GridFunc->GetVectorValues(T, ir, M);
184 }
185 
187  GridFunction *gf)
188  : VectorCoefficient((gf) ?
189  gf -> FESpace() -> GetMesh() -> SpaceDimension() : 0)
190 {
191  GridFunc = gf;
192 }
193 
195 {
196  GridFunc = gf; vdim = (gf) ?
197  gf -> FESpace() -> GetMesh() -> SpaceDimension() : 0;
198 }
199 
201  const IntegrationPoint &ip)
202 {
203  GridFunc->GetGradient(T, V);
204 }
205 
208 {
209  GridFunc->GetGradients(T, ir, M);
210 }
211 
213  GridFunction *gf)
214  : VectorCoefficient ((gf) ?
215  gf -> FESpace() -> GetMesh() -> SpaceDimension() : 0)
216 {
217  GridFunc = gf;
218 }
219 
221 {
222  GridFunc = gf; vdim = (gf) ?
223  gf -> FESpace() -> GetMesh() -> SpaceDimension() : 0;
224 }
225 
227  const IntegrationPoint &ip)
228 {
229  GridFunc->GetCurl(T, V);
230 }
231 
233  GridFunction *gf) : Coefficient()
234 {
235  GridFunc = gf;
236 }
237 
239  const IntegrationPoint &ip)
240 {
241  return GridFunc->GetDivergence(T);
242 }
243 
245 {
246  dir = _d;
247  (*this).vdim = dir.Size();
248 }
249 
252 {
253  V = dir;
254  d.SetTime(GetTime());
255  V *= d.EvalDelta(T, ip);
256 }
257 
259  const IntegrationPoint &ip)
260 {
261  V.SetSize(vdim);
262  if (active_attr[T.Attribute-1])
263  {
264  c->SetTime(GetTime());
265  c->Eval(V, T, ip);
266  }
267  else
268  {
269  V = 0.0;
270  }
271 }
272 
275 {
276  if (active_attr[T.Attribute-1])
277  {
278  c->SetTime(GetTime());
279  c->Eval(M, T, ir);
280  }
281  else
282  {
283  M.SetSize(vdim, ir.GetNPoints());
284  M = 0.0;
285  }
286 }
287 
289  const IntegrationPoint &ip)
290 {
291  double x[3];
292  Vector transip(x, 3);
293 
294  T.Transform(ip, transip);
295 
296  K.SetSize(height, width);
297 
298  if (Function)
299  {
300  (*Function)(transip, K);
301  }
302  else if (TDFunction)
303  {
304  (*TDFunction)(transip, GetTime(), K);
305  }
306  else
307  {
308  K = mat;
309  }
310  if (Q)
311  {
312  K *= Q->Eval(T, ip, GetTime());
313  }
314 }
315 
317  : MatrixCoefficient (dim)
318 {
319  Coeff.SetSize(height*width);
320  ownCoeff.SetSize(height*width);
321  for (int i = 0; i < (height*width); i++)
322  {
323  Coeff[i] = NULL;
324  ownCoeff[i] = true;
325  }
326 }
327 
328 void MatrixArrayCoefficient::Set(int i, int j, Coefficient * c, bool own)
329 {
330  if (ownCoeff[i*width+j]) { delete Coeff[i*width+j]; }
331  Coeff[i*width+j] = c;
332  ownCoeff[i*width+j] = own;
333 }
334 
336 {
337  for (int i=0; i < height*width; i++)
338  {
339  if (ownCoeff[i]) { delete Coeff[i]; }
340  }
341 }
342 
344  const IntegrationPoint &ip)
345 {
346  for (int i = 0; i < height; i++)
347  {
348  for (int j = 0; j < width; j++)
349  {
350  K(i,j) = this->Eval(i, j, T, ip);
351  }
352  }
353 }
354 
356  const IntegrationPoint &ip)
357 {
358  if (active_attr[T.Attribute-1])
359  {
360  c->SetTime(GetTime());
361  c->Eval(K, T, ip);
362  }
363  else
364  {
365  K.SetSize(height, width);
366  K = 0.0;
367  }
368 }
369 
372  : a(&A), b(&B)
373 {
374  MFEM_ASSERT(A.GetVDim() == B.GetVDim(),
375  "InnerProductCoefficient: "
376  "Arguments have incompatible dimensions.");
377 }
378 
380  const IntegrationPoint &ip)
381 {
382  a->Eval(va, T, ip);
383  b->Eval(vb, T, ip);
384  return va * vb;
385 }
386 
389  : a(&A), b(&B), va(A.GetVDim()), vb(B.GetVDim())
390 {
391  MFEM_ASSERT(A.GetVDim() == 2 && B.GetVDim() == 2,
392  "VectorRotProductCoefficient: "
393  "Arguments must have dimension equal to two.");
394 }
395 
397  const IntegrationPoint &ip)
398 {
399  a->Eval(va, T, ip);
400  b->Eval(vb, T, ip);
401  return va[0] * vb[1] - va[1] * vb[0];
402 }
403 
405  : a(&A), ma(A.GetHeight(), A.GetWidth())
406 {
407  MFEM_ASSERT(A.GetHeight() == A.GetWidth(),
408  "DeterminantCoefficient: "
409  "Argument must be a square matrix.");
410 }
411 
413  const IntegrationPoint &ip)
414 {
415  a->Eval(ma, T, ip);
416  return ma.Det();
417 }
418 
421  double _alpha, double _beta)
422  : VectorCoefficient(A.GetVDim()), a(&A), b(&B), alpha(_alpha), beta(_beta),
423  va(A.GetVDim())
424 {
425  MFEM_ASSERT(A.GetVDim() == B.GetVDim(),
426  "VectorSumCoefficient: "
427  "Arguments must have the same dimension.");
428 }
429 
431  const IntegrationPoint &ip)
432 {
433  b->Eval(V, T, ip);
434  if ( beta != 1.0 ) { V *= beta; }
435  a->Eval(va, T, ip);
436  V.Add(alpha, va);
437 }
438 
440  Coefficient &A,
442  : VectorCoefficient(B.GetVDim()), a(&A), b(&B)
443 {}
444 
446  const IntegrationPoint &ip)
447 {
448  double sa = a->Eval(T, ip);
449  b->Eval(V, T, ip);
450  V *= sa;
451 }
452 
456  : VectorCoefficient(3), a(&A), b(&B), va(A.GetVDim()), vb(B.GetVDim())
457 {
458  MFEM_ASSERT(A.GetVDim() == 3 && B.GetVDim() == 3,
459  "VectorCrossProductCoefficient: "
460  "Arguments must have dimension equal to three.");
461 }
462 
464  const IntegrationPoint &ip)
465 {
466  a->Eval(va, T, ip);
467  b->Eval(vb, T, ip);
468  V.SetSize(3);
469  V[0] = va[1] * vb[2] - va[2] * vb[1];
470  V[1] = va[2] * vb[0] - va[0] * vb[2];
471  V[2] = va[0] * vb[1] - va[1] * vb[0];
472 }
473 
476  : VectorCoefficient(A.GetHeight()), a(&A), b(&B),
477  ma(A.GetHeight(), A.GetWidth()), vb(B.GetVDim())
478 {
479  MFEM_ASSERT(A.GetWidth() == B.GetVDim(),
480  "MatVecCoefficient: Arguments have incompatible dimensions.");
481 }
482 
484  const IntegrationPoint &ip)
485 {
486  a->Eval(ma, T, ip);
487  b->Eval(vb, T, ip);
488  ma.Mult(vb, V);
489 }
490 
492  const IntegrationPoint &ip)
493 {
494  M.SetSize(dim);
495  M = 0.0;
496  for (int d=0; d<dim; d++) { M(d,d) = 1.0; }
497 }
498 
501  double _alpha, double _beta)
502  : MatrixCoefficient(A.GetHeight(), A.GetWidth()),
503  a(&A), b(&B), alpha(_alpha), beta(_beta),
504  ma(A.GetHeight(), A.GetWidth())
505 {
506  MFEM_ASSERT(A.GetHeight() == B.GetHeight() && A.GetWidth() == B.GetWidth(),
507  "MatrixSumCoefficient: "
508  "Arguments must have the same dimensions.");
509 }
510 
512  const IntegrationPoint &ip)
513 {
514  b->Eval(M, T, ip);
515  if ( beta != 1.0 ) { M *= beta; }
516  a->Eval(ma, T, ip);
517  M.Add(alpha, ma);
518 }
519 
521  Coefficient &A,
523  : MatrixCoefficient(B.GetHeight(), B.GetWidth()), a(&A), b(&B)
524 {}
525 
528  const IntegrationPoint &ip)
529 {
530  double sa = a->Eval(T, ip);
531  b->Eval(M, T, ip);
532  M *= sa;
533 }
534 
536  : MatrixCoefficient(A.GetWidth(), A.GetHeight()), a(&A)
537 {}
538 
541  const IntegrationPoint &ip)
542 {
543  a->Eval(M, T, ip);
544  M.Transpose();
545 }
546 
548  : MatrixCoefficient(A.GetHeight(), A.GetWidth()), a(&A)
549 {
550  MFEM_ASSERT(A.GetHeight() == A.GetWidth(),
551  "InverseMatrixCoefficient: "
552  "Argument must be a square matrix.");
553 }
554 
557  const IntegrationPoint &ip)
558 {
559  a->Eval(M, T, ip);
560  M.Invert();
561 }
562 
565  : MatrixCoefficient(A.GetVDim(), B.GetVDim()), a(&A), b(&B),
566  va(A.GetVDim()), vb(B.GetVDim())
567 {}
568 
570  const IntegrationPoint &ip)
571 {
572  a->Eval(va, T, ip);
573  b->Eval(vb, T, ip);
574  M.SetSize(va.Size(), vb.Size());
575  for (int i=0; i<va.Size(); i++)
576  {
577  for (int j=0; j<vb.Size(); j++)
578  {
579  M(i, j) = va[i] * vb[j];
580  }
581  }
582 }
583 
584 double LpNormLoop(double p, Coefficient &coeff, Mesh &mesh,
585  const IntegrationRule *irs[])
586 {
587  double norm = 0.0;
589 
590  for (int i = 0; i < mesh.GetNE(); i++)
591  {
592  tr = mesh.GetElementTransformation(i);
593  const IntegrationRule &ir = *irs[mesh.GetElementType(i)];
594  for (int j = 0; j < ir.GetNPoints(); j++)
595  {
596  const IntegrationPoint &ip = ir.IntPoint(j);
597  tr->SetIntPoint(&ip);
598  double val = fabs(coeff.Eval(*tr, ip));
599  if (p < infinity())
600  {
601  norm += ip.weight * tr->Weight() * pow(val, p);
602  }
603  else
604  {
605  if (norm < val)
606  {
607  norm = val;
608  }
609  }
610  }
611  }
612  return norm;
613 }
614 
615 double LpNormLoop(double p, VectorCoefficient &coeff, Mesh &mesh,
616  const IntegrationRule *irs[])
617 {
618  double norm = 0.0;
620  int vdim = coeff.GetVDim();
621  Vector vval(vdim);
622  double val;
623 
624  for (int i = 0; i < mesh.GetNE(); i++)
625  {
626  tr = mesh.GetElementTransformation(i);
627  const IntegrationRule &ir = *irs[mesh.GetElementType(i)];
628  for (int j = 0; j < ir.GetNPoints(); j++)
629  {
630  const IntegrationPoint &ip = ir.IntPoint(j);
631  tr->SetIntPoint(&ip);
632  coeff.Eval(vval, *tr, ip);
633  if (p < infinity())
634  {
635  for (int idim(0); idim < vdim; ++idim)
636  {
637  norm += ip.weight * tr->Weight() * pow(fabs( vval(idim) ), p);
638  }
639  }
640  else
641  {
642  for (int idim(0); idim < vdim; ++idim)
643  {
644  val = fabs(vval(idim));
645  if (norm < val)
646  {
647  norm = val;
648  }
649  }
650  }
651  }
652  }
653 
654  return norm;
655 }
656 
657 double ComputeLpNorm(double p, Coefficient &coeff, Mesh &mesh,
658  const IntegrationRule *irs[])
659 {
660  double norm = LpNormLoop(p, coeff, mesh, irs);
661 
662  if (p < infinity())
663  {
664  // negative quadrature weights may cause norm to be negative
665  if (norm < 0.0)
666  {
667  norm = -pow(-norm, 1.0/p);
668  }
669  else
670  {
671  norm = pow(norm, 1.0/p);
672  }
673  }
674 
675  return norm;
676 }
677 
678 double ComputeLpNorm(double p, VectorCoefficient &coeff, Mesh &mesh,
679  const IntegrationRule *irs[])
680 {
681  double norm = LpNormLoop(p, coeff, mesh, irs);
682 
683  if (p < infinity())
684  {
685  // negative quadrature weights may cause norm to be negative
686  if (norm < 0.0)
687  {
688  norm = -pow(-norm, 1.0/p);
689  }
690  else
691  {
692  norm = pow(norm, 1.0/p);
693  }
694  }
695 
696  return norm;
697 }
698 
699 #ifdef MFEM_USE_MPI
700 double ComputeGlobalLpNorm(double p, Coefficient &coeff, ParMesh &pmesh,
701  const IntegrationRule *irs[])
702 {
703  double loc_norm = LpNormLoop(p, coeff, pmesh, irs);
704  double glob_norm = 0;
705 
706  MPI_Comm comm = pmesh.GetComm();
707 
708  if (p < infinity())
709  {
710  MPI_Allreduce(&loc_norm, &glob_norm, 1, MPI_DOUBLE, MPI_SUM, comm);
711 
712  // negative quadrature weights may cause norm to be negative
713  if (glob_norm < 0.0)
714  {
715  glob_norm = -pow(-glob_norm, 1.0/p);
716  }
717  else
718  {
719  glob_norm = pow(glob_norm, 1.0/p);
720  }
721  }
722  else
723  {
724  MPI_Allreduce(&loc_norm, &glob_norm, 1, MPI_DOUBLE, MPI_MAX, comm);
725  }
726 
727  return glob_norm;
728 }
729 
730 double ComputeGlobalLpNorm(double p, VectorCoefficient &coeff, ParMesh &pmesh,
731  const IntegrationRule *irs[])
732 {
733  double loc_norm = LpNormLoop(p, coeff, pmesh, irs);
734  double glob_norm = 0;
735 
736  MPI_Comm comm = pmesh.GetComm();
737 
738  if (p < infinity())
739  {
740  MPI_Allreduce(&loc_norm, &glob_norm, 1, MPI_DOUBLE, MPI_SUM, comm);
741 
742  // negative quadrature weights may cause norm to be negative
743  if (glob_norm < 0.0)
744  {
745  glob_norm = -pow(-glob_norm, 1.0/p);
746  }
747  else
748  {
749  glob_norm = pow(glob_norm, 1.0/p);
750  }
751  }
752  else
753  {
754  MPI_Allreduce(&loc_norm, &glob_norm, 1, MPI_DOUBLE, MPI_MAX, comm);
755  }
756 
757  return glob_norm;
758 }
759 #endif
760 
761 }
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition: intrules.hpp:245
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 ...
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:27
void Set(int i, int j, Coefficient *c, bool own=true)
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
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 ...
TransposeMatrixCoefficient(MatrixCoefficient &A)
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)
void GetVectorValue(int i, const IntegrationPoint &ip, Vector &val) const
Definition: gridfunc.cpp:408
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:407
double Eval(int i, ElementTransformation &T, const IntegrationPoint &ip)
Evaluates i&#39;th component of the vector.
double Det() const
Definition: densemat.cpp:449
void SetDeltaCenter(const Vector &center)
Definition: coefficient.cpp:69
Element::Type GetElementType(int i) const
Returns the type of element i.
Definition: mesh.cpp:4656
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 ...
void SetIntPoint(const IntegrationPoint *ip)
Definition: eltrans.hpp:56
GradientGridFunctionCoefficient(GridFunction *gf)
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int Size() const
Returns the size of the vector.
Definition: vector.hpp:157
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 ...
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:693
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
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 ...
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
Definition: coefficient.cpp:55
void SetTime(double t)
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate coefficient.
Definition: coefficient.cpp:31
DeterminantCoefficient(MatrixCoefficient &A)
void SetGridFunction(GridFunction *gf)
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
double Eval(int i, int j, ElementTransformation &T, const IntegrationPoint &ip)
InnerProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
Definition: intrules.hpp:248
double ComputeLpNorm(double p, Coefficient &coeff, Mesh &mesh, const IntegrationRule *irs[])
double LpNormLoop(double p, Coefficient &coeff, Mesh &mesh, const IntegrationRule *irs[])
void Add(const double c, const DenseMatrix &A)
Adds the matrix A multiplied by the number c to the matrix.
Definition: densemat.cpp:542
ScalarVectorProductCoefficient(Coefficient &A, VectorCoefficient &B)
double b
Definition: lissajous.cpp:42
void SetDirection(const Vector &_d)
OuterProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
void Invert()
Replaces the current matrix with its inverse.
Definition: densemat.cpp:635
void Set(int i, Coefficient *c, bool own=true)
Sets coefficient in the vector.
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 ...
virtual ~VectorArrayCoefficient()
Destroys vector coefficient.
int GetVDim()
Returns dimension of the vector.
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
VectorSumCoefficient(VectorCoefficient &A, VectorCoefficient &B, double _alpha=1.0, double _beta=1.0)
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
void SetTime(double t)
Definition: coefficient.hpp:39
void SetGridFunction(GridFunction *gf)
void GetDeltaCenter(Vector &center)
Definition: coefficient.cpp:77
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient.
void Transpose()
(*this) = (*this)^t
Definition: densemat.cpp:1376
MPI_Comm GetComm() const
Definition: pmesh.hpp:230
MatVecCoefficient(MatrixCoefficient &A, VectorCoefficient &B)
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
void SetSize(int nsize)
Change logical size of the array, keep existing entries.
Definition: array.hpp:635
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient function.
Definition: coefficient.cpp:24
void GetColumnReference(int c, Vector &col)
Definition: densemat.hpp:270
void SetTime(double t)
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 ...
virtual void EvalDelta(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Return the specified direction vector multiplied by the value returned by DeltaCoefficient::EvalDelta...
double ComputeGlobalLpNorm(double p, Coefficient &coeff, ParMesh &pmesh, const IntegrationRule *irs[])
double GetDivergence(ElementTransformation &tr) const
Definition: gridfunc.cpp:1010
virtual double EvalDelta(ElementTransformation &T, const IntegrationPoint &ip)
Return the Scale() multiplied by the weight Coefficient, if any.
Definition: coefficient.cpp:83
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
Definition: coefficient.cpp:49
void SetGridFunction(GridFunction *gf)
double a
Definition: lissajous.cpp:41
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 ...
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 ...
Vector & Add(const double a, const Vector &Va)
(*this) += a * Va
Definition: vector.cpp:190
Class for integration point with weight.
Definition: intrules.hpp:25
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the coefficient in the element described by T at the point ip.
void GetElementTransformation(int i, IsoparametricTransformation *ElTr)
Definition: mesh.cpp:338
MatrixSumCoefficient(MatrixCoefficient &A, MatrixCoefficient &B, double _alpha=1.0, double _beta=1.0)
int dim
Definition: ex24.cpp:43
void GetGradients(ElementTransformation &tr, const IntegrationRule &ir, DenseMatrix &grad) const
Definition: gridfunc.cpp:1113
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
void GetVectorValues(ElementTransformation &T, const IntegrationRule &ir, DenseMatrix &vals) const
Definition: gridfunc.cpp:620
double infinity()
Define a shortcut for std::numeric_limits&lt;double&gt;::infinity()
Definition: vector.hpp:42
DivergenceGridFunctionCoefficient(GridFunction *gf)
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.
const double alpha
Definition: ex15.cpp:336
void GetCurl(ElementTransformation &tr, Vector &curl) const
Definition: gridfunc.cpp:1044
Vector data type.
Definition: vector.hpp:48
void Mult(const double *x, double *y) const
Matrix vector multiplication.
Definition: densemat.cpp:173
virtual void Transform(const IntegrationPoint &, Vector &)=0
VectorRotProductCoefficient(VectorCoefficient &A, VectorCoefficient &B)
virtual void Eval(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
virtual void Eval(DenseMatrix &M, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the matrix coefficient in the element described by T at the point ip, storing the result in ...
CurlGridFunctionCoefficient(GridFunction *gf)
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
Definition: densemat.hpp:90
Class for parallel meshes.
Definition: pmesh.hpp:32
InverseMatrixCoefficient(MatrixCoefficient &A)
ScalarMatrixProductCoefficient(Coefficient &A, MatrixCoefficient &B)
void GetGradient(ElementTransformation &tr, Vector &grad) const
Definition: gridfunc.cpp:1095