23 void BilinearFormIntegrator::AssembleElementMatrix (
27 mfem_error (
"BilinearFormIntegrator::AssembleElementMatrix (...)\n"
28 " is not implemented fot this class.");
31 void BilinearFormIntegrator::AssembleElementMatrix2 (
35 mfem_error (
"BilinearFormIntegrator::AssembleElementMatrix2 (...)\n"
36 " is not implemented fot this class.");
39 void BilinearFormIntegrator::AssembleFaceMatrix (
43 mfem_error (
"BilinearFormIntegrator::AssembleFaceMatrix (...)\n"
44 " is not implemented fot this class.");
47 void BilinearFormIntegrator::AssembleFaceMatrix(
52 MFEM_ABORT(
"AssembleFaceMatrix (mixed form) is not implemented for this"
53 " Integrator class.");
56 void BilinearFormIntegrator::AssembleElementVector(
60 mfem_error(
"BilinearFormIntegrator::AssembleElementVector\n"
61 " is not implemented fot this class.");
65 void TransposeIntegrator::AssembleElementMatrix (
68 bfi -> AssembleElementMatrix (el, Trans, bfi_elmat);
73 void TransposeIntegrator::AssembleElementMatrix2 (
77 bfi -> AssembleElementMatrix2 (test_fe, trial_fe, Trans, bfi_elmat);
82 void TransposeIntegrator::AssembleFaceMatrix (
86 bfi -> AssembleFaceMatrix (el1, el2, Trans, bfi_elmat);
91 void LumpedIntegrator::AssembleElementMatrix (
94 bfi -> AssembleElementMatrix (el, Trans, elmat);
98 void InverseIntegrator::AssembleElementMatrix(
101 integrator->AssembleElementMatrix(el, Trans, elmat);
105 void SumIntegrator::AssembleElementMatrix(
108 MFEM_ASSERT(integrators.Size() > 0,
"empty SumIntegrator.");
110 integrators[0]->AssembleElementMatrix(el, Trans, elmat);
111 for (
int i = 1; i < integrators.Size(); i++)
113 integrators[i]->AssembleElementMatrix(el, Trans, elem_mat);
118 SumIntegrator::~SumIntegrator()
122 for (
int i = 0; i < integrators.Size(); i++)
123 delete integrators[i];
128 void DiffusionIntegrator::AssembleElementMatrix
135 bool square = (dim == spaceDim);
138 #ifdef MFEM_THREAD_SAFE
139 DenseMatrix dshape(nd,dim), dshapedxt(nd,spaceDim), invdfdx(dim,spaceDim);
141 dshape.SetSize(nd,dim);
142 dshapedxt.SetSize(nd,spaceDim);
151 if (el.
Space() == FunctionSpace::Pk)
157 if (el.
Space() == FunctionSpace::rQk)
174 w = ip.
weight / (square ? w : w*w*w);
175 Mult(dshape, invdfdx, dshapedxt);
179 w *= Q->Eval(Trans, ip);
184 MQ->Eval(invdfdx, Trans, ip);
186 Mult(dshapedxt, invdfdx, dshape);
192 void DiffusionIntegrator::AssembleElementMatrix2(
196 int tr_nd = trial_fe.
GetDof();
197 int te_nd = test_fe.
GetDof();
198 int dim = trial_fe.
GetDim();
200 bool square = (dim == spaceDim);
203 #ifdef MFEM_THREAD_SAFE
204 DenseMatrix dshape(tr_nd, dim), dshapedxt(tr_nd, spaceDim);
205 DenseMatrix te_dshape(te_nd, dim), te_dshapedxt(te_nd, spaceDim);
208 dshape.SetSize(tr_nd, dim);
209 dshapedxt.
SetSize(tr_nd, spaceDim);
210 te_dshape.SetSize(te_nd, dim);
211 te_dshapedxt.
SetSize(te_nd, spaceDim);
212 invdfdx.
SetSize(dim, spaceDim);
220 if (trial_fe.
Space() == FunctionSpace::Pk)
225 if (trial_fe.
Space() == FunctionSpace::rQk)
241 w = ip.
weight / (square ? w : w*w*w);
242 Mult(dshape, invdfdx, dshapedxt);
243 Mult(te_dshape, invdfdx, te_dshapedxt);
248 w *= Q->Eval(Trans, ip);
254 MQ->Eval(invdfdx, Trans, ip);
256 Mult(te_dshapedxt, invdfdx, te_dshape);
262 void DiffusionIntegrator::AssembleElementVector(
270 #ifdef MFEM_THREAD_SAFE
273 dshape.SetSize(nd,dim);
274 invdfdx.SetSize(dim);
278 pointflux.SetSize(dim);
286 if (el.
Space() == FunctionSpace::Pk)
292 if (el.
Space() == FunctionSpace::rQk)
310 dshape.MultTranspose(elfun, vec);
311 invdfdx.MultTranspose(vec, pointflux);
313 w *= Q->Eval(Tr, ip);
318 dshape.MultTranspose(elfun, pointflux);
319 invdfdx.MultTranspose(pointflux, vec);
320 MQ->Eval(mq, Tr, ip);
321 mq.
Mult(vec, pointflux);
324 invdfdx.Mult(pointflux, vec);
325 dshape.AddMult(vec, elvect);
329 void DiffusionIntegrator::ComputeElementFlux
333 int i, j, nd, dim, spaceDim, fnd;
339 #ifdef MFEM_THREAD_SAFE
340 DenseMatrix dshape(nd,dim), invdfdx(dim, spaceDim);
342 dshape.SetSize(nd,dim);
343 invdfdx.
SetSize(dim, spaceDim);
346 pointflux.SetSize(spaceDim);
352 for (i = 0; i < fnd; i++)
356 dshape.MultTranspose(u, vec);
364 for (j = 0; j < dim; j++)
365 flux(fnd*j+i) = pointflux(j);
370 pointflux *= Q->Eval(Trans,ip);
371 for (j = 0; j < dim; j++)
372 flux(fnd*j+i) = pointflux(j);
376 MQ->Eval(invdfdx, Trans, ip);
377 invdfdx.
Mult(pointflux, vec);
378 for (j = 0; j < dim; j++)
379 flux(fnd*j+i) = vec(j);
384 double DiffusionIntegrator::ComputeFluxEnergy
388 int i, j, k, nd, dim, order;
394 #ifdef MFEM_THREAD_SAFE
399 pointflux.SetSize(dim);
416 for (k = 0; k < dim; k++)
417 for (j = 0; j < nd; j++)
418 pointflux(k) += flux(k*nd+j)*shape(j);
425 co *= ( pointflux * pointflux );
427 co *= Q->Eval(Trans, ip);
431 MQ->Eval(invdfdx, Trans, ip);
442 void MassIntegrator::AssembleElementMatrix
459 if (el.
Space() == FunctionSpace::rQk)
474 w *= Q -> Eval(Trans, ip);
480 void MassIntegrator::AssembleElementMatrix2(
484 int tr_nd = trial_fe.
GetDof();
485 int te_nd = test_fe.
GetDof();
490 shape.SetSize (tr_nd);
491 te_shape.SetSize (te_nd);
511 w *= Q -> Eval(Trans, ip);
519 void ConvectionIntegrator::AssembleElementMatrix(
526 dshape.SetSize(nd,dim);
541 Q.Eval(Q_ir, Trans, *ir);
552 Q_ir.GetColumnReference(i, vec1);
553 vec1 *= alpha * ip.
weight;
555 adjJ.Mult(vec1, vec2);
556 dshape.Mult(vec2, BdFidxT);
563 void GroupConvectionIntegrator::AssembleElementMatrix(
570 dshape.SetSize(nd,dim);
573 grad.SetSize(nd,dim);
582 Q.Eval(Q_nodal, Trans, el.
GetNodes());
594 Mult(dshape, adjJ, grad);
596 double w = alpha * ip.
weight;
599 for (
int k = 0; k < nd; k++)
601 double wsk = w*shape(k);
602 for (
int l = 0; l < nd; l++)
605 for (
int s = 0; s < dim; s++)
606 a += Q_nodal(s,k)*grad(l,s);
614 void VectorMassIntegrator::AssembleElementMatrix
625 vdim = (VQ) ? (VQ -> GetVDim()) : ((MQ) ? (MQ -> GetVDim()) : (dim));
629 partelmat.SetSize(nd);
633 mcoeff.SetSize(vdim);
640 if (el.
Space() == FunctionSpace::rQk)
659 VQ->Eval(vec, Trans, ip);
660 for (
int k = 0; k < vdim; k++)
661 elmat.
AddMatrix(norm*vec(k), partelmat, nd*k, nd*k);
665 MQ->Eval(mcoeff, Trans, ip);
666 for (
int i = 0; i < vdim; i++)
667 for (
int j = 0; j < vdim; j++)
668 elmat.
AddMatrix(norm*mcoeff(i,j), partelmat, nd*i, nd*j);
673 norm *= Q->Eval(Trans, ip);
675 for (
int k = 0; k < vdim; k++)
681 void VectorMassIntegrator::AssembleElementMatrix2(
685 int tr_nd = trial_fe.
GetDof();
686 int te_nd = test_fe.
GetDof();
687 int dim = trial_fe.
GetDim();
693 vdim = (VQ) ? (VQ -> GetVDim()) : ((MQ) ? (MQ -> GetVDim()) : (dim));
695 elmat.
SetSize(te_nd*vdim, tr_nd*vdim);
696 shape.SetSize(tr_nd);
697 te_shape.SetSize(te_nd);
698 partelmat.SetSize(te_nd, tr_nd);
702 mcoeff.SetSize(vdim);
708 Trans.
OrderW() + Q_order);
710 if (trial_fe.
Space() == FunctionSpace::rQk)
726 MultVWt(te_shape, shape, partelmat);
730 VQ->Eval(vec, Trans, ip);
731 for (
int k = 0; k < vdim; k++)
732 elmat.
AddMatrix(norm*vec(k), partelmat, te_nd*k, tr_nd*k);
736 MQ->Eval(mcoeff, Trans, ip);
737 for (
int i = 0; i < vdim; i++)
738 for (
int j = 0; j < vdim; j++)
739 elmat.
AddMatrix(norm*mcoeff(i,j), partelmat, te_nd*i, tr_nd*j);
744 norm *= Q->Eval(Trans, ip);
746 for (
int k = 0; k < vdim; k++)
747 elmat.
AddMatrix(partelmat, te_nd*k, tr_nd*k);
752 void VectorFEDivergenceIntegrator::AssembleElementMatrix2(
756 int trial_nd = trial_fe.
GetDof(), test_nd = test_fe.
GetDof(), i;
758 #ifdef MFEM_THREAD_SAFE
759 Vector divshape(trial_nd), shape(test_nd);
761 divshape.SetSize(trial_nd);
765 elmat.
SetSize(test_nd, trial_nd);
784 w *= Q->Eval(Trans, ip);
791 void VectorFECurlIntegrator::AssembleElementMatrix2(
795 int trial_nd = trial_fe.
GetDof(), test_nd = test_fe.
GetDof(), i;
796 int dim = trial_fe.
GetDim();
798 #ifdef MFEM_THREAD_SAFE
803 curlshapeTrial.
SetSize(trial_nd, dim);
804 curlshapeTrial_dFT.
SetSize(trial_nd, dim);
805 vshapeTest.
SetSize(test_nd, dim);
808 elmat.
SetSize(test_nd, trial_nd);
827 w *= Q->Eval(Trans, ip);
829 AddMultABt(vshapeTest, curlshapeTrial_dFT, elmat);
833 void DerivativeIntegrator::AssembleElementMatrix2 (
839 int dim = trial_fe.
GetDim();
840 int trial_nd = trial_fe.
GetDof();
841 int test_nd = test_fe.
GetDof();
846 elmat.
SetSize (test_nd,trial_nd);
847 dshape.SetSize (trial_nd,dim);
848 dshapedxt.SetSize(trial_nd,dim);
849 dshapedxi.SetSize(trial_nd);
850 invdfdx.SetSize(dim);
851 shape.SetSize (test_nd);
857 if (trial_fe.
Space() == FunctionSpace::Pk)
862 if (trial_fe.
Space() == FunctionSpace::rQk)
877 Mult (dshape, invdfdx, dshapedxt);
881 for (l = 0; l < trial_nd; l++)
882 dshapedxi(l) = dshapedxt(l,xi);
884 shape *= Q.Eval(Trans,ip) * det * ip.
weight;
889 void CurlCurlIntegrator::AssembleElementMatrix
897 #ifdef MFEM_THREAD_SAFE
898 DenseMatrix Curlshape(nd,dim), Curlshape_dFt(nd,dim);
900 Curlshape.SetSize(nd,dim);
909 if (el.
Space() == FunctionSpace::Pk)
930 w *= Q->Eval(Trans, ip);
937 void VectorCurlCurlIntegrator::AssembleElementMatrix(
942 int cld = (dim*(dim-1))/2;
944 #ifdef MFEM_THREAD_SAFE
945 DenseMatrix dshape_hat(dof, dim), dshape(dof, dim);
948 dshape_hat.SetSize(dof, dim);
950 curlshape.SetSize(dim*dof, cld);
972 Mult(dshape_hat, Jadj, dshape);
976 w *= Q->Eval(Trans, ip);
982 double VectorCurlCurlIntegrator::GetElementEnergy(
988 #ifdef MFEM_THREAD_SAFE
989 DenseMatrix dshape_hat(dof, dim), Jadj(dim), grad_hat(dim), grad(dim);
991 dshape_hat.SetSize(dof, dim);
993 grad_hat.SetSize(dim);
1007 for (
int i = 0; i < ir->GetNPoints(); i++)
1012 MultAtB(elfun_mat, dshape_hat, grad_hat);
1018 Mult(grad_hat, Jadj, grad);
1022 double curl = grad(0,1) - grad(1,0);
1027 double curl_x = grad(2,1) - grad(1,2);
1028 double curl_y = grad(0,2) - grad(2,0);
1029 double curl_z = grad(1,0) - grad(0,1);
1030 w *= curl_x * curl_x + curl_y * curl_y + curl_z * curl_z;
1034 w *= Q->Eval(Tr, ip);
1039 elfun_mat.ClearExternalData();
1041 return 0.5 * energy;
1045 void VectorFEMassIntegrator::AssembleElementMatrix(
1055 #ifdef MFEM_THREAD_SAFE
1056 Vector D(VQ ? VQ->GetVDim() : 0);
1058 DenseMatrix K(MQ ? MQ->GetVDim() : 0, MQ ? MQ->GetVDim() : 0);
1061 D.SetSize(VQ ? VQ->GetVDim() : 0);
1062 K.SetSize(MQ ? MQ->GetVDim() : 0, MQ ? MQ->GetVDim() : 0);
1088 MQ->Eval(K, Trans, ip);
1095 VQ->Eval(D, Trans, ip);
1102 w *= Q -> Eval (Trans, ip);
1108 void VectorFEMassIntegrator::AssembleElementMatrix2(
1113 int dim = test_fe.
GetDim();
1114 int trial_dof = trial_fe.
GetDof();
1115 int test_dof = test_fe.
GetDof();
1119 mfem_error(
"VectorFEMassIntegrator::AssembleElementMatrix2(...)\n"
1120 " is not implemented for vector/tensor permeability");
1122 #ifdef MFEM_THREAD_SAFE
1126 vshape.
SetSize(trial_dof, dim);
1130 elmat.
SetSize (dim*test_dof, trial_dof);
1151 w *= Q -> Eval (Trans, ip);
1153 for (
int d = 0; d < dim; d++)
1155 for (
int j = 0; j < test_dof; j++)
1157 for (
int k = 0; k < trial_dof; k++)
1159 elmat(d * test_dof + j, k) += w * shape(j) * vshape(k, d);
1167 void VectorDivergenceIntegrator::AssembleElementMatrix2(
1173 int dim = trial_fe.
GetDim();
1174 int trial_dof = trial_fe.
GetDof();
1175 int test_dof = test_fe.
GetDof();
1178 dshape.SetSize (trial_dof, dim);
1179 gshape.SetSize (trial_dof, dim);
1181 divshape.SetSize (dim*trial_dof);
1182 shape.SetSize (test_dof);
1184 elmat.
SetSize (test_dof, dim*trial_dof);
1195 for (
int i = 0; i < ir -> GetNPoints(); i++)
1205 Mult (dshape, Jadj, gshape);
1207 gshape.GradToDiv (divshape);
1211 c *= Q -> Eval (Trans, ip);
1220 void DivDivIntegrator::AssembleElementMatrix(
1228 #ifdef MFEM_THREAD_SAFE
1244 for (
int i = 0; i < ir -> GetNPoints(); i++)
1254 c *= Q -> Eval (Trans, ip);
1262 void VectorDiffusionIntegrator::AssembleElementMatrix(
1274 Jinv. SetSize (dim);
1275 dshape.SetSize (dof, dim);
1276 gshape.SetSize (dof, dim);
1277 pelmat.SetSize (dof);
1284 if (el.
Space() == FunctionSpace::rQk)
1292 for (
int i = 0; i < ir -> GetNPoints(); i++)
1302 Mult (dshape, Jinv, gshape);
1307 norm *= Q -> Eval (Trans, ip);
1311 for (
int d = 0; d < dim; d++)
1313 for (
int k = 0; k < dof; k++)
1314 for (
int l = 0; l < dof; l++)
1315 elmat (dof*d+k, dof*d+l) += pelmat (k, l);
1321 void ElasticityIntegrator::AssembleElementMatrix(
1328 #ifdef MFEM_THREAD_SAFE
1329 DenseMatrix dshape(dof, dim), Jinv(dim), gshape(dof, dim), pelmat(dof);
1330 Vector divshape(dim*dof);
1333 dshape.SetSize(dof, dim);
1334 gshape.SetSize(dof, dim);
1350 for (
int i = 0; i < ir -> GetNPoints(); i++)
1359 Mult(dshape, Jinv, gshape);
1361 gshape.GradToDiv (divshape);
1363 M = mu->Eval(Trans, ip);
1365 L = lambda->Eval(Trans, ip);
1377 for (
int d = 0; d < dim; d++)
1379 for (
int k = 0; k < dof; k++)
1380 for (
int l = 0; l < dof; l++)
1381 elmat (dof*d+k, dof*d+l) += (M * w) * pelmat(k, l);
1383 for (
int i = 0; i < dim; i++)
1384 for (
int j = 0; j < dim; j++)
1386 for (
int k = 0; k < dof; k++)
1387 for (
int l = 0; l < dof; l++)
1388 elmat(dof*i+k, dof*j+l) +=
1389 (M * w) * gshape(k, j) * gshape(l, i);
1401 int dim, ndof1, ndof2;
1407 Vector vu(dim), nor(dim);
1414 shape1.SetSize(ndof1);
1415 shape2.SetSize(ndof2);
1429 if (el1.
Space() == FunctionSpace::Pk)
1446 u->Eval(vu, *Trans.
Elem1, eip1);
1449 nor(0) = 2*eip1.
x - 1.0;
1454 a = 0.5 * alpha * un;
1455 b = beta * fabs(un);
1463 if (un >= 0.0 && ndof2)
1466 rho_p = rho->Eval(*Trans.
Elem2, eip2);
1470 rho_p = rho->Eval(*Trans.
Elem1, eip1);
1479 for (
int i = 0; i < ndof1; i++)
1480 for (
int j = 0; j < ndof1; j++)
1481 elmat(i, j) += w * shape1(i) * shape1(j);
1489 for (
int i = 0; i < ndof2; i++)
1490 for (
int j = 0; j < ndof1; j++)
1491 elmat(ndof1+i, j) -= w * shape2(i) * shape1(j);
1496 for (
int i = 0; i < ndof2; i++)
1497 for (
int j = 0; j < ndof2; j++)
1498 elmat(ndof1+i, ndof1+j) += w * shape2(i) * shape2(j);
1500 for (
int i = 0; i < ndof1; i++)
1501 for (
int j = 0; j < ndof2; j++)
1502 elmat(i, ndof1+j) -= w * shape1(i) * shape2(j);
1508 void DGDiffusionIntegrator::AssembleFaceMatrix(
1512 int dim, ndof1, ndof2, ndofs;
1513 bool kappa_is_nonzero = (
kappa != 0.);
1526 shape1.SetSize(ndof1);
1527 dshape1.SetSize(ndof1, dim);
1528 dshape1dn.SetSize(ndof1);
1532 shape2.SetSize(ndof2);
1533 dshape2.SetSize(ndof2, dim);
1534 dshape2dn.SetSize(ndof2);
1539 ndofs = ndof1 + ndof2;
1542 if (kappa_is_nonzero)
1570 nor(0) = 2*eip1.
x - 1.0;
1583 w *= Q->Eval(*Trans.
Elem1, eip1);
1589 MQ->Eval(mq, *Trans.
Elem1, eip1);
1590 mq.MultTranspose(nh, ni);
1594 if (kappa_is_nonzero)
1607 dshape1.Mult(nh, dshape1dn);
1608 for (
int i = 0; i < ndof1; i++)
1609 for (
int j = 0; j < ndof1; j++)
1610 elmat(i, j) += shape1(i) * dshape1dn(j);
1622 w *= Q->Eval(*Trans.
Elem2, eip2);
1628 MQ->Eval(mq, *Trans.
Elem2, eip2);
1629 mq.MultTranspose(nh, ni);
1633 if (kappa_is_nonzero)
1636 dshape2.Mult(nh, dshape2dn);
1638 for (
int i = 0; i < ndof1; i++)
1639 for (
int j = 0; j < ndof2; j++)
1640 elmat(i, ndof1 + j) += shape1(i) * dshape2dn(j);
1642 for (
int i = 0; i < ndof2; i++)
1643 for (
int j = 0; j < ndof1; j++)
1644 elmat(ndof1 + i, j) -= shape2(i) * dshape1dn(j);
1646 for (
int i = 0; i < ndof2; i++)
1647 for (
int j = 0; j < ndof2; j++)
1648 elmat(ndof1 + i, ndof1 + j) -= shape2(i) * dshape2dn(j);
1651 if (kappa_is_nonzero)
1655 for (
int i = 0; i < ndof1; i++)
1657 const double wsi = wq*shape1(i);
1658 for (
int j = 0; j <= i; j++)
1659 jmat(i, j) += wsi * shape1(j);
1663 for (
int i = 0; i < ndof2; i++)
1665 const int i2 = ndof1 + i;
1666 const double wsi = wq*shape2(i);
1667 for (
int j = 0; j < ndof1; j++)
1668 jmat(i2, j) -= wsi * shape1(j);
1669 for (
int j = 0; j <= i; j++)
1670 jmat(i2, ndof1 + j) += wsi * shape2(j);
1677 if (kappa_is_nonzero)
1679 for (
int i = 0; i < ndofs; i++)
1681 for (
int j = 0; j < i; j++)
1683 double aij = elmat(i,j), aji = elmat(j,i), mij = jmat(i,j);
1684 elmat(i,j) = sigma*aji - aij + mij;
1685 elmat(j,i) = sigma*aij - aji + mij;
1687 elmat(i,i) = (sigma - 1.)*elmat(i,i) + jmat(i,i);
1692 for (
int i = 0; i < ndofs; i++)
1694 for (
int j = 0; j < i; j++)
1696 double aij = elmat(i,j), aji = elmat(j,i);
1697 elmat(i,j) = sigma*aji - aij;
1698 elmat(j,i) = sigma*aij - aji;
1700 elmat(i,i) *= (sigma - 1.);
1705 void TraceJumpIntegrator::AssembleFaceMatrix(
1710 int i, j, face_ndof, ndof1, ndof2;
1715 face_ndof = trial_face_fe.
GetDof();
1716 ndof1 = test_fe1.
GetDof();
1718 face_shape.SetSize(face_ndof);
1719 shape1.SetSize(ndof1);
1723 ndof2 = test_fe2.
GetDof();
1724 shape2.SetSize(ndof2);
1729 elmat.
SetSize(ndof1 + ndof2, face_ndof);
1740 if (trial_face_fe.
GetMapType() == FiniteElement::VALUE)
1751 trial_face_fe.
CalcShape(ip, face_shape);
1764 if (trial_face_fe.
GetMapType() == FiniteElement::VALUE)
1767 for (i = 0; i < ndof1; i++)
1768 for (j = 0; j < face_ndof; j++)
1769 elmat(i, j) += shape1(i) * face_shape(j);
1773 for (i = 0; i < ndof2; i++)
1774 for (j = 0; j < face_ndof; j++)
1775 elmat(ndof1+i, j) -= shape2(i) * face_shape(j);
int GetNPoints() const
Returns the number of the points in the integration rule.
Abstract class for Finite Elements.
void MultABt(const DenseMatrix &A, const DenseMatrix &B, DenseMatrix &ABt)
Multiply a matrix A with the transpose of a matrix B: A*Bt.
int GetDim() const
Returns the space dimension for the finite element.
void AddMultVWt(const Vector &v, const Vector &w, DenseMatrix &VWt)
VWt += v w^t.
Class for integration rule.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
void MultVWt(const Vector &v, const Vector &w, DenseMatrix &VWt)
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
void SetSize(int s)
Resizes the vector if the new size is different.
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
double InnerProduct(const double *x, const double *y) const
Compute y^t A x.
void CalcAdjugate(const DenseMatrix &a, DenseMatrix &adja)
int GetOrder() const
Returns the order of the finite element.
int Space() const
Returns the type of space on each element.
void CalcOrtho(const DenseMatrix &J, Vector &n)
friend void Mult(const DenseMatrix &b, const DenseMatrix &c, DenseMatrix &a)
Matrix matrix multiplication. A = B * C.
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
void MultTranspose(const double *x, double *y) const
Multiply a vector with the transpose matrix.
IntegrationRules IntRules(0)
A global object with all integration rules (defined in intrules.cpp)
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const =0
void Invert()
Replaces the current matrix with its inverse.
const IntegrationRule & GetNodes() const
void AddMult_a_VVt(const double a, const Vector &v, DenseMatrix &VVt)
VVt += a * v v^t.
void CalcInverse(const DenseMatrix &a, DenseMatrix &inva)
int GetGeomType() const
Returns the geometry type:
void Transpose()
(*this) = (*this)^t
void MultVVt(const Vector &v, DenseMatrix &vvt)
Make a matrix from a vector V.Vt.
void AddMultABt(const DenseMatrix &A, const DenseMatrix &B, DenseMatrix &ABt)
ABt += A * B^t.
int GetDof() const
Returns the degrees of freedom in the FE space.
void mfem_error(const char *msg)
void MultAAt(const DenseMatrix &a, DenseMatrix &aat)
Calculate the matrix A.At.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
void AddMatrix(DenseMatrix &A, int ro, int co)
Perform (ro+i,co+j)+=A(i,j) for 0<=i<A.Height, 0<=j<A.Width.
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Class for integration point with weight.
void GradToCurl(DenseMatrix &curl)
void MultAtB(const DenseMatrix &A, const DenseMatrix &B, DenseMatrix &AtB)
Multiply the transpose of a matrix A with a matrix B: At*B.
IntegrationRules RefinedIntRules(1)
A global object with all refined integration rules.
void AddMultADAt(const DenseMatrix &A, const Vector &D, DenseMatrix &ADAt)
ADAt += A D A^t, where D is diagonal.
void SetSize(int s)
If the matrix is not a square matrix of size s then recreate it.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const =0
void AddMult_a_AAt(double a, const DenseMatrix &A, DenseMatrix &AAt)
AAt += a * A * A^t.