Nektar++
Expansion.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: Expansion.cpp
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description: File for Expansion routines
32//
33///////////////////////////////////////////////////////////////////////////////
34
38
39using namespace std;
40
42{
44 : m_indexMapManager(
45 std::bind(&Expansion::CreateIndexMap, this, std::placeholders::_1),
46 std::string("ExpansionIndexMap")),
47 m_geom(pGeom), m_metricinfo(m_geom->GetGeomFactors()),
48 m_elementTraceLeft(-1), m_elementTraceRight(-1)
49{
50 if (!m_metricinfo)
51 {
52 return;
53 }
54
55 if (!m_metricinfo->IsValid())
56 {
57 int nDim = m_base.size();
58 string type = "regular";
59 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
60 {
61 type = "deformed";
62 }
63
64 stringstream err;
65 err << nDim << "D " << type << " Jacobian not positive "
66 << "(element ID = " << m_geom->GetGlobalID() << ") "
67 << "(first vertex ID = " << m_geom->GetVid(0) << ")";
68 NEKERROR(ErrorUtil::ewarning, err.str());
69 }
70
71 m_traceExp.clear();
72}
73
75 : StdExpansion(pSrc), m_indexMapManager(pSrc.m_indexMapManager),
76 m_geom(pSrc.m_geom), m_metricinfo(pSrc.m_metricinfo)
77{
78}
79
81{
82}
83
85 const LocalRegions::MatrixKey &mkey)
86{
87 return v_GetLocMatrix(mkey);
88}
89
91{
92 return v_DropLocMatrix(mkey);
93}
94
96 const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
97{
98 return v_BuildTransformationMatrix(r_bnd, matrixType);
99}
100
102{
103 return v_BuildVertexMatrix(r_bnd);
104}
105
107 const NekDouble *data, const std::vector<unsigned int> &nummodes,
108 const int nmodes_offset, NekDouble *coeffs,
109 std::vector<LibUtilities::BasisType> &fromType)
110{
111 v_ExtractDataToCoeffs(data, nummodes, nmodes_offset, coeffs, fromType);
112}
113
115 const int edge, const std::shared_ptr<Expansion> &EdgeExp,
118{
119 v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fx, Fy, outarray);
120}
121
123 const int edge, const std::shared_ptr<Expansion> &EdgeExp,
125{
126 v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fn, outarray);
127}
128
130 const int face, const std::shared_ptr<Expansion> &FaceExp,
132{
133 v_AddFaceNormBoundaryInt(face, FaceExp, Fn, outarray);
134}
135
136void Expansion::DGDeriv(const int dir,
137 const Array<OneD, const NekDouble> &inarray,
140 Array<OneD, NekDouble> &outarray)
141{
142 v_DGDeriv(dir, inarray, EdgeExp, coeffs, outarray);
143}
144
146{
147 return v_VectorFlux(vec);
148}
149
152 Array<OneD, Array<OneD, NekDouble>> &d0factors,
153 Array<OneD, Array<OneD, NekDouble>> &d1factors)
154{
155 return v_NormalTraceDerivFactors(factors, d0factors, d1factors);
156}
157
159 const StdRegions::MatrixType mtype,
161 const StdRegions::VarCoeffMap &varcoeffs)
162{
163 MatrixKey mkey(mtype, DetShapeType(), *this, factors, varcoeffs);
164 return GetLocMatrix(mkey);
165}
166
168{
169 return m_geom;
170}
171
173{
174 // Clear metrics
175 m_metrics.clear();
176
177 // Regenerate geometry factors
178 m_metricinfo = m_geom->GetGeomFactors();
179}
180
182{
183 IndexMapValuesSharedPtr returnval;
184
185 IndexMapType itype = ikey.GetIndexMapType();
186
187 int entity = ikey.GetIndexEntity();
188
190
193
194 switch (itype)
195 {
196 case eEdgeToElement:
197 {
198 GetTraceToElementMap(entity, map, sign, orient);
199 }
200 break;
201 case eFaceToElement:
202 {
203 GetTraceToElementMap(entity, map, sign, orient);
204 }
205 break;
206 case eEdgeInterior:
207 {
208 ASSERTL0(false, "Boundary Index Map not implemented yet.");
209 // v_GetEdgeInteriorMap(entity,orient,map,sign);
210 }
211 break;
212 case eFaceInterior:
213 {
214 ASSERTL0(false, "Boundary Index Map not implemented yet.");
215 // v_GetFaceInteriorMap(entity,orient,map,sign);
216 }
217 break;
218 case eBoundary:
219 {
220 ASSERTL0(false, "Boundary Index Map not implemented yet.");
221 }
222 break;
223 case eVertex:
224 {
225 ASSERTL0(false, "Vertex Index Map not implemented yet.");
226 }
227 break;
228 default:
229 {
230 ASSERTL0(false, "The Index Map you are requiring "
231 "is not between the possible options.");
232 }
233 }
234
236
237 for (int i = 0; i < map.size(); i++)
238 {
239 (*returnval)[i].index = map[i];
240 (*returnval)[i].sign = sign[i];
241 }
242
243 return returnval;
244}
245
247{
248 return m_metricinfo;
249}
250
252{
253 std::map<int, NormalVector>::const_iterator x;
254 x = m_traceNormals.find(id);
255
256 // if edge normal not defined compute it
257 if (x == m_traceNormals.end())
258 {
260 x = m_traceNormals.find(id);
261 }
262 return x->second;
263}
264
266 [[maybe_unused]] const LocalRegions::MatrixKey &mkey)
267{
268 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
270}
271
273{
274 DNekScalBlkMatSharedPtr returnval;
275
277 "Geometric information is not set up");
278
279 // set up block matrix system
280 unsigned int nbdry = NumBndryCoeffs();
281 unsigned int nint = (unsigned int)(m_ncoeffs - nbdry);
282 unsigned int exp_size[] = {nbdry, nint};
283 unsigned int nblks = 2;
285 nblks, nblks, exp_size, exp_size);
286 // Really need a constructor which takes Arrays
287 NekDouble factor = 1.0;
288
289 switch (mkey.GetMatrixType())
290 {
291 // this can only use stdregions statically condensed system
292 // for mass matrix
294 if ((m_metricinfo->GetGtype() == SpatialDomains::eDeformed) ||
295 (mkey.GetNVarCoeff()))
296 {
297 factor = 1.0;
298 goto UseLocRegionsMatrix;
299 }
300 else
301 {
302 factor = (m_metricinfo->GetJac(GetPointsKeys()))[0];
303 goto UseStdRegionsMatrix;
304 }
305 break;
306 default: // use Deformed case for both
307 // regular and deformed geometries
308 factor = 1.0;
309 goto UseLocRegionsMatrix;
310 break;
311 UseStdRegionsMatrix:
312 {
313 NekDouble invfactor = 1.0 / factor;
314 NekDouble one = 1.0;
317 DNekMatSharedPtr Asubmat;
318
319 returnval->SetBlock(
320 0, 0,
322 factor, Asubmat = mat->GetBlock(0, 0)));
323 returnval->SetBlock(
324 0, 1,
326 one, Asubmat = mat->GetBlock(0, 1)));
327 returnval->SetBlock(
328 1, 0,
330 factor, Asubmat = mat->GetBlock(1, 0)));
331 returnval->SetBlock(
332 1, 1,
334 invfactor, Asubmat = mat->GetBlock(1, 1)));
335 }
336 break;
337 UseLocRegionsMatrix:
338 {
339 int i, j;
340 NekDouble invfactor = 1.0 / factor;
341 NekDouble one = 1.0;
342 DNekScalMat &mat = *GetLocMatrix(mkey);
351
352 Array<OneD, unsigned int> bmap(nbdry);
353 Array<OneD, unsigned int> imap(nint);
354 GetBoundaryMap(bmap);
355 GetInteriorMap(imap);
356
357 for (i = 0; i < nbdry; ++i)
358 {
359 for (j = 0; j < nbdry; ++j)
360 {
361 (*A)(i, j) = mat(bmap[i], bmap[j]);
362 }
363
364 for (j = 0; j < nint; ++j)
365 {
366 (*B)(i, j) = mat(bmap[i], imap[j]);
367 }
368 }
369
370 for (i = 0; i < nint; ++i)
371 {
372 for (j = 0; j < nbdry; ++j)
373 {
374 (*C)(i, j) = mat(imap[i], bmap[j]);
375 }
376
377 for (j = 0; j < nint; ++j)
378 {
379 (*D)(i, j) = mat(imap[i], imap[j]);
380 }
381 }
382
383 // Calculate static condensed system
384 if (nint)
385 {
386 D->Invert();
387 (*B) = (*B) * (*D);
388 (*A) = (*A) - (*B) * (*C);
389 }
390
392
393 returnval->SetBlock(
394 0, 0,
395 Atmp =
397 returnval->SetBlock(
398 0, 1,
400 returnval->SetBlock(
401 1, 0,
402 Atmp =
404 returnval->SetBlock(
405 1, 1,
407 D));
408
409 // Remove the local matrix from manager if using this option since
410 // we assume it is only created to generate static condensed system
411 v_DropLocMatrix(mkey);
412 }
413 }
414 return returnval;
415}
416
418 [[maybe_unused]] const LocalRegions::MatrixKey &mkey)
419{
420 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
421}
422
424 const Array<OneD, const NekDouble> &inarray,
425 Array<OneD, NekDouble> &outarray)
426{
427 const int nqtot = GetTotPoints();
428
429 if (m_metrics.count(eMetricQuadrature) == 0)
430 {
432 }
433
434 Vmath::Vmul(nqtot, m_metrics[eMetricQuadrature], 1, inarray, 1, outarray,
435 1);
436}
437
439 const Array<OneD, const NekDouble> &inarray,
440 Array<OneD, NekDouble> &outarray)
441{
442 const int nqtot = GetTotPoints();
443
444 if (m_metrics.count(eMetricQuadrature) == 0)
445 {
447 }
448
449 Vmath::Vdiv(nqtot, inarray, 1, m_metrics[eMetricQuadrature], 1, outarray,
450 1);
451}
452
454{
456}
457
459{
460 unsigned int nqtot = GetTotPoints();
461 SpatialDomains::GeomType type = m_metricinfo->GetGtype();
463 if (type == SpatialDomains::eRegular ||
465 {
467 Array<OneD, NekDouble>(nqtot, m_metricinfo->GetJac(p)[0]);
468 }
469 else
470 {
472 }
473
476}
477
479{
480 int nquad = GetTotPoints();
481 int ntraces = GetNtraces();
482 int ndir = m_base.size();
483
485 Array<OneD, NekDouble> phys(nquad);
486
487 Array<OneD, Array<OneD, int>> traceids(ntraces);
488
489 int tottracepts = 0;
490 for (int i = 0; i < ntraces; ++i)
491 {
492 GetTracePhysMap(i, traceids[i]);
493 tottracepts += GetTraceNumPoints(i);
494 }
495
496 // initialise array to null so can call for
497 // differnt dimensions
499
500 DerivMat = Array<OneD, DNekMatSharedPtr>(ndir);
501
502 for (int i = 0; i < ndir; ++i)
503 {
504 Deriv[i] = Array<OneD, NekDouble>(nquad);
505 DerivMat[i] =
507 }
508
509 for (int i = 0; i < m_ncoeffs; ++i)
510 {
511 Vmath::Zero(m_ncoeffs, coeffs, 1);
512 coeffs[i] = 1.0;
513 BwdTrans(coeffs, phys);
514
515 // dphi_i/d\xi_1, dphi_i/d\xi_2 dphi_i/d\xi_3
516 StdPhysDeriv(phys, Deriv[0], Deriv[1], Deriv[2]);
517
518 int cnt = 0;
519 for (int j = 0; j < ntraces; ++j)
520 {
521 int nTracePts = GetTraceNumPoints(j);
522 for (int k = 0; k < nTracePts; ++k)
523 {
524 for (int d = 0; d < ndir; ++d)
525 {
526 (*DerivMat[d])(i, cnt + k) = Deriv[d][traceids[j][k]];
527 }
528 }
529 cnt += nTracePts;
530 }
531 }
532}
533
535 Array<OneD, NekDouble> &coords_1,
536 Array<OneD, NekDouble> &coords_2)
537{
538 ASSERTL1(m_geom, "m_geom not defined");
539
540 // get physical points defined in Geom
541 m_geom->FillGeom();
542
543 const int expDim = m_base.size();
544 int nqGeom = 1;
545 bool doCopy = true;
546
549
550 for (int i = 0; i < expDim; ++i)
551 {
552 CBasis[i] = m_geom->GetXmap()->GetBasis(i);
553 nqGeom *= CBasis[i]->GetNumPoints();
554 doCopy = doCopy &&
555 m_base[i]->GetBasisKey().SamePoints(CBasis[i]->GetBasisKey());
556 }
557
558 tmp[0] = coords_0;
559 tmp[1] = coords_1;
560 tmp[2] = coords_2;
561
562 if (doCopy)
563 {
564 for (int i = 0; i < m_geom->GetCoordim(); ++i)
565 {
566 m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmp[i]);
567 }
568 }
569 else
570 {
571 for (int i = 0; i < m_geom->GetCoordim(); ++i)
572 {
573 Array<OneD, NekDouble> tmpGeom(nqGeom);
574 m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmpGeom);
575
576 switch (expDim)
577 {
578 case 1:
579 {
581 CBasis[0]->GetPointsKey(), &tmpGeom[0],
582 m_base[0]->GetPointsKey(), &tmp[i][0]);
583 break;
584 }
585 case 2:
586 {
588 CBasis[0]->GetPointsKey(), CBasis[1]->GetPointsKey(),
589 &tmpGeom[0], m_base[0]->GetPointsKey(),
590 m_base[1]->GetPointsKey(), &tmp[i][0]);
591 break;
592 }
593 case 3:
594 {
596 CBasis[0]->GetPointsKey(), CBasis[1]->GetPointsKey(),
597 CBasis[2]->GetPointsKey(), &tmpGeom[0],
598 m_base[0]->GetPointsKey(), m_base[1]->GetPointsKey(),
599 m_base[2]->GetPointsKey(), &tmp[i][0]);
600 break;
601 }
602 }
603 }
604 }
605}
606
608 const Array<OneD, const NekDouble> &direction,
610{
611 int shapedim = dfdir.size();
612 int coordim = GetCoordim();
613 int nqtot = direction.size() / coordim;
614
615 for (int j = 0; j < shapedim; j++)
616 {
617 dfdir[j] = Array<OneD, NekDouble>(nqtot, 0.0);
618 for (int k = 0; k < coordim; k++)
619 {
620 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
621 {
622 Vmath::Vvtvp(nqtot, &df[shapedim * k + j][0], 1,
623 &direction[k * nqtot], 1, &dfdir[j][0], 1,
624 &dfdir[j][0], 1);
625 }
626 else
627 {
628 Vmath::Svtvp(nqtot, df[shapedim * k + j][0],
629 &direction[k * nqtot], 1, &dfdir[j][0], 1,
630 &dfdir[j][0], 1);
631 }
632 }
633 }
634}
635
636// Get Moving frames
638 const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
639{
640 int coordim = GetCoordim();
641
642 int nquad0, nquad1, nquad2;
643 int nqtot = 1;
644 switch (shapedim)
645 {
646 case 2:
647 {
648 nquad0 = m_base[0]->GetNumPoints();
649 nquad1 = m_base[1]->GetNumPoints();
650 nqtot = nquad0 * nquad1;
651 break;
652 }
653 case 3:
654 {
655 nquad0 = m_base[0]->GetNumPoints();
656 nquad1 = m_base[1]->GetNumPoints();
657 nquad2 = m_base[2]->GetNumPoints();
658 nqtot = nquad0 * nquad1 * nquad2;
659 break;
660 }
661 default:
662 break;
663 }
664
665 StdRegions::VarCoeffType MMFCoeffs[15] = {
674
675 Array<OneD, NekDouble> MF(coordim * nqtot);
676 Array<OneD, NekDouble> tmp(nqtot);
677 for (int k = 0; k < coordim; k++)
678 {
679 StdRegions::VarCoeffMap::const_iterator MFdir =
680 varcoeffs.find(MMFCoeffs[5 * dir + k]);
681 tmp = MFdir->second.GetValue();
682
683 Vmath::Vcopy(nqtot, &tmp[0], 1, &MF[k * nqtot], 1);
684 }
685
686 return MF;
687}
688
689// Get magnitude of MF
691 const int dir, const StdRegions::VarCoeffMap &varcoeffs)
692{
693 int indxDiv = 3;
694
695 StdRegions::VarCoeffType MMFCoeffs[15] = {
704
705 StdRegions::VarCoeffMap::const_iterator MFdir =
706 varcoeffs.find(MMFCoeffs[5 * dir + indxDiv]);
707 Array<OneD, NekDouble> MFDiv = MFdir->second.GetValue();
708
709 return MFDiv;
710}
711
712// Get magnitude of MF
714 const int dir, const StdRegions::VarCoeffMap &varcoeffs)
715{
716 int indxMag = 4;
717
718 StdRegions::VarCoeffType MMFCoeffs[15] = {
727
728 StdRegions::VarCoeffMap::const_iterator MFdir =
729 varcoeffs.find(MMFCoeffs[5 * dir + indxMag]);
730 Array<OneD, NekDouble> MFmag = MFdir->second.GetValue();
731
732 return MFmag;
733}
734
736 [[maybe_unused]] const DNekScalMatSharedPtr &r_bnd,
737 [[maybe_unused]] const StdRegions::MatrixType matrixType)
738{
739 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
741}
742
744 [[maybe_unused]] const DNekScalMatSharedPtr &r_bnd)
745{
746 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
748}
749
751 [[maybe_unused]] const NekDouble *data,
752 [[maybe_unused]] const std::vector<unsigned int> &nummodes,
753 [[maybe_unused]] const int nmodes_offset,
754 [[maybe_unused]] NekDouble *coeffs,
755 [[maybe_unused]] std::vector<LibUtilities::BasisType> &fromType)
756{
757 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
758}
759
761 [[maybe_unused]] const int edge,
762 [[maybe_unused]] const std::shared_ptr<Expansion> &EdgeExp,
763 [[maybe_unused]] const Array<OneD, const NekDouble> &Fx,
764 [[maybe_unused]] const Array<OneD, const NekDouble> &Fy,
765 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
766{
767 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
768}
769
771 [[maybe_unused]] const int edge,
772 [[maybe_unused]] const std::shared_ptr<Expansion> &EdgeExp,
773 [[maybe_unused]] const Array<OneD, const NekDouble> &Fn,
774 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
775{
776 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
777}
778
780 [[maybe_unused]] const int face,
781 [[maybe_unused]] const std::shared_ptr<Expansion> &FaceExp,
782 [[maybe_unused]] const Array<OneD, const NekDouble> &Fn,
783 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
784{
785 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
786}
787
789 [[maybe_unused]] const int dir,
790 [[maybe_unused]] const Array<OneD, const NekDouble> &inarray,
791 [[maybe_unused]] Array<OneD, ExpansionSharedPtr> &EdgeExp,
792 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &coeffs,
793 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
794{
795 NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
796}
797
799 [[maybe_unused]] const Array<OneD, Array<OneD, NekDouble>> &vec)
800{
802 "This function is only valid for "
803 "shape expansion in LocalRegions, not parant class");
804 return 0.0;
805}
806
808 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &factors,
809 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &d0factors,
810 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &d1factors)
811{
813 "This function is only valid for "
814 "shape expansion in LocalRegions, not parant class");
815}
816
818{
820}
821
823 [[maybe_unused]] StdRegions::Orientation dir,
824 [[maybe_unused]] Array<OneD, const NekDouble> &inarray,
825 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
826{
827 NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
828}
829
831 [[maybe_unused]] const int trace,
832 [[maybe_unused]] Array<OneD, NekDouble> &outarray)
833{
835 "Method does not exist for this shape or library");
836}
837
839 [[maybe_unused]] const int trace,
840 [[maybe_unused]] const StdRegions::StdExpansionSharedPtr &TraceExp,
841 [[maybe_unused]] const Array<OneD, const NekDouble> &inarray,
842 [[maybe_unused]] Array<OneD, NekDouble> &outarray,
843 [[maybe_unused]] StdRegions::Orientation orient)
844{
846 "Method does not exist for this shape or library");
847}
848
849void Expansion::v_GetTracePhysMap([[maybe_unused]] const int edge,
850 [[maybe_unused]] Array<OneD, int> &outarray)
851{
853 "Method does not exist for this shape or library");
854}
855
857 [[maybe_unused]] const StdRegions::Orientation orient,
858 [[maybe_unused]] Array<OneD, int> &idmap, [[maybe_unused]] const int nq0,
859 [[maybe_unused]] const int nq1)
860{
862 "Method does not exist for this shape or library");
863}
864
865void Expansion::v_GenTraceExp([[maybe_unused]] const int traceid,
866 [[maybe_unused]] ExpansionSharedPtr &exp)
867{
869 "Method does not exist for this shape or library");
870}
871
872void Expansion::v_ComputeTraceNormal([[maybe_unused]] const int id)
873{
874 ASSERTL0(false, "Cannot compute trace normal for this expansion.");
875}
876
878{
879 NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
881}
882
884 [[maybe_unused]] Array<OneD, const NekDouble> &normal)
885{
886 NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
887}
888
889void Expansion::v_SetUpPhysNormals([[maybe_unused]] const int edge)
890{
891 NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
892}
893
895 [[maybe_unused]] const int trace,
896 [[maybe_unused]] const Array<OneD, const NekDouble> &primCoeffs,
897 [[maybe_unused]] DNekMatSharedPtr &inoutmat)
898{
899 NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
900}
901
903 [[maybe_unused]] const int traceid,
904 [[maybe_unused]] const Array<OneD, const NekDouble> &primCoeffs,
905 [[maybe_unused]] const Array<OneD, NekDouble> &incoeffs,
906 [[maybe_unused]] Array<OneD, NekDouble> &coeffs)
907{
908 NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
909}
910
911void Expansion::v_TraceNormLen([[maybe_unused]] const int traceid,
912 [[maybe_unused]] NekDouble &h,
913 [[maybe_unused]] NekDouble &p)
914{
915 NEKERROR(ErrorUtil::efatal, "This method has not been defined");
916}
917
919 const int nbnd) const
920{
921 auto x = m_elmtBndNormDirElmtLen.find(nbnd);
923 "m_elmtBndNormDirElmtLen normal not computed.");
924 return x->second;
925}
926
928 [[maybe_unused]] const int dir,
929 [[maybe_unused]] const Array<OneD, const NekDouble> &inarray,
930 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &outarray)
931{
932 NEKERROR(ErrorUtil::efatal, "v_AlignVectorToCollapsedDir is not defined");
933}
934} // namespace Nektar::LocalRegions
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
Definition: ErrorUtil.hpp:265
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:47
std::map< int, NormalVector > m_traceNormals
Definition: Expansion.h:276
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
Definition: Expansion.h:286
virtual void v_SetPhysNormals(Array< OneD, const NekDouble > &normal)
Definition: Expansion.cpp:883
SpatialDomains::GeometrySharedPtr GetGeom() const
Definition: Expansion.cpp:167
virtual void v_TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
Definition: Expansion.cpp:911
IndexMapValuesSharedPtr CreateIndexMap(const IndexMapKey &ikey)
Definition: Expansion.cpp:181
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:273
void GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.h:209
void DropLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:90
const SpatialDomains::GeomFactorsSharedPtr & GetMetricInfo() const
Definition: Expansion.cpp:246
DNekMatSharedPtr BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:101
virtual void v_DivideByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:438
virtual void v_GenTraceExp(const int traceid, ExpansionSharedPtr &exp)
Definition: Expansion.cpp:865
void v_SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Definition: Expansion.cpp:822
virtual void v_AddRobinTraceContribution(const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
Definition: Expansion.cpp:902
void ComputeGmatcdotMF(const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble > > &dfdir)
Definition: Expansion.cpp:607
virtual void v_GetTraceQFactors(const int trace, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:830
virtual void v_ComputeTraceNormal(const int id)
Definition: Expansion.cpp:872
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
Definition: Expansion.cpp:534
Array< OneD, NekDouble > GetMFMag(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:713
virtual DNekMatSharedPtr v_BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:743
virtual DNekScalMatSharedPtr v_GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:265
virtual void v_AddFaceNormBoundaryInt(const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:779
virtual void v_DGDeriv(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:788
virtual DNekMatSharedPtr v_BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:735
const Array< OneD, const NekDouble > & GetElmtBndNormDirElmtLen(const int nbnd) const
Definition: Expansion.cpp:918
virtual void v_AddRobinMassMatrix(const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
Definition: Expansion.cpp:894
std::map< int, ExpansionWeakPtr > m_traceExp
Definition: Expansion.h:272
virtual void v_NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
Definition: Expansion.cpp:807
void AddEdgeNormBoundaryInt(const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:114
void AddFaceNormBoundaryInt(const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:129
virtual void v_DropLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:417
virtual NekDouble v_VectorFlux(const Array< OneD, Array< OneD, NekDouble > > &vec)
Definition: Expansion.cpp:798
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:274
virtual StdRegions::Orientation v_GetTraceOrient(int trace)
Definition: Expansion.cpp:817
void DGDeriv(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:136
void StdDerivBaseOnTraceMat(Array< OneD, DNekMatSharedPtr > &DerivMat)
Definition: Expansion.cpp:478
void ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
Definition: Expansion.cpp:106
void v_MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Definition: Expansion.cpp:423
virtual void v_ReOrientTracePhysMap(const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
Definition: Expansion.cpp:856
void NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
Definition: Expansion.cpp:150
virtual void v_ComputeLaplacianMetric()
Definition: Expansion.h:311
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals()
Definition: Expansion.cpp:877
virtual void v_ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
Definition: Expansion.cpp:750
virtual void v_AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Expansion.cpp:927
virtual void v_AddEdgeNormBoundaryInt(const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:760
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:84
virtual void v_SetUpPhysNormals(const int id)
Definition: Expansion.cpp:889
virtual void v_GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.cpp:849
Array< OneD, NekDouble > GetMFDiv(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:690
const NormalVector & GetTraceNormal(const int id)
Definition: Expansion.cpp:251
Array< OneD, NekDouble > GetMF(const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:637
virtual void v_GetTracePhysVals(const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
Definition: Expansion.cpp:838
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:43
DNekMatSharedPtr BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:95
NekDouble VectorFlux(const Array< OneD, Array< OneD, NekDouble > > &vec)
Definition: Expansion.cpp:145
StdRegions::Orientation GetIndexOrientation() const
Definition: IndexMapKey.h:94
IndexMapType GetIndexMapType() const
Definition: IndexMapKey.h:89
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
void GetBoundaryMap(Array< OneD, unsigned int > &outarray)
Definition: StdExpansion.h:678
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:134
int GetTraceNumPoints(const int i) const
This function returns the number of quadrature points belonging to the i-th trace.
Definition: StdExpansion.h:281
const LibUtilities::PointsKeyVector GetPointsKeys() const
DNekBlkMatSharedPtr GetStdStaticCondMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:617
int GetNtraces() const
Returns the number of trace elements connected to this element.
Definition: StdExpansion.h:354
virtual void v_MultiplyByStdQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void GetTraceToElementMap(const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
Definition: StdExpansion.h:693
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:370
void GetInteriorMap(Array< OneD, unsigned int > &outarray)
Definition: StdExpansion.h:683
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Backward transformation from coefficient space to physical space.
Definition: StdExpansion.h:433
void StdPhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Definition: StdExpansion.h:891
Array< OneD, LibUtilities::BasisSharedPtr > m_base
MatrixType GetMatrixType() const
Definition: StdMatrixKey.h:83
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 ...
Definition: Interp.cpp:47
void Interp3D(const BasisKey &fbasis0, const BasisKey &fbasis1, const BasisKey &fbasis2, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, const BasisKey &tbasis2, Array< OneD, NekDouble > &to)
this function interpolates a 3D function evaluated at the quadrature points of the 3D basis,...
Definition: Interp.cpp:162
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis,...
Definition: Interp.cpp:101
std::vector< PointsKey > PointsKeyVector
Definition: Points.h:231
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:66
std::shared_ptr< IndexMapValues > IndexMapValuesSharedPtr
Definition: IndexMapKey.h:126
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:60
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eNoGeomType
No type defined.
@ eMovingRegular
Currently unused.
@ eDeformed
Geometry is curved or has non-constant factors.
std::shared_ptr< Geometry > GeometrySharedPtr
Definition: Geometry.h:51
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:430
std::map< StdRegions::VarCoeffType, VarCoeffEntry > VarCoeffMap
Definition: StdRegions.hpp:375
std::vector< double > d(NPUPPER *NPUPPER)
StdRegions::ConstFactorMap factors
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekBlkMat > DNekBlkMatSharedPtr
Definition: NekTypeDefs.hpp:77
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:79
static DNekMatSharedPtr NullDNekMatSharedPtr
Definition: NekTypeDefs.hpp:83
static Array< OneD, NekDouble > NullNekDouble1DArray
static DNekScalMatSharedPtr NullDNekScalMatSharedPtr
Definition: NekTypeDefs.hpp:84
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble
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.
Definition: Vmath.hpp:72
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.
Definition: Vmath.hpp:396
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
Definition: Vmath.hpp:366
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.
Definition: Vmath.hpp:126
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
STL namespace.