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 
35 #include <boost/core/ignore_unused.hpp>
36 
38 #include <LocalRegions/Expansion.h>
39 #include <LocalRegions/MatrixKey.h>
40 
41 using namespace std;
42 
43 namespace Nektar
44 {
45 namespace LocalRegions
46 {
47 Expansion::Expansion(SpatialDomains::GeometrySharedPtr pGeom)
48  : m_indexMapManager(
49  std::bind(&Expansion::CreateIndexMap, this, std::placeholders::_1),
50  std::string("ExpansionIndexMap")),
51  m_geom(pGeom), m_metricinfo(m_geom->GetGeomFactors()),
52  m_elementTraceLeft(-1), m_elementTraceRight(-1)
53 {
54  if (!m_metricinfo)
55  {
56  return;
57  }
58 
59  if (!m_metricinfo->IsValid())
60  {
61  int nDim = m_base.size();
62  string type = "regular";
63  if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
64  {
65  type = "deformed";
66  }
67 
68  stringstream err;
69  err << nDim << "D " << type << " Jacobian not positive "
70  << "(element ID = " << m_geom->GetGlobalID() << ") "
71  << "(first vertex ID = " << m_geom->GetVid(0) << ")";
72  NEKERROR(ErrorUtil::ewarning, err.str());
73  }
74 
75  m_traceExp.empty();
76 }
77 
79  : StdExpansion(pSrc), m_indexMapManager(pSrc.m_indexMapManager),
80  m_geom(pSrc.m_geom), m_metricinfo(pSrc.m_metricinfo)
81 {
82 }
83 
85 {
86 }
87 
89  const LocalRegions::MatrixKey &mkey)
90 {
91  return v_GetLocMatrix(mkey);
92 }
93 
95  const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
96 {
97  return v_BuildTransformationMatrix(r_bnd, matrixType);
98 }
99 
101 {
102  return v_BuildVertexMatrix(r_bnd);
103 }
104 
106  const NekDouble *data, const std::vector<unsigned int> &nummodes,
107  const int nmodes_offset, NekDouble *coeffs,
108  std::vector<LibUtilities::BasisType> &fromType)
109 {
110  v_ExtractDataToCoeffs(data, nummodes, nmodes_offset, coeffs, fromType);
111 }
112 
114  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
117 {
118  v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fx, Fy, outarray);
119 }
120 
122  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
124 {
125  v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fn, outarray);
126 }
127 
129  const int face, const std::shared_ptr<Expansion> &FaceExp,
131 {
132  v_AddFaceNormBoundaryInt(face, FaceExp, Fn, outarray);
133 }
134 
135 void Expansion::DGDeriv(const int dir,
136  const Array<OneD, const NekDouble> &inarray,
139  Array<OneD, NekDouble> &outarray)
140 {
141  v_DGDeriv(dir, inarray, EdgeExp, coeffs, outarray);
142 }
143 
145 {
146  return v_VectorFlux(vec);
147 }
148 
150  Array<OneD, Array<OneD, NekDouble>> &factors,
151  Array<OneD, Array<OneD, NekDouble>> &d0factors,
152  Array<OneD, Array<OneD, NekDouble>> &d1factors)
153 {
154  return v_NormalTraceDerivFactors(factors, d0factors, d1factors);
155 }
156 
158  const StdRegions::MatrixType mtype,
159  const StdRegions::ConstFactorMap &factors,
160  const StdRegions::VarCoeffMap &varcoeffs)
161 {
162  MatrixKey mkey(mtype, DetShapeType(), *this, factors, varcoeffs);
163  return GetLocMatrix(mkey);
164 }
165 
167 {
168  return m_geom;
169 }
170 
172 {
173  // Clear metrics
174  m_metrics.clear();
175 
176  // Regenerate geometry factors
177  m_metricinfo = m_geom->GetGeomFactors();
178 }
179 
181 {
182  IndexMapValuesSharedPtr returnval;
183 
184  IndexMapType itype = ikey.GetIndexMapType();
185 
186  int entity = ikey.GetIndexEntity();
187 
189 
192 
193  switch (itype)
194  {
195  case eEdgeToElement:
196  {
197  GetTraceToElementMap(entity, map, sign, orient);
198  }
199  break;
200  case eFaceToElement:
201  {
202  GetTraceToElementMap(entity, map, sign, orient);
203  }
204  break;
205  case eEdgeInterior:
206  {
207  ASSERTL0(false, "Boundary Index Map not implemented yet.");
208  // v_GetEdgeInteriorMap(entity,orient,map,sign);
209  }
210  break;
211  case eFaceInterior:
212  {
213  ASSERTL0(false, "Boundary Index Map not implemented yet.");
214  // v_GetFaceInteriorMap(entity,orient,map,sign);
215  }
216  break;
217  case eBoundary:
218  {
219  ASSERTL0(false, "Boundary Index Map not implemented yet.");
220  }
221  break;
222  case eVertex:
223  {
224  ASSERTL0(false, "Vertex Index Map not implemented yet.");
225  }
226  break;
227  default:
228  {
229  ASSERTL0(false, "The Index Map you are requiring "
230  "is not between the possible options.");
231  }
232  }
233 
234  returnval = MemoryManager<IndexMapValues>::AllocateSharedPtr(map.size());
235 
236  for (int i = 0; i < map.size(); i++)
237  {
238  (*returnval)[i].index = map[i];
239  (*returnval)[i].sign = sign[i];
240  }
241 
242  return returnval;
243 }
244 
246 {
247  return m_metricinfo;
248 }
249 
251 {
252  std::map<int, NormalVector>::const_iterator x;
253  x = m_traceNormals.find(id);
254 
255  // if edge normal not defined compute it
256  if (x == m_traceNormals.end())
257  {
259  x = m_traceNormals.find(id);
260  }
261  return x->second;
262 }
263 
265  const LocalRegions::MatrixKey &mkey)
266 {
267  boost::ignore_unused(mkey);
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
293  case StdRegions::eMass:
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);
345  DNekMatSharedPtr B =
347  DNekMatSharedPtr C =
349  DNekMatSharedPtr D =
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  }
410  return returnval;
411 }
412 
414  const Array<OneD, const NekDouble> &inarray,
415  Array<OneD, NekDouble> &outarray)
416 {
417  const int nqtot = GetTotPoints();
418 
419  if (m_metrics.count(eMetricQuadrature) == 0)
420  {
422  }
423 
424  Vmath::Vmul(nqtot, m_metrics[eMetricQuadrature], 1, inarray, 1, outarray,
425  1);
426 }
427 
429  const Array<OneD, const NekDouble> &inarray,
430  Array<OneD, NekDouble> &outarray)
431 {
432  const int nqtot = GetTotPoints();
433 
434  if (m_metrics.count(eMetricQuadrature) == 0)
435  {
437  }
438 
439  Vmath::Vdiv(nqtot, inarray, 1, m_metrics[eMetricQuadrature], 1, outarray,
440  1);
441 }
442 
444 {
446 }
447 
449 {
450  unsigned int nqtot = GetTotPoints();
451  SpatialDomains::GeomType type = m_metricinfo->GetGtype();
453  if (type == SpatialDomains::eRegular ||
455  {
457  Array<OneD, NekDouble>(nqtot, m_metricinfo->GetJac(p)[0]);
458  }
459  else
460  {
462  }
463 
466 }
467 
469 {
470  int nquad = GetTotPoints();
471  int ntraces = GetNtraces();
472  int ndir = m_base.size();
473 
475  Array<OneD, NekDouble> phys(nquad);
476 
477  Array<OneD, Array<OneD, int>> traceids(ntraces);
478 
479  int tottracepts = 0;
480  for (int i = 0; i < ntraces; ++i)
481  {
482  GetTracePhysMap(i, traceids[i]);
483  tottracepts += GetTraceNumPoints(i);
484  }
485 
486  // initialise array to null so can call for
487  // differnt dimensions
489 
490  DerivMat = Array<OneD, DNekMatSharedPtr>(ndir);
491 
492  for (int i = 0; i < ndir; ++i)
493  {
494  Deriv[i] = Array<OneD, NekDouble>(nquad);
495  DerivMat[i] =
497  }
498 
499  for (int i = 0; i < m_ncoeffs; ++i)
500  {
501  Vmath::Zero(m_ncoeffs, coeffs, 1);
502  coeffs[i] = 1.0;
503  BwdTrans(coeffs, phys);
504 
505  // dphi_i/d\xi_1, dphi_i/d\xi_2 dphi_i/d\xi_3
506  StdPhysDeriv(phys, Deriv[0], Deriv[1], Deriv[2]);
507 
508  int cnt = 0;
509  for (int j = 0; j < ntraces; ++j)
510  {
511  int nTracePts = GetTraceNumPoints(j);
512  for (int k = 0; k < nTracePts; ++k)
513  {
514  for (int d = 0; d < ndir; ++d)
515  {
516  (*DerivMat[d])(i, cnt + k) = Deriv[d][traceids[j][k]];
517  }
518  }
519  cnt += nTracePts;
520  }
521  }
522 }
523 
525  Array<OneD, NekDouble> &coords_1,
526  Array<OneD, NekDouble> &coords_2)
527 {
528  ASSERTL1(m_geom, "m_geom not defined");
529 
530  // get physical points defined in Geom
531  m_geom->FillGeom();
532 
533  const int expDim = m_base.size();
534  int nqGeom = 1;
535  bool doCopy = true;
536 
539 
540  for (int i = 0; i < expDim; ++i)
541  {
542  CBasis[i] = m_geom->GetXmap()->GetBasis(i);
543  nqGeom *= CBasis[i]->GetNumPoints();
544  doCopy = doCopy &&
545  m_base[i]->GetBasisKey().SamePoints(CBasis[i]->GetBasisKey());
546  }
547 
548  tmp[0] = coords_0;
549  tmp[1] = coords_1;
550  tmp[2] = coords_2;
551 
552  if (doCopy)
553  {
554  for (int i = 0; i < m_geom->GetCoordim(); ++i)
555  {
556  m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmp[i]);
557  }
558  }
559  else
560  {
561  for (int i = 0; i < m_geom->GetCoordim(); ++i)
562  {
563  Array<OneD, NekDouble> tmpGeom(nqGeom);
564  m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmpGeom);
565 
566  switch (expDim)
567  {
568  case 1:
569  {
571  CBasis[0]->GetPointsKey(), &tmpGeom[0],
572  m_base[0]->GetPointsKey(), &tmp[i][0]);
573  break;
574  }
575  case 2:
576  {
578  CBasis[0]->GetPointsKey(), CBasis[1]->GetPointsKey(),
579  &tmpGeom[0], m_base[0]->GetPointsKey(),
580  m_base[1]->GetPointsKey(), &tmp[i][0]);
581  break;
582  }
583  case 3:
584  {
586  CBasis[0]->GetPointsKey(), CBasis[1]->GetPointsKey(),
587  CBasis[2]->GetPointsKey(), &tmpGeom[0],
588  m_base[0]->GetPointsKey(), m_base[1]->GetPointsKey(),
589  m_base[2]->GetPointsKey(), &tmp[i][0]);
590  break;
591  }
592  }
593  }
594  }
595 }
596 
598  const Array<OneD, const NekDouble> &direction,
600 {
601  int shapedim = dfdir.size();
602  int coordim = GetCoordim();
603  int nqtot = direction.size() / coordim;
604 
605  for (int j = 0; j < shapedim; j++)
606  {
607  dfdir[j] = Array<OneD, NekDouble>(nqtot, 0.0);
608  for (int k = 0; k < coordim; k++)
609  {
610  if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
611  {
612  Vmath::Vvtvp(nqtot, &df[shapedim * k + j][0], 1,
613  &direction[k * nqtot], 1, &dfdir[j][0], 1,
614  &dfdir[j][0], 1);
615  }
616  else
617  {
618  Vmath::Svtvp(nqtot, df[shapedim * k + j][0],
619  &direction[k * nqtot], 1, &dfdir[j][0], 1,
620  &dfdir[j][0], 1);
621  }
622  }
623  }
624 }
625 
626 // Get Moving frames
628  const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
629 {
630  int coordim = GetCoordim();
631 
632  int nquad0, nquad1, nquad2;
633  int nqtot = 1;
634  switch (shapedim)
635  {
636  case 2:
637  {
638  nquad0 = m_base[0]->GetNumPoints();
639  nquad1 = m_base[1]->GetNumPoints();
640  nqtot = nquad0 * nquad1;
641  break;
642  }
643  case 3:
644  {
645  nquad0 = m_base[0]->GetNumPoints();
646  nquad1 = m_base[1]->GetNumPoints();
647  nquad2 = m_base[2]->GetNumPoints();
648  nqtot = nquad0 * nquad1 * nquad2;
649  break;
650  }
651  default:
652  break;
653  }
654 
655  StdRegions::VarCoeffType MMFCoeffs[15] = {
664 
665  Array<OneD, NekDouble> MF(coordim * nqtot);
666  Array<OneD, NekDouble> tmp(nqtot);
667  for (int k = 0; k < coordim; k++)
668  {
669  StdRegions::VarCoeffMap::const_iterator MFdir =
670  varcoeffs.find(MMFCoeffs[5 * dir + k]);
671  tmp = MFdir->second;
672 
673  Vmath::Vcopy(nqtot, &tmp[0], 1, &MF[k * nqtot], 1);
674  }
675 
676  return MF;
677 }
678 
679 // Get magnitude of MF
681  const int dir, const StdRegions::VarCoeffMap &varcoeffs)
682 {
683  int indxDiv = 3;
684 
685  StdRegions::VarCoeffType MMFCoeffs[15] = {
694 
695  StdRegions::VarCoeffMap::const_iterator MFdir =
696  varcoeffs.find(MMFCoeffs[5 * dir + indxDiv]);
697  Array<OneD, NekDouble> MFDiv = MFdir->second;
698 
699  return MFDiv;
700 }
701 
702 // Get magnitude of MF
704  const int dir, const StdRegions::VarCoeffMap &varcoeffs)
705 {
706  int indxMag = 4;
707 
708  StdRegions::VarCoeffType MMFCoeffs[15] = {
717 
718  StdRegions::VarCoeffMap::const_iterator MFdir =
719  varcoeffs.find(MMFCoeffs[5 * dir + indxMag]);
720  Array<OneD, NekDouble> MFmag = MFdir->second;
721 
722  return MFmag;
723 }
724 
726  const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
727 {
728  boost::ignore_unused(r_bnd, matrixType);
729  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
730  return NullDNekMatSharedPtr;
731 }
732 
734  const DNekScalMatSharedPtr &r_bnd)
735 {
736  boost::ignore_unused(r_bnd);
737  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
738  return NullDNekMatSharedPtr;
739 }
740 
742  const NekDouble *data, const std::vector<unsigned int> &nummodes,
743  const int nmodes_offset, NekDouble *coeffs,
744  std::vector<LibUtilities::BasisType> &fromType)
745 {
746  boost::ignore_unused(data, nummodes, nmodes_offset, coeffs, fromType);
747  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
748 }
749 
751  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
754 {
755  boost::ignore_unused(edge, EdgeExp, Fx, Fy, outarray);
756  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
757 }
758 
760  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
762 {
763  boost::ignore_unused(edge, EdgeExp, Fn, outarray);
764  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
765 }
766 
768  const int face, const std::shared_ptr<Expansion> &FaceExp,
770 {
771  boost::ignore_unused(face, FaceExp, Fn, outarray);
772  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
773 }
774 
775 void Expansion::v_DGDeriv(const int dir,
776  const Array<OneD, const NekDouble> &inarray,
779  Array<OneD, NekDouble> &outarray)
780 {
781  boost::ignore_unused(dir, inarray, EdgeExp, coeffs, outarray);
782  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
783 }
784 
786  const Array<OneD, Array<OneD, NekDouble>> &vec)
787 {
788  boost::ignore_unused(vec);
790  "This function is only valid for "
791  "shape expansion in LocalRegions, not parant class");
792  return 0.0;
793 }
794 
796  Array<OneD, Array<OneD, NekDouble>> &factors,
797  Array<OneD, Array<OneD, NekDouble>> &d0factors,
798  Array<OneD, Array<OneD, NekDouble>> &d1factors)
799 {
800  boost::ignore_unused(factors, d0factors, d1factors);
802  "This function is only valid for "
803  "shape expansion in LocalRegions, not parant class");
804 }
805 
807 {
808  boost::ignore_unused(trace);
809  return StdRegions::eForwards;
810 }
811 
814  Array<OneD, NekDouble> &outarray)
815 {
816  boost::ignore_unused(dir, inarray, outarray);
817  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
818 }
819 
820 void Expansion::v_GetTraceQFactors(const int trace,
821  Array<OneD, NekDouble> &outarray)
822 {
823  boost::ignore_unused(trace, outarray);
825  "Method does not exist for this shape or library");
826 }
827 
829  const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp,
830  const Array<OneD, const NekDouble> &inarray,
832 {
833  boost::ignore_unused(trace, TraceExp, inarray, outarray, orient);
835  "Method does not exist for this shape or library");
836 }
837 
838 void Expansion::v_GetTracePhysMap(const int edge, Array<OneD, int> &outarray)
839 {
840  boost::ignore_unused(edge, outarray);
842  "Method does not exist for this shape or library");
843 }
844 
846  Array<OneD, int> &idmap, const int nq0,
847  const int nq1)
848 {
849  boost::ignore_unused(orient, idmap, nq0, nq1);
851  "Method does not exist for this shape or library");
852 }
853 
854 void Expansion::v_GenTraceExp(const int traceid, ExpansionSharedPtr &exp)
855 {
856  boost::ignore_unused(traceid, exp);
858  "Method does not exist for this shape or library");
859 }
860 
862 {
863  boost::ignore_unused(id);
864  ASSERTL0(false, "Cannot compute trace normal for this expansion.");
865 }
866 
868 {
869  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
870  return NullNekDouble1DArray;
871 }
872 
874 {
875  boost::ignore_unused(normal);
876  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
877 }
878 
879 void Expansion::v_SetUpPhysNormals(const int edge)
880 {
881  boost::ignore_unused(edge);
882  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
883 }
884 
886  const int trace, const Array<OneD, const NekDouble> &primCoeffs,
887  DNekMatSharedPtr &inoutmat)
888 {
889  boost::ignore_unused(trace, primCoeffs, inoutmat);
890  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
891 }
892 
894  const int traceid, const Array<OneD, const NekDouble> &primCoeffs,
895  const Array<OneD, NekDouble> &incoeffs, Array<OneD, NekDouble> &coeffs)
896 {
897  boost::ignore_unused(traceid, primCoeffs, incoeffs, coeffs);
898  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
899 }
900 
901 void Expansion::v_TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
902 {
903  boost::ignore_unused(traceid, h, p);
904  NEKERROR(ErrorUtil::efatal, "This method has not been defined");
905 }
906 
908  const int nbnd) const
909 {
910  auto x = m_elmtBndNormDirElmtLen.find(nbnd);
912  "m_elmtBndNormDirElmtLen normal not computed.");
913  return x->second;
914 }
915 
917  const int dir, const Array<OneD, const NekDouble> &inarray,
918  Array<OneD, Array<OneD, NekDouble>> &outarray)
919 {
920  boost::ignore_unused(dir, inarray, outarray);
921  NEKERROR(ErrorUtil::efatal, "v_AlignVectorToCollapsedDir is not defined");
922 }
923 } // namespace LocalRegions
924 } // namespace Nektar
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
Definition: ErrorUtil.hpp:272
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:15
std::map< int, NormalVector > m_traceNormals
Definition: Expansion.h:275
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:285
virtual void v_SetPhysNormals(Array< OneD, const NekDouble > &normal)
Definition: Expansion.cpp:873
SpatialDomains::GeometrySharedPtr GetGeom() const
Definition: Expansion.cpp:166
virtual void v_TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
Definition: Expansion.cpp:901
IndexMapValuesSharedPtr CreateIndexMap(const IndexMapKey &ikey)
Definition: Expansion.cpp:180
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:272
Array< OneD, NekDouble > v_GetMF(const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:627
void GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.h:208
const SpatialDomains::GeomFactorsSharedPtr & GetMetricInfo() const
Definition: Expansion.cpp:245
DNekMatSharedPtr BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:100
void ComputeGmatcdotMF(const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble >> &dfdir)
Definition: Expansion.cpp:597
virtual void v_DivideByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:428
virtual void v_GenTraceExp(const int traceid, ExpansionSharedPtr &exp)
Definition: Expansion.cpp:854
NekDouble VectorFlux(const Array< OneD, Array< OneD, NekDouble >> &vec)
Definition: Expansion.cpp:144
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:893
virtual void v_AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Expansion.cpp:916
virtual void v_GetTraceQFactors(const int trace, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:820
virtual void v_ComputeTraceNormal(const int id)
Definition: Expansion.cpp:861
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
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:775
virtual DNekMatSharedPtr v_BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:733
virtual void v_SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:812
virtual DNekScalMatSharedPtr v_GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:264
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:767
void NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble >> &factors, Array< OneD, Array< OneD, NekDouble >> &d0factors, Array< OneD, Array< OneD, NekDouble >> &d1factors)
Definition: Expansion.cpp:149
Array< OneD, NekDouble > v_GetMFMag(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:703
virtual DNekMatSharedPtr v_BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:725
const Array< OneD, const NekDouble > & GetElmtBndNormDirElmtLen(const int nbnd) const
Definition: Expansion.cpp:907
virtual void v_AddRobinMassMatrix(const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
Definition: Expansion.cpp:885
std::map< int, ExpansionWeakPtr > m_traceExp
Definition: Expansion.h:271
Array< OneD, NekDouble > v_GetMFDiv(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:680
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:113
virtual NekDouble v_VectorFlux(const Array< OneD, Array< OneD, NekDouble >> &vec)
Definition: Expansion.cpp:785
void AddFaceNormBoundaryInt(const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:128
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:135
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals(void)
Definition: Expansion.cpp:867
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:273
virtual StdRegions::Orientation v_GetTraceOrient(int trace)
Definition: Expansion.cpp:806
void StdDerivBaseOnTraceMat(Array< OneD, DNekMatSharedPtr > &DerivMat)
Definition: Expansion.cpp:468
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:105
virtual void v_ReOrientTracePhysMap(const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
Definition: Expansion.cpp:845
virtual void v_MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:413
virtual void v_ComputeLaplacianMetric()
Definition: Expansion.h:301
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
Definition: Expansion.cpp:524
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:741
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:750
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:88
virtual void v_SetUpPhysNormals(const int id)
Definition: Expansion.cpp:879
virtual void v_GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.cpp:838
const NormalVector & GetTraceNormal(const int id)
Definition: Expansion.cpp:250
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:795
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:828
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:47
DNekMatSharedPtr BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:94
StdRegions::Orientation GetIndexOrientation() const
Definition: IndexMapKey.h:96
IndexMapType GetIndexMapType() const
Definition: IndexMapKey.h:91
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:677
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:140
int GetTraceNumPoints(const int i) const
This function returns the number of quadrature points belonging to the i-th trace.
Definition: StdExpansion.h:289
const LibUtilities::PointsKeyVector GetPointsKeys() const
DNekBlkMatSharedPtr GetStdStaticCondMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:616
int GetNtraces() const
Returns the number of trace elements connected to this element.
Definition: StdExpansion.h:359
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:692
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:375
void GetInteriorMap(Array< OneD, unsigned int > &outarray)
Definition: StdExpansion.h:682
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:432
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:883
Array< OneD, LibUtilities::BasisSharedPtr > m_base
MatrixType GetMatrixType() const
Definition: StdMatrixKey.h:85
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:52
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:167
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:106
std::vector< PointsKey > PointsKeyVector
Definition: Points.h:250
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
std::shared_ptr< IndexMapValues > IndexMapValuesSharedPtr
Definition: IndexMapKey.h:128
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:62
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:53
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:240
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:282
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
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.cpp:209
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.cpp:622
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.cpp:574
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.cpp:284
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255