34 #include <boost/core/ignore_unused.hpp>
47 namespace LocalRegions
53 Ba.GetNumModes(), (Bb.GetNumModes())),
56 Ba.GetNumModes(), (Bb.GetNumModes())),
60 std::bind(&
Expansion2D::CreateMatrix, this, std::placeholders::_1),
61 std::string(
"TriExpMatrix")),
62 m_staticCondMatrixManager(std::bind(&
Expansion::CreateStaticCondMatrix,
63 this, std::placeholders::_1),
64 std::string(
"TriExpStaticCondMatrix"))
69 : StdExpansion(T), StdExpansion2D(T), StdTriExp(T),
Expansion(T),
71 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
81 int nquad0 =
m_base[0]->GetNumPoints();
82 int nquad1 =
m_base[1]->GetNumPoints();
90 Vmath::Vmul(nquad0 * nquad1, jac, 1, inarray, 1, tmp, 1);
94 Vmath::Smul(nquad0 * nquad1, jac[0], inarray, 1, tmp, 1);
98 ival = StdTriExp::v_Integral(tmp);
107 int nquad0 =
m_base[0]->GetNumPoints();
108 int nquad1 =
m_base[1]->GetNumPoints();
109 int nqtot = nquad0 * nquad1;
116 StdTriExp::v_PhysDeriv(inarray, diff0, diff1);
123 Vmath::Vvtvp(nqtot, df[1], 1, diff1, 1, out_d0, 1, out_d0, 1);
129 Vmath::Vvtvp(nqtot, df[3], 1, diff1, 1, out_d1, 1, out_d1, 1);
135 Vmath::Vvtvp(nqtot, df[5], 1, diff1, 1, out_d2, 1, out_d2, 1);
186 ASSERTL1(
false,
"input dir is out of range");
201 int nquad0 =
m_base[0]->GetNumPoints();
202 int nquad1 =
m_base[1]->GetNumPoints();
203 int nqtot = nquad0 * nquad1;
212 StdTriExp::v_PhysDeriv(inarray, diff0, diff1);
220 for (
int i = 0; i < 2; ++i)
223 for (
int k = 0; k < (
m_geom->GetCoordim()); ++k)
225 Vmath::Vvtvp(nqtot, &df[2 * k + i][0], 1, &direction[k * nqtot],
226 1, &tangmat[i][0], 1, &tangmat[i][0], 1);
231 Vmath::Vmul(nqtot, &tangmat[0][0], 1, &diff0[0], 1, &out[0], 1);
232 Vmath::Vvtvp(nqtot, &tangmat[1][0], 1, &diff1[0], 1, &out[0], 1,
239 for (
int i = 0; i < 2; ++i)
242 for (
int k = 0; k < (
m_geom->GetCoordim()); ++k)
244 Vmath::Svtvp(nqtot, df[2 * k + i][0], &direction[k * nqtot], 1,
245 &tangmat[i][0], 1, &tangmat[i][0], 1);
250 Vmath::Vmul(nqtot, &tangmat[0][0], 1, &diff0[0], 1, &out[0], 1);
252 Vmath::Vvtvp(nqtot, &tangmat[1][0], 1, &diff1[0], 1, &out[0], 1,
270 out = (*matsys) * in;
278 int npoints[2] = {
m_base[0]->GetNumPoints(),
m_base[1]->GetNumPoints()};
279 int nmodes[2] = {
m_base[0]->GetNumModes(),
m_base[1]->GetNumModes()};
281 fill(outarray.get(), outarray.get() +
m_ncoeffs, 0.0);
283 if (nmodes[0] == 1 && nmodes[1] == 1)
285 outarray[0] = inarray[0];
291 for (i = 0; i < 3; i++)
299 for (i = 0; i < npoints[0]; i++)
301 physEdge[0][i] = inarray[i];
305 for (i = 0; i < npoints[1]; i++)
307 physEdge[1][i] = inarray[npoints[0] - 1 + i * npoints[0]];
308 physEdge[2][i] = inarray[i * npoints[0]];
317 for (i = 1; i < 3; i++)
325 for (i = 1; i < 3; i++)
331 m_base[0]->GetPointsKey(), physEdge[i]);
333 npoints[1] = npoints[0];
343 for (i = 0; i < 3; i++)
345 segexp[i]->FwdTransBndConstrained(physEdge[i], coeffEdge[i]);
350 for (j = 0; j < nmodes[i != 0]; j++)
353 outarray[mapArray[j]] =
sign * coeffEdge[i][j];
358 int nInteriorDofs =
m_ncoeffs - nBoundaryDofs;
360 if (nInteriorDofs > 0)
384 for (i = 0; i < nInteriorDofs; i++)
386 rhs[i] = tmp1[mapArray[i]];
389 Blas::Dgemv(
'N', nInteriorDofs, nInteriorDofs, matsys->Scale(),
390 &((matsys->GetOwnedMatrix())->GetPtr())[0], nInteriorDofs,
391 rhs.get(), 1, 0.0, result.get(), 1);
393 for (i = 0; i < nInteriorDofs; i++)
395 outarray[mapArray[i]] = result[i];
417 int nquad0 =
m_base[0]->GetNumPoints();
418 int nquad1 =
m_base[1]->GetNumPoints();
419 int order0 =
m_base[0]->GetNumModes();
421 if (multiplybyweights)
428 m_base[1]->GetBdata(), tmp, outarray, wsp);
435 m_base[1]->GetBdata(), inarray, outarray,
449 (iprodmat->GetOwnedMatrix())->GetPtr().get(),
m_ncoeffs,
450 inarray.get(), 1, 0.0, outarray.get(), 1);
457 int nquad0 =
m_base[0]->GetNumPoints();
458 int nquad1 =
m_base[1]->GetNumPoints();
459 int nqtot = nquad0 * nquad1;
460 int nmodes0 =
m_base[0]->GetNumModes();
461 int wspsize = max(max(nqtot,
m_ncoeffs), nquad1 * nmodes0);
480 tmp2, outarray, tmp0);
488 ASSERTL1((dir == 0) || (dir == 1) || (dir == 2),
"Invalid direction.");
490 "Invalid direction.");
492 int nquad0 =
m_base[0]->GetNumPoints();
493 int nquad1 =
m_base[1]->GetNumPoints();
494 int nqtot = nquad0 * nquad1;
495 int nmodes0 =
m_base[0]->GetNumModes();
496 int wspsize = max(max(nqtot,
m_ncoeffs), nquad1 * nmodes0);
513 for (
int i = 0; i < nquad1; ++i)
515 gfac0[i] = 2.0 / (1 - z1[i]);
517 for (
int i = 0; i < nquad0; ++i)
519 gfac1[i] = 0.5 * (1 + z0[i]);
522 for (
int i = 0; i < nquad1; ++i)
524 Vmath::Smul(nquad0, gfac0[i], &inarray[0] + i * nquad0, 1,
525 &tmp0[0] + i * nquad0, 1);
528 for (
int i = 0; i < nquad1; ++i)
530 Vmath::Vmul(nquad0, &gfac1[0], 1, &tmp0[0] + i * nquad0, 1,
531 &tmp1[0] + i * nquad0, 1);
536 Vmath::Vmul(nqtot, &df[2 * dir][0], 1, &tmp0[0], 1, &tmp0[0], 1);
537 Vmath::Vmul(nqtot, &df[2 * dir + 1][0], 1, &tmp1[0], 1, &tmp1[0], 1);
538 Vmath::Vmul(nqtot, &df[2 * dir + 1][0], 1, &inarray[0], 1, &tmp2[0], 1);
542 Vmath::Smul(nqtot, df[2 * dir][0], tmp0, 1, tmp0, 1);
543 Vmath::Smul(nqtot, df[2 * dir + 1][0], tmp1, 1, tmp1, 1);
544 Vmath::Smul(nqtot, df[2 * dir + 1][0], inarray, 1, tmp2, 1);
575 ASSERTL1(
false,
"input dir is out of range");
584 (iprodmat->GetOwnedMatrix())->GetPtr().get(),
m_ncoeffs,
585 inarray.get(), 1, 0.0, outarray.get(), 1);
607 int nquad0 =
m_base[0]->GetNumPoints();
608 int nquad1 =
m_base[1]->GetNumPoints();
609 int nqtot = nquad0 * nquad1;
610 int nmodes0 =
m_base[0]->GetNumModes();
611 int wspsize = max(max(nqtot,
m_ncoeffs), nquad1 * nmodes0);
627 for (i = 0; i < nquad1; ++i)
629 gfac0[i] = 2.0 / (1 - z1[i]);
631 for (i = 0; i < nquad0; ++i)
633 gfac1[i] = 0.5 * (1 + z0[i]);
635 for (i = 0; i < nquad1; ++i)
637 Vmath::Smul(nquad0, gfac0[i], &inarray[0] + i * nquad0, 1,
638 &tmp0[0] + i * nquad0, 1);
640 for (i = 0; i < nquad1; ++i)
642 Vmath::Vmul(nquad0, &gfac1[0], 1, &tmp0[0] + i * nquad0, 1,
643 &tmp1[0] + i * nquad0, 1);
650 Vmath::Vmul(nqtot, &dfdir[0][0], 1, &tmp0[0], 1, &tmp0[0], 1);
651 Vmath::Vmul(nqtot, &dfdir[1][0], 1, &tmp1[0], 1, &tmp1[0], 1);
652 Vmath::Vmul(nqtot, &dfdir[1][0], 1, &inarray[0], 1, &tmp2[0], 1);
654 Vmath::Vadd(nqtot, &tmp0[0], 1, &tmp1[0], 1, &tmp1[0], 1);
662 tmp2, outarray, tmp0);
671 int nq =
m_base[0]->GetNumPoints() *
m_base[1]->GetNumPoints();
680 Vmath::Vvtvvtp(nq, &normals[0][0], 1, &Fx[0], 1, &normals[1][0], 1,
681 &Fy[0], 1, &Fn[0], 1);
682 Vmath::Vvtvp(nq, &normals[2][0], 1, &Fz[0], 1, &Fn[0], 1, &Fn[0], 1);
686 Vmath::Svtsvtp(nq, normals[0][0], &Fx[0], 1, normals[1][0], &Fy[0], 1,
688 Vmath::Svtvp(nq, normals[2][0], &Fz[0], 1, &Fn[0], 1, &Fn[0], 1);
711 m_base[0]->GetPointsKey());
713 m_base[1]->GetPointsKey());
724 ASSERTL1(Lcoords[0] >= -1.0 && Lcoords[1] <= 1.0 && Lcoords[1] >= -1.0 &&
726 "Local coordinates are not in region [-1,1]");
730 for (i = 0; i <
m_geom->GetCoordim(); ++i)
732 coords[i] =
m_geom->GetCoord(i, Lcoords);
753 return StdTriExp::v_PhysEvaluate(Lcoord, physvals);
762 m_geom->GetLocCoords(coord, Lcoord);
764 return StdTriExp::v_PhysEvaluate(Lcoord, physvals);
772 int nquad0 =
m_base[0]->GetNumPoints();
773 int nquad1 =
m_base[1]->GetNumPoints();
780 Vmath::Vcopy(nquad0, &(inarray[0]), 1, &(outarray[0]), 1);
784 Vmath::Vcopy(nquad1, &(inarray[0]) + (nquad0 - 1), nquad0,
789 Vmath::Vcopy(nquad1, &(inarray[0]), nquad0, &(outarray[0]), 1);
793 ASSERTL0(
false,
"edge value (< 3) is out of range");
797 ASSERTL1(EdgeExp->GetBasis(0)->GetPointsType() ==
799 "Edge expansion should be GLL");
802 if (
m_base[edge ? 1 : 0]->GetPointsKey() !=
803 EdgeExp->GetBasis(0)->GetPointsKey())
810 EdgeExp->GetBasis(0)->GetPointsKey(), outarray);
821 Vmath::Reverse(EdgeExp->GetNumPoints(0), &outarray[0], 1, &outarray[0],
830 boost::ignore_unused(edge, inarray, outarray);
831 ASSERTL0(
false,
"Routine not implemented for triangular elements");
837 boost::ignore_unused(edge, outarray);
838 ASSERTL0(
false,
"Routine not implemented for triangular elements");
843 int nquad0 =
m_base[0]->GetNumPoints();
844 int nquad1 =
m_base[1]->GetNumPoints();
851 for (
int i = 0; i < nquad0; ++i)
858 for (
int i = 0; i < nquad1; ++i)
860 outarray[i] = (nquad0 - 1) + i * nquad0;
865 for (
int i = 0; i < nquad1; ++i)
867 outarray[i] = i * nquad0;
871 ASSERTL0(
false,
"edge value (< 3) is out of range");
883 for (i = 0; i < ptsKeys.size(); ++i)
895 geomFactors->GetDerivFactors(ptsKeys);
897 int nqe =
m_base[0]->GetNumPoints();
902 for (i = 0; i < dim; ++i)
929 Vmath::Fill(nqe, df[2 * i + 1][0] + df[2 * i][0], normal[i],
940 ASSERTL0(
false,
"Edge is out of range (edge < 3)");
947 fac += normal[i][0] * normal[i][0];
949 fac = 1.0 /
sqrt(fac);
962 int nquad0 = ptsKeys[0].GetNumPoints();
963 int nquad1 = ptsKeys[1].GetNumPoints();
976 for (j = 0; j < nquad0; ++j)
981 normals[i * nquad0 + j] =
982 -df[2 * i + 1][j] * edgejac[j];
985 from_key = ptsKeys[0];
988 for (j = 0; j < nquad1; ++j)
990 edgejac[j] = jac[nquad0 * j + nquad0 - 1];
993 normals[i * nquad1 + j] =
994 (df[2 * i][nquad0 * j + nquad0 - 1] +
995 df[2 * i + 1][nquad0 * j + nquad0 - 1]) *
999 from_key = ptsKeys[1];
1002 for (j = 0; j < nquad1; ++j)
1004 edgejac[j] = jac[nquad0 * j];
1007 normals[i * nquad1 + j] =
1008 -df[2 * i][nquad0 * j] * edgejac[j];
1011 from_key = ptsKeys[1];
1014 ASSERTL0(
false,
"edge is out of range (edge < 3)");
1028 m_base[0]->GetPointsKey(), &normal[i][0]);
1029 Vmath::Vmul(nqe, work, 1, normal[i], 1, normal[i], 1);
1036 Vmath::Vvtvp(nqe, normal[i], 1, normal[i], 1, work, 1, work, 1);
1046 Vmath::Vmul(nqe, normal[i], 1, work, 1, normal[i], 1);
1064 return m_geom->GetCoordim();
1068 const NekDouble *data,
const std::vector<unsigned int> &nummodes,
1069 const int mode_offset,
NekDouble *coeffs,
1070 std::vector<LibUtilities::BasisType> &fromType)
1072 boost::ignore_unused(fromType);
1074 int data_order0 = nummodes[mode_offset];
1075 int fillorder0 = min(
m_base[0]->GetNumModes(), data_order0);
1076 int data_order1 = nummodes[mode_offset + 1];
1077 int order1 =
m_base[1]->GetNumModes();
1078 int fillorder1 = min(order1, data_order1);
1091 "Extraction routine not set up for this basis");
1094 for (i = 0; i < fillorder0; ++i)
1096 Vmath::Vcopy(fillorder1 - i, &data[cnt], 1, &coeffs[cnt1], 1);
1097 cnt += data_order1 - i;
1103 ASSERTL0(
false,
"basis is either not set up or not hierarchicial");
1114 ASSERTL1(dir >= 0 && dir <= 1,
"input dir is out of range");
1138 returnval = StdTriExp::v_GenMatrix(mkey);
1152 return tmp->GetStdMatrix(mkey);
1174 StdExpansion::MassMatrixOp_MatFree(inarray, outarray, mkey);
1189 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
1197 StdExpansion::WeakDerivMatrixOp_MatFree(i, inarray, outarray, mkey);
1204 StdExpansion::WeakDirectionalDerivMatrixOp_MatFree(inarray, outarray, mkey);
1211 StdExpansion::MassLevelCurvatureMatrixOp_MatFree(inarray, outarray, mkey);
1227 if (inarray.get() == outarray.get())
1233 (mat->GetOwnedMatrix())->GetPtr().get(),
m_ncoeffs,
1234 tmp.get(), 1, 0.0, outarray.get(), 1);
1239 (mat->GetOwnedMatrix())->GetPtr().get(),
m_ncoeffs,
1240 inarray.get(), 1, 0.0, outarray.get(), 1);
1253 int nquad0 =
m_base[0]->GetNumPoints();
1254 int nquad1 =
m_base[1]->GetNumPoints();
1255 int nqtot = nquad0 * nquad1;
1256 int nmodes0 =
m_base[0]->GetNumModes();
1257 int nmodes1 =
m_base[1]->GetNumModes();
1259 max(max(max(nqtot,
m_ncoeffs), nquad1 * nmodes0), nquad0 * nmodes1);
1261 ASSERTL1(wsp.size() >= 3 * wspsize,
"Workspace is of insufficient size.");
1279 StdExpansion2D::PhysTensorDeriv(inarray, wsp1, wsp2);
1285 Vmath::Vvtvvtp(nqtot, &metric00[0], 1, &wsp1[0], 1, &metric01[0], 1,
1286 &wsp2[0], 1, &wsp0[0], 1);
1287 Vmath::Vvtvvtp(nqtot, &metric01[0], 1, &wsp1[0], 1, &metric11[0], 1,
1288 &wsp2[0], 1, &wsp2[0], 1);
1310 const unsigned int dim = 2;
1319 for (i = 0; i < dim; ++i)
1321 for (j = i; j < dim; ++j)
1329 const unsigned int nquad0 =
m_base[0]->GetNumPoints();
1330 const unsigned int nquad1 =
m_base[1]->GetNumPoints();
1334 for (i = 0; i < nquad1; i++)
1336 Blas::Dscal(nquad0, 2.0 / (1 - z1[i]), &dEta_dXi[0][0] + i * nquad0, 1);
1337 Blas::Dscal(nquad0, 2.0 / (1 - z1[i]), &dEta_dXi[1][0] + i * nquad0, 1);
1339 for (i = 0; i < nquad0; i++)
1341 Blas::Dscal(nquad1, 0.5 * (1 + z0[i]), &dEta_dXi[1][0] + i, nquad0);
1348 Vmath::Smul(nqtot, df[0][0], &dEta_dXi[0][0], 1, &tmp[0], 1);
1349 Vmath::Svtvp(nqtot, df[1][0], &dEta_dXi[1][0], 1, &tmp[0], 1, &tmp[0],
1357 Vmath::Smul(nqtot, df[2][0], &dEta_dXi[0][0], 1, &tmp[0], 1);
1358 Vmath::Svtvp(nqtot, df[3][0], &dEta_dXi[1][0], 1, &tmp[0], 1, &tmp[0],
1370 Vmath::Smul(nqtot, df[4][0], &dEta_dXi[0][0], 1, &tmp[0], 1);
1371 Vmath::Svtvp(nqtot, df[5][0], &dEta_dXi[1][0], 1, &tmp[0], 1,
1382 NekDouble g2 = df[1][0] * df[1][0] + df[3][0] * df[3][0];
1385 g2 += df[5][0] * df[5][0];
1392 Vmath::Vmul(nqtot, &df[0][0], 1, &dEta_dXi[0][0], 1, &tmp[0], 1);
1393 Vmath::Vvtvp(nqtot, &df[1][0], 1, &dEta_dXi[1][0], 1, &tmp[0], 1,
1403 Vmath::Vmul(nqtot, &df[2][0], 1, &dEta_dXi[0][0], 1, &tmp[0], 1);
1404 Vmath::Vvtvp(nqtot, &df[3][0], 1, &dEta_dXi[1][0], 1, &tmp[0], 1,
1419 Vmath::Vmul(nqtot, &df[4][0], 1, &dEta_dXi[0][0], 1, &tmp[0], 1);
1420 Vmath::Vvtvp(nqtot, &df[5][0], 1, &dEta_dXi[1][0], 1, &tmp[0], 1,
1435 for (
unsigned int i = 0; i < dim; ++i)
1437 for (
unsigned int j = i; j < dim; ++j)
1460 int n_coeffs = inarray.size();
1461 int nquad0 =
m_base[0]->GetNumPoints();
1462 int nquad1 =
m_base[1]->GetNumPoints();
1463 int nqtot = nquad0 * nquad1;
1464 int nmodes0 =
m_base[0]->GetNumModes();
1465 int nmodes1 =
m_base[1]->GetNumModes();
1466 int numMin2 = nmodes0, i;
1495 m_TriExp->BwdTrans(inarray, phys_tmp);
1496 m_OrthoTriExp->FwdTrans(phys_tmp, coeff);
1498 for (i = 0; i < n_coeffs; i++)
1503 numMin += numMin2 - 1;
1508 m_OrthoTriExp->BwdTrans(coeff, phys_tmp);
1509 m_TriExp->FwdTrans(phys_tmp, outarray);
1533 StdTriExp::v_SVVLaplacianFilter(array, mkey);
1549 boost::ignore_unused(d2factors);
1556 if (d0factors.size() != 3)
1562 if (d0factors[0].size() != nquad0)
1568 if (d0factors[1].size() != nquad1)
1584 int ncoords = normal_0.size();
1590 for (
int i = 0; i < nquad0; ++i)
1592 d1factors[0][i] = df[1][i] * normal_0[0][i];
1596 for (
int i = 0; i < nquad1; ++i)
1598 d0factors[1][i] = df[0][(i + 1) * nquad0 - 1] * normal_1[0][i];
1599 d0factors[2][i] = df[0][i * nquad0] * normal_2[0][i];
1602 for (
int n = 1; n < ncoords; ++n)
1606 for (
int i = 0; i < nquad0; ++i)
1608 d1factors[0][i] += df[2 * n + 1][i] * normal_0[n][i];
1613 for (
int i = 0; i < nquad1; ++i)
1616 df[2 * n][(i + 1) * nquad0 - 1] * normal_1[n][i];
1617 d0factors[2][i] += df[2 * n][i * nquad0] * normal_2[n][i];
1623 for (
int i = 0; i < nquad0; ++i)
1625 d0factors[0][i] = df[0][i] * normal_0[0][i];
1629 for (
int i = 0; i < nquad1; ++i)
1631 d1factors[1][i] = df[1][(i + 1) * nquad0 - 1] * normal_1[0][i];
1632 d1factors[2][i] = df[1][i * nquad0] * normal_2[0][i];
1635 for (
int n = 1; n < ncoords; ++n)
1639 for (
int i = 0; i < nquad0; ++i)
1641 d0factors[0][i] += df[2 * n][i] * normal_0[n][i];
1646 for (
int i = 0; i < nquad1; ++i)
1649 df[2 * n + 1][(i + 1) * nquad0 - 1] * normal_1[n][i];
1650 d1factors[2][i] += df[2 * n + 1][i * nquad0] * normal_2[n][i];
1657 for (
int i = 0; i < nquad0; ++i)
1659 d1factors[0][i] = df[1][0] * normal_0[0][i];
1663 for (
int i = 0; i < nquad1; ++i)
1665 d0factors[1][i] = df[0][0] * normal_1[0][i];
1666 d0factors[2][i] = df[0][0] * normal_2[0][i];
1669 for (
int n = 1; n < ncoords; ++n)
1673 for (
int i = 0; i < nquad0; ++i)
1675 d1factors[0][i] += df[2 * n + 1][0] * normal_0[n][i];
1680 for (
int i = 0; i < nquad1; ++i)
1682 d0factors[1][i] += df[2 * n][0] * normal_1[n][i];
1683 d0factors[2][i] += df[2 * n][0] * normal_2[n][i];
1689 for (
int i = 0; i < nquad0; ++i)
1691 d0factors[0][i] = df[0][0] * normal_0[0][i];
1695 for (
int i = 0; i < nquad1; ++i)
1697 d1factors[1][i] = df[1][0] * normal_1[0][i];
1698 d1factors[2][i] = df[1][0] * normal_2[0][i];
1701 for (
int n = 1; n < ncoords; ++n)
1705 for (
int i = 0; i < nquad0; ++i)
1707 d0factors[0][i] += df[2 * n][0] * normal_0[n][i];
1712 for (
int i = 0; i < nquad1; ++i)
1714 d1factors[1][i] += df[2 * n + 1][0] * normal_1[n][i];
1715 d1factors[2][i] += df[2 * n + 1][0] * normal_2[n][i];
#define ASSERTL0(condition, msg)
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
#define sign(a, b)
return the sign(b)*a
Describes the specification for a Basis.
Defines a specification for a set of points.
unsigned int GetNumPoints() const
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
SpatialDomains::Geometry2DSharedPtr GetGeom2D() const
std::map< int, NormalVector > m_traceNormals
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
SpatialDomains::GeometrySharedPtr GetGeom() const
SpatialDomains::GeometrySharedPtr m_geom
ExpansionSharedPtr GetLeftAdjacentElementExp() const
void ComputeGmatcdotMF(const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble >> &dfdir)
void ComputeLaplacianMetric()
int GetLeftAdjacentElementTrace() const
void ComputeQuadratureMetric()
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
StdRegions::Orientation GetTraceOrient(int trace)
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
const NormalVector & GetTraceNormal(const int id)
virtual void v_ComputeTraceNormal(const int edge)
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
virtual void v_PhysDirectionalDeriv(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &direction, Array< OneD, NekDouble > &out)
Physical derivative along a direction vector.
virtual void v_NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, const Array< OneD, const NekDouble > &Fz, Array< OneD, NekDouble > &outarray)
virtual void v_GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
virtual void v_GetCoord(const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords)
virtual void v_LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
virtual void v_WeakDerivMatrixOp(const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual DNekScalMatSharedPtr v_GetLocMatrix(const MatrixKey &mkey)
virtual void v_FwdTransBndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_ReduceOrderCoeffs(int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual const LibUtilities::BasisSharedPtr & v_GetBasis(int dir) const
virtual void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix(const MatrixKey &mkey)
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Calculate the derivative of the physical points.
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
virtual void v_AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
virtual int v_GetCoordim()
virtual int v_GetNumPoints(const int dir) const
virtual void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
virtual void v_GetEdgeInterpVals(const int edge, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual StdRegions::Orientation v_GetTraceOrient(int edge)
virtual StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const
virtual void v_GetTraceQFactors(const int edge, Array< OneD, NekDouble > &outarray)
void v_DropLocStaticCondMatrix(const MatrixKey &mkey)
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Integrates the specified function over the domain.
virtual void v_MassLevelCurvatureMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_ComputeLaplacianMetric()
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
virtual void v_SVVLaplacianFilter(Array< OneD, NekDouble > &array, const StdRegions::StdMatrixKey &mkey)
virtual void v_LaplacianMatrixOp_MatFree_Kernel(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
virtual DNekMatSharedPtr v_CreateStdMatrix(const StdRegions::StdMatrixKey &mkey)
virtual void v_IProductWRTBase_MatOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Transform a given function from physical quadrature space to coefficient space.
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Calculate the inner product of inarray with respect to the basis B=base0[p]*base1[pq] and put into ou...
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coord, const Array< OneD, const NekDouble > &physvals)
This function evaluates the expansion at a single (arbitrary) point of the domain.
virtual NekDouble v_StdPhysEvaluate(const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
virtual void v_NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble >> &factors, Array< OneD, Array< OneD, NekDouble >> &d0factors, Array< OneD, Array< OneD, NekDouble >> &d1factors)
: This method gets all of the factors which are required as part of the Gradient Jump Penalty stabili...
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
virtual void v_HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_GetTracePhysVals(const int edge, const StdRegions::StdExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
virtual void v_IProductWRTDirectionalDerivBase(const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const
virtual void v_IProductWRTDirectionalDerivBase_SumFac(const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Directinoal Derivative in the modal space in the dir direction of varcoeffs.
virtual void v_GeneralMatrixOp_MatOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_WeakDirectionalDerivMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
virtual void v_IProductWRTDerivBase_MatOp(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0=true, bool doCheckCollDir1=true)
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
void HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
int NumBndryCoeffs(void) const
const LibUtilities::PointsKeyVector GetPointsKeys() const
void MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
void IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
this function calculates the inner product of a given function f with the different modes of the expa...
void LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void GetTraceToElementMap(const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
void GetInteriorMap(Array< OneD, unsigned int > &outarray)
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
void IProductWRTDirectionalDerivBase_SumFac(const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Array< OneD, LibUtilities::BasisSharedPtr > m_base
MatrixType GetMatrixType() const
static void Dgemv(const char &trans, const int &m, const int &n, const double &alpha, const double *a, const int &lda, const double *x, const int &incx, const double &beta, double *y, const int &incy)
BLAS level 2: Matrix vector multiply y = A x where A[m x n].
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
int getNumberOfCoefficients(int Na)
std::shared_ptr< Basis > BasisSharedPtr
void Interp1D(const BasisKey &fbasis0, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, Array< OneD, NekDouble > &to)
this function interpolates a 1D function evaluated at the quadrature points of the basis fbasis0 to ...
std::vector< PointsKey > PointsKeyVector
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
@ eModified_B
Principle Modified Functions .
@ eOrtho_A
Principle Orthogonal Functions .
@ eOrtho_B
Principle Orthogonal Functions .
@ eModified_A
Principle Modified Functions .
std::shared_ptr< SegExp > SegExpSharedPtr
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eMovingRegular
Currently unused.
@ eDeformed
Geometry is curved or has non-constant factors.
std::shared_ptr< Geometry2D > Geometry2DSharedPtr
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
@ eInvLaplacianWithUnityMean
std::shared_ptr< StdTriExp > StdTriExpSharedPtr
The above copyright notice and this permission notice shall be included.
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
static Array< OneD, NekDouble > NullNekDouble1DArray
std::shared_ptr< DNekMat > DNekMatSharedPtr
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
vvtvvtp (scalar times vector plus scalar times vector):
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
svtvp (scalar times vector plus vector): z = alpha*x + y
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Sdiv(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha/y.
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
void Zero(int n, T *x, const int incx)
Zero vector.
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
void Reverse(int n, const T *x, const int incx, T *y, const int incy)
void Vvtvvtp(int n, const T *v, int incv, const T *w, int incw, const T *x, int incx, const T *y, int incy, T *z, int incz)
vvtvvtp (vector times vector plus vector times vector):
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
scalarT< T > sqrt(scalarT< T > in)