MFEM  v4.1.0 Finite element discretization library
lininteg.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, 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
13 #include <cmath>
14 #include "fem.hpp"
15
16 namespace mfem
17 {
18
20  const FiniteElement &el, FaceElementTransformations &Tr, Vector &elvect)
21 {
22  mfem_error("LinearFormIntegrator::AssembleRHSElementVect(...)");
23 }
24
25
28  Vector &elvect)
29 {
30  int dof = el.GetDof();
31
32  shape.SetSize(dof); // vector of size dof
33  elvect.SetSize(dof);
34  elvect = 0.0;
35
36  const IntegrationRule *ir = IntRule;
37  if (ir == NULL)
38  {
39  // ir = &IntRules.Get(el.GetGeomType(),
40  // oa * el.GetOrder() + ob + Tr.OrderW());
41  ir = &IntRules.Get(el.GetGeomType(), oa * el.GetOrder() + ob);
42  }
43
44  for (int i = 0; i < ir->GetNPoints(); i++)
45  {
46  const IntegrationPoint &ip = ir->IntPoint(i);
47
48  Tr.SetIntPoint (&ip);
49  double val = Tr.Weight() * Q.Eval(Tr, ip);
50
51  el.CalcShape(ip, shape);
52
53  add(elvect, ip.weight * val, shape, elvect);
54  }
55 }
56
58  const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
59 {
60  MFEM_ASSERT(delta != NULL, "coefficient must be DeltaCoefficient");
61  elvect.SetSize(fe.GetDof());
62  fe.CalcPhysShape(Trans, elvect);
63  elvect *= delta->EvalDelta(Trans, Trans.GetIntPoint());
64 }
65
66
68  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
69 {
70  int dof = el.GetDof();
71
72  shape.SetSize(dof); // vector of size dof
73  elvect.SetSize(dof);
74  elvect = 0.0;
75
76  const IntegrationRule *ir = IntRule;
77  if (ir == NULL)
78  {
79  int intorder = oa * el.GetOrder() + ob; // <----------
80  ir = &IntRules.Get(el.GetGeomType(), intorder);
81  }
82
83  for (int i = 0; i < ir->GetNPoints(); i++)
84  {
85  const IntegrationPoint &ip = ir->IntPoint(i);
86
87  Tr.SetIntPoint (&ip);
88  double val = Tr.Weight() * Q.Eval(Tr, ip);
89
90  el.CalcShape(ip, shape);
91
92  add(elvect, ip.weight * val, shape, elvect);
93  }
94 }
95
97  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
98 {
99  int dim = el.GetDim()+1;
100  int dof = el.GetDof();
101  Vector nor(dim), Qvec;
102
103  shape.SetSize(dof);
104  elvect.SetSize(dof);
105  elvect = 0.0;
106
107  const IntegrationRule *ir = IntRule;
108  if (ir == NULL)
109  {
110  int intorder = oa * el.GetOrder() + ob; // <----------
111  ir = &IntRules.Get(el.GetGeomType(), intorder);
112  }
113
114  for (int i = 0; i < ir->GetNPoints(); i++)
115  {
116  const IntegrationPoint &ip = ir->IntPoint(i);
117
118  Tr.SetIntPoint(&ip);
119  CalcOrtho(Tr.Jacobian(), nor);
120  Q.Eval(Qvec, Tr, ip);
121
122  el.CalcShape(ip, shape);
123
125  }
126 }
127
129  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
130 {
131  int dim = el.GetDim()+1;
132  int dof = el.GetDof();
133  Vector tangent(dim), Qvec;
134
135  shape.SetSize(dof);
136  elvect.SetSize(dof);
137  elvect = 0.0;
138
139  if (dim != 2)
140  {
141  mfem_error("These methods make sense only in 2D problems.");
142  }
143
144  const IntegrationRule *ir = IntRule;
145  if (ir == NULL)
146  {
147  int intorder = oa * el.GetOrder() + ob; // <----------
148  ir = &IntRules.Get(el.GetGeomType(), intorder);
149  }
150
151  for (int i = 0; i < ir->GetNPoints(); i++)
152  {
153  const IntegrationPoint &ip = ir->IntPoint(i);
154
155  Tr.SetIntPoint(&ip);
156  const DenseMatrix &Jac = Tr.Jacobian();
157  tangent(0) = Jac(0,0);
158  tangent(1) = Jac(1,0);
159
160  Q.Eval(Qvec, Tr, ip);
161
162  el.CalcShape(ip, shape);
163
165  }
166 }
167
169  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
170 {
171  int vdim = Q.GetVDim();
172  int dof = el.GetDof();
173
174  double val,cf;
175
176  shape.SetSize(dof); // vector of size dof
177
178  elvect.SetSize(dof * vdim);
179  elvect = 0.0;
180
181  const IntegrationRule *ir = IntRule;
182  if (ir == NULL)
183  {
184  int intorder = 2*el.GetOrder();
185  ir = &IntRules.Get(el.GetGeomType(), intorder);
186  }
187
188  for (int i = 0; i < ir->GetNPoints(); i++)
189  {
190  const IntegrationPoint &ip = ir->IntPoint(i);
191
192  Tr.SetIntPoint (&ip);
193  val = Tr.Weight();
194
195  el.CalcShape(ip, shape);
196  Q.Eval (Qvec, Tr, ip);
197
198  for (int k = 0; k < vdim; k++)
199  {
200  cf = val * Qvec(k);
201
202  for (int s = 0; s < dof; s++)
203  {
204  elvect(dof*k+s) += ip.weight * cf * shape(s);
205  }
206  }
207  }
208 }
209
211  const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
212 {
213  MFEM_ASSERT(vec_delta != NULL, "coefficient must be VectorDeltaCoefficient");
214  int vdim = Q.GetVDim();
215  int dof = fe.GetDof();
216
217  shape.SetSize(dof);
218  fe.CalcPhysShape(Trans, shape);
219
220  vec_delta->EvalDelta(Qvec, Trans, Trans.GetIntPoint());
221
222  elvect.SetSize(dof*vdim);
223  DenseMatrix elvec_as_mat(elvect.GetData(), dof, vdim);
224  MultVWt(shape, Qvec, elvec_as_mat);
225 }
226
227
229  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
230 {
231  int vdim = Q.GetVDim();
232  int dof = el.GetDof();
233
234  shape.SetSize(dof);
235  vec.SetSize(vdim);
236
237  elvect.SetSize(dof * vdim);
238  elvect = 0.0;
239
240  const IntegrationRule *ir = IntRule;
241  if (ir == NULL)
242  {
243  int intorder = 2*el.GetOrder();
244  ir = &IntRules.Get(el.GetGeomType(), intorder);
245  }
246
247  for (int i = 0; i < ir->GetNPoints(); i++)
248  {
249  const IntegrationPoint &ip = ir->IntPoint(i);
250
251  Q.Eval(vec, Tr, ip);
252  Tr.SetIntPoint (&ip);
253  vec *= Tr.Weight() * ip.weight;
254  el.CalcShape(ip, shape);
255  for (int k = 0; k < vdim; k++)
256  for (int s = 0; s < dof; s++)
257  {
258  elvect(dof*k+s) += vec(k) * shape(s);
259  }
260  }
261 }
262
264  const FiniteElement &el, FaceElementTransformations &Tr, Vector &elvect)
265 {
266  int vdim = Q.GetVDim();
267  int dof = el.GetDof();
268
269  shape.SetSize(dof);
270  vec.SetSize(vdim);
271
272  elvect.SetSize(dof * vdim);
273  elvect = 0.0;
274
275  const IntegrationRule *ir = IntRule;
276  if (ir == NULL)
277  {
278  int intorder = 2*el.GetOrder();
279  ir = &IntRules.Get(Tr.FaceGeom, intorder);
280  }
281
282  for (int i = 0; i < ir->GetNPoints(); i++)
283  {
284  const IntegrationPoint &ip = ir->IntPoint(i);
285  IntegrationPoint eip;
286  Tr.Loc1.Transform(ip, eip);
287
288  Tr.Face->SetIntPoint(&ip);
289  Q.Eval(vec, *Tr.Face, ip);
290  vec *= Tr.Face->Weight() * ip.weight;
291  el.CalcShape(eip, shape);
292  for (int k = 0; k < vdim; k++)
293  {
294  for (int s = 0; s < dof; s++)
295  {
296  elvect(dof*k+s) += vec(k) * shape(s);
297  }
298  }
299  }
300 }
301
302
304  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
305 {
306  int dof = el.GetDof();
307  int spaceDim = Tr.GetSpaceDim();
308
309  vshape.SetSize(dof,spaceDim);
310  vec.SetSize(spaceDim);
311
312  elvect.SetSize(dof);
313  elvect = 0.0;
314
315  const IntegrationRule *ir = IntRule;
316  if (ir == NULL)
317  {
318  // int intorder = 2*el.GetOrder() - 1; // ok for O(h^{k+1}) conv. in L2
319  int intorder = 2*el.GetOrder();
320  ir = &IntRules.Get(el.GetGeomType(), intorder);
321  }
322
323  for (int i = 0; i < ir->GetNPoints(); i++)
324  {
325  const IntegrationPoint &ip = ir->IntPoint(i);
326
327  Tr.SetIntPoint (&ip);
328  el.CalcVShape(Tr, vshape);
329
330  QF.Eval (vec, Tr, ip);
331  vec *= ip.weight * Tr.Weight();
332
334  }
335 }
336
338  const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
339 {
340  MFEM_ASSERT(vec_delta != NULL, "coefficient must be VectorDeltaCoefficient");
341  int dof = fe.GetDof();
342  int spaceDim = Trans.GetSpaceDim();
343
344  vshape.SetSize(dof, spaceDim);
345  fe.CalcPhysVShape(Trans, vshape);
346
347  vec_delta->EvalDelta(vec, Trans, Trans.GetIntPoint());
348
349  elvect.SetSize(dof);
350  vshape.Mult(vec, elvect);
351 }
352
354  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
355 {
356  int dim = el.GetDim()+1;
357  int dof = el.GetDof();
358
359  shape.SetSize (dof);
360  nor.SetSize (dim);
361  elvect.SetSize (dim*dof);
362
363  const IntegrationRule *ir = IntRule;
364  if (ir == NULL)
365  {
366  ir = &IntRules.Get(el.GetGeomType(), el.GetOrder() + 1);
367  }
368
369  elvect = 0.0;
370  for (int i = 0; i < ir->GetNPoints(); i++)
371  {
372  const IntegrationPoint &ip = ir->IntPoint(i);
373  Tr.SetIntPoint (&ip);
374  CalcOrtho(Tr.Jacobian(), nor);
375  el.CalcShape (ip, shape);
376  nor *= Sign * ip.weight * F -> Eval (Tr, ip);
377  for (int j = 0; j < dof; j++)
378  for (int k = 0; k < dim; k++)
379  {
380  elvect(dof*k+j) += nor(k) * shape(j);
381  }
382  }
383 }
384
385
387  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
388 {
389  int dof = el.GetDof();
390
391  shape.SetSize(dof);
392  elvect.SetSize(dof);
393  elvect = 0.0;
394
395  const IntegrationRule *ir = IntRule;
396  if (ir == NULL)
397  {
398  int intorder = 2*el.GetOrder(); // <----------
399  if (F == NULL)
400  {
401  intorder -= el.GetOrder() + 1;
402  }
403  ir = &IntRules.Get(el.GetGeomType(), intorder);
404  }
405
406  for (int i = 0; i < ir->GetNPoints(); i++)
407  {
408  const IntegrationPoint &ip = ir->IntPoint(i);
409  el.CalcShape(ip, shape);
410
411  double val = ip.weight;
412  if (F)
413  {
414  Tr.SetIntPoint (&ip);
415  val *= F->Eval(Tr, ip);
416  }
417
419  }
420 }
421
422
424  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
425 {
426  int dof = el.GetDof();
427  DenseMatrix vshape(dof, 2);
428  Vector f_loc(3);
429  Vector f_hat(2);
430
431  elvect.SetSize(dof);
432  elvect = 0.0;
433
434  const IntegrationRule *ir = IntRule;
435  if (ir == NULL)
436  {
437  int intorder = 2*el.GetOrder(); // <----------
438  ir = &IntRules.Get(el.GetGeomType(), intorder);
439  }
440
441  for (int i = 0; i < ir->GetNPoints(); i++)
442  {
443  const IntegrationPoint &ip = ir->IntPoint(i);
444
445  Tr.SetIntPoint(&ip);
446  f.Eval(f_loc, Tr, ip);
447  Tr.Jacobian().MultTranspose(f_loc, f_hat);
448  el.CalcVShape(ip, vshape);
449
450  Swap<double>(f_hat(0), f_hat(1));
451  f_hat(0) = -f_hat(0);
452  f_hat *= ip.weight;
454  }
455 }
456
457
459  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
460 {
461  mfem_error("BoundaryFlowIntegrator::AssembleRHSElementVect\n"
462  " is not implemented as boundary integrator!\n"
465 }
466
468  const FiniteElement &el, FaceElementTransformations &Tr, Vector &elvect)
469 {
470  int dim, ndof, order;
471  double un, w, vu_data[3], nor_data[3];
472
473  dim = el.GetDim();
474  ndof = el.GetDof();
475  Vector vu(vu_data, dim), nor(nor_data, dim);
476
477  const IntegrationRule *ir = IntRule;
478  if (ir == NULL)
479  {
480  // Assuming order(u)==order(mesh)
481  order = Tr.Elem1->OrderW() + 2*el.GetOrder();
482  if (el.Space() == FunctionSpace::Pk)
483  {
484  order++;
485  }
486  ir = &IntRules.Get(Tr.FaceGeom, order);
487  }
488
489  shape.SetSize(ndof);
490  elvect.SetSize(ndof);
491  elvect = 0.0;
492
493  for (int p = 0; p < ir->GetNPoints(); p++)
494  {
495  const IntegrationPoint &ip = ir->IntPoint(p);
496  IntegrationPoint eip;
497  Tr.Loc1.Transform(ip, eip);
498  el.CalcShape(eip, shape);
499
500  Tr.Face->SetIntPoint(&ip);
501
502  u->Eval(vu, *Tr.Elem1, eip);
503
504  if (dim == 1)
505  {
506  nor(0) = 2*eip.x - 1.0;
507  }
508  else
509  {
510  CalcOrtho(Tr.Face->Jacobian(), nor);
511  }
512
513  un = vu * nor;
514  w = 0.5*alpha*un - beta*fabs(un);
515  w *= ip.weight*f->Eval(*Tr.Elem1, eip);
517  }
518 }
519
520
522  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
523 {
524  mfem_error("DGDirichletLFIntegrator::AssembleRHSElementVect");
525 }
526
528  const FiniteElement &el, FaceElementTransformations &Tr, Vector &elvect)
529 {
530  int dim, ndof;
531  bool kappa_is_nonzero = (kappa != 0.);
532  double w;
533
534  dim = el.GetDim();
535  ndof = el.GetDof();
536
537  nor.SetSize(dim);
538  nh.SetSize(dim);
539  ni.SetSize(dim);
541  if (MQ)
542  {
543  mq.SetSize(dim);
544  }
545
546  shape.SetSize(ndof);
547  dshape.SetSize(ndof, dim);
548  dshape_dn.SetSize(ndof);
549
550  elvect.SetSize(ndof);
551  elvect = 0.0;
552
553  const IntegrationRule *ir = IntRule;
554  if (ir == NULL)
555  {
556  // a simple choice for the integration order; is this OK?
557  int order = 2*el.GetOrder();
558  ir = &IntRules.Get(Tr.FaceGeom, order);
559  }
560
561  for (int p = 0; p < ir->GetNPoints(); p++)
562  {
563  const IntegrationPoint &ip = ir->IntPoint(p);
564  IntegrationPoint eip;
565
566  Tr.Loc1.Transform(ip, eip);
567  Tr.Face->SetIntPoint(&ip);
568  if (dim == 1)
569  {
570  nor(0) = 2*eip.x - 1.0;
571  }
572  else
573  {
574  CalcOrtho(Tr.Face->Jacobian(), nor);
575  }
576
577  el.CalcShape(eip, shape);
578  el.CalcDShape(eip, dshape);
579  Tr.Elem1->SetIntPoint(&eip);
580  // compute uD through the face transformation
581  w = ip.weight * uD->Eval(*Tr.Face, ip) / Tr.Elem1->Weight();
582  if (!MQ)
583  {
584  if (Q)
585  {
586  w *= Q->Eval(*Tr.Elem1, eip);
587  }
588  ni.Set(w, nor);
589  }
590  else
591  {
592  nh.Set(w, nor);
593  MQ->Eval(mq, *Tr.Elem1, eip);
594  mq.MultTranspose(nh, ni);
595  }
598
601
602  if (kappa_is_nonzero)
603  {
605  }
606  }
607 }
608
609
611  const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
612 {
613  mfem_error("DGElasticityDirichletLFIntegrator::AssembleRHSElementVect");
614 }
615
617  const FiniteElement &el, FaceElementTransformations &Tr, Vector &elvect)
618 {
619  MFEM_ASSERT(Tr.Elem2No < 0, "interior boundary is not supported");
620
622  Vector shape;
626  Vector nor;
629  Vector u_dir;
630 #endif
631
632  const int dim = el.GetDim();
633  const int ndofs = el.GetDof();
634  const int nvdofs = dim*ndofs;
635
636  elvect.SetSize(nvdofs);
637  elvect = 0.0;
638
640  shape.SetSize(ndofs);
641  dshape.SetSize(ndofs, dim);
642  dshape_ps.SetSize(ndofs, dim);
643  nor.SetSize(dim);
644  dshape_dn.SetSize(ndofs);
645  dshape_du.SetSize(ndofs);
646  u_dir.SetSize(dim);
647
648  const IntegrationRule *ir = IntRule;
649  if (ir == NULL)
650  {
651  const int order = 2*el.GetOrder(); // <-----
652  ir = &IntRules.Get(Tr.FaceGeom, order);
653  }
654
655  for (int pi = 0; pi < ir->GetNPoints(); ++pi)
656  {
657  const IntegrationPoint &ip = ir->IntPoint(pi);
658  IntegrationPoint eip;
659  Tr.Loc1.Transform(ip, eip);
660  Tr.Face->SetIntPoint(&ip);
661  Tr.Elem1->SetIntPoint(&eip);
662
663  // Evaluate the Dirichlet b.c. using the face transformation.
664  uD.Eval(u_dir, *Tr.Face, ip);
665
666  el.CalcShape(eip, shape);
667  el.CalcDShape(eip, dshape);
668
671
672  if (dim == 1)
673  {
674  nor(0) = 2*eip.x - 1.0;
675  }
676  else
677  {
678  CalcOrtho(Tr.Face->Jacobian(), nor);
679  }
680
681  double wL, wM, jcoef;
682  {
683  const double w = ip.weight / Tr.Elem1->Weight();
684  wL = w * lambda->Eval(*Tr.Elem1, eip);
685  wM = w * mu->Eval(*Tr.Elem1, eip);
686  jcoef = kappa * (wL + 2.0*wM) * (nor*nor);
687  dshape_ps.Mult(nor, dshape_dn);
688  dshape_ps.Mult(u_dir, dshape_du);
689  }
690
691  // alpha < uD, (lambda div(v) I + mu (grad(v) + grad(v)^T)) . n > +
692  // + kappa < h^{-1} (lambda + 2 mu) uD, v >
693
694  // i = idof + ndofs * im
695  // v_phi(i,d) = delta(im,d) phi(idof)
696  // div(v_phi(i)) = dphi(idof,im)
697  // (grad(v_phi(i)))(k,l) = delta(im,k) dphi(idof,l)
698  //
699  // term 1:
700  // alpha < uD, lambda div(v_phi(i)) n >
701  // alpha lambda div(v_phi(i)) (uD.n) =
702  // alpha lambda dphi(idof,im) (uD.n) --> quadrature -->
703  // ip.weight/det(J1) alpha lambda (uD.nor) dshape_ps(idof,im) =
704  // alpha * wL * (u_dir*nor) * dshape_ps(idof,im)
705  // term 2:
706  // < alpha uD, mu grad(v_phi(i)).n > =
707  // alpha mu uD^T grad(v_phi(i)) n =
708  // alpha mu uD(k) delta(im,k) dphi(idof,l) n(l) =
709  // alpha mu uD(im) dphi(idof,l) n(l) --> quadrature -->
710  // ip.weight/det(J1) alpha mu uD(im) dshape_ps(idof,l) nor(l) =
711  // alpha * wM * u_dir(im) * dshape_dn(idof)
712  // term 3:
713  // < alpha uD, mu (grad(v_phi(i)))^T n > =
714  // alpha mu n^T grad(v_phi(i)) uD =
715  // alpha mu n(k) delta(im,k) dphi(idof,l) uD(l) =
716  // alpha mu n(im) dphi(idof,l) uD(l) --> quadrature -->
717  // ip.weight/det(J1) alpha mu nor(im) dshape_ps(idof,l) uD(l) =
718  // alpha * wM * nor(im) * dshape_du(idof)
719  // term j:
720  // < kappa h^{-1} (lambda + 2 mu) uD, v_phi(i) > =
721  // kappa/h (lambda + 2 mu) uD(k) v_phi(i,k) =
722  // kappa/h (lambda + 2 mu) uD(k) delta(im,k) phi(idof) =
723  // kappa/h (lambda + 2 mu) uD(im) phi(idof) --> quadrature -->
724  // [ 1/h = |nor|/det(J1) ]
725  // ip.weight/det(J1) |nor|^2 kappa (lambda + 2 mu) uD(im) phi(idof) =
726  // jcoef * u_dir(im) * shape(idof)
727
728  wM *= alpha;
729  const double t1 = alpha * wL * (u_dir*nor);
730  for (int im = 0, i = 0; im < dim; ++im)
731  {
732  const double t2 = wM * u_dir(im);
733  const double t3 = wM * nor(im);
734  const double tj = jcoef * u_dir(im);
735  for (int idof = 0; idof < ndofs; ++idof, ++i)
736  {
737  elvect(i) += (t1*dshape_ps(idof,im) + t2*dshape_dn(idof) +
738  t3*dshape_du(idof) + tj*shape(idof));
739  }
740  }
741  }
742 }
743
744 }
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition: intrules.hpp:245
Abstract class for Finite Elements.
Definition: fe.hpp:232
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:386
int GetDim() const
Returns the reference space dimension for the finite element.
Definition: fe.hpp:308
ElementTransformation * Face
Definition: eltrans.hpp:356
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
Definition: fe.cpp:40
void MultVWt(const Vector &v, const Vector &w, DenseMatrix &VWt)
Definition: densemat.cpp:2757
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
Definition: intrules.cpp:890
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 ...
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:407
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:458
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:610
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
Definition: table.cpp:476
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 ...
virtual void AssembleDeltaElementVect(const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
Assemble the delta coefficient at the IntegrationPoint set in Trans which is assumed to map to the de...
Definition: lininteg.cpp:57
void SetIntPoint(const IntegrationPoint *ip)
Definition: eltrans.hpp:56
Definition: densemat.cpp:2089
int GetOrder() const
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order...
Definition: fe.hpp:318
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:168
int Space() const
Returns the type of space on each element.
Definition: fe.hpp:328
const IntegrationRule * IntRule
Definition: lininteg.hpp:25
const IntegrationPoint & GetIntPoint()
Definition: eltrans.hpp:58
Polynomials of order k.
Definition: fe.hpp:215
double * GetData() const
Return a pointer to the beginning of the Vector data.
Definition: vector.hpp:166
void CalcOrtho(const DenseMatrix &J, Vector &n)
Definition: densemat.cpp:2303
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:67
void add(const Vector &v1, const Vector &v2, Vector &v)
Definition: vector.cpp:238
Geometry::Type GetGeomType() const
Returns the Geometry::Type of the reference element.
Definition: fe.hpp:311
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
Definition: intrules.hpp:248
void MultTranspose(const double *x, double *y) const
Multiply a vector with the transpose matrix.
Definition: densemat.cpp:201
IntegrationPointTransformation Loc1
Definition: eltrans.hpp:357
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:153
const DenseMatrix & Jacobian()
Return the Jacobian matrix of the transformation at the currently set IntegrationPoint, using the method SetIntPoint().
Definition: eltrans.hpp:71
int GetSpaceDim() const
Get the dimension of the target (physical) space.
Definition: eltrans.hpp:100
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:228
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const =0
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
VectorDeltaCoefficient * vec_delta
Definition: lininteg.hpp:51
MatrixCoefficient * MQ
Definition: lininteg.hpp:353
int GetVDim()
Returns dimension of the vector.
void AddMult(const Vector &x, Vector &y) const
y += A.x
Definition: densemat.cpp:224
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:128
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:521
void CalcPhysVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Equivalent to the CalcVShape() method with the same arguments.
Definition: fe.hpp:389
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)=0
int GetDof() const
Returns the number of degrees of freedom in the finite element.
Definition: fe.hpp:314
void CalcPhysShape(ElementTransformation &Trans, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in physical space at the point ...
Definition: fe.cpp:185
DeltaCoefficient * delta
Definition: lininteg.hpp:50
void Transform(const IntegrationPoint &, IntegrationPoint &)
Definition: eltrans.cpp:537
virtual void EvalDelta(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Return the specified direction vector multiplied by the value returned by DeltaCoefficient::EvalDelta...
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:353
virtual double EvalDelta(ElementTransformation &T, const IntegrationPoint &ip)
Return the Scale() multiplied by the weight Coefficient, if any.
Definition: coefficient.cpp:83
Vector & Set(const double a, const Vector &x)
(*this) = a * x
Definition: vector.cpp:205
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:423
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 void AssembleDeltaElementVect(const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
Assemble the delta coefficient at the IntegrationPoint set in Trans which is assumed to map to the de...
Definition: lininteg.cpp:337
int dim
Definition: ex24.cpp:43
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:96
ElementTransformation * Elem1
Definition: eltrans.hpp:356
virtual double Eval(ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the coefficient in the element described by T at the point ip.
Vector data type.
Definition: vector.hpp:48
void Mult(const double *x, double *y) const
Matrix vector multiplication.
Definition: densemat.cpp:173
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:303
virtual void AssembleDeltaElementVect(const FiniteElement &fe, ElementTransformation &Trans, Vector &elvect)
Assemble the delta coefficient at the IntegrationPoint set in Trans which is assumed to map to the de...
Definition: lininteg.cpp:210
virtual void AssembleRHSElementVect(const FiniteElement &el, ElementTransformation &Tr, Vector &elvect)
Definition: lininteg.cpp:26
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
Definition: densemat.hpp:90
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const =0
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...