Nektar++
StdExpansion.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File StdExpansion.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: Definition of methods in class StdExpansion which is
32 // the base class to all expansion shapes
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <boost/core/ignore_unused.hpp>
37 
38 #include <LibUtilities/Foundations/ManagerAccess.h> // for BasisManager, etc
40 
41 namespace Nektar
42 {
43 namespace StdRegions
44 {
45 StdExpansion::StdExpansion(void) : m_elmt_id(0), m_ncoeffs(0)
46 {
47 }
48 
49 StdExpansion::StdExpansion(const int numcoeffs, const int numbases,
50  const LibUtilities::BasisKey &Ba,
51  const LibUtilities::BasisKey &Bb,
52  const LibUtilities::BasisKey &Bc)
53  : m_base(numbases), m_elmt_id(0), m_ncoeffs(numcoeffs),
54  m_stdMatrixManager(std::bind(&StdExpansion::CreateStdMatrix, this,
55  std::placeholders::_1),
56  std::string("StdExpansionStdMatrix")),
57  m_stdStaticCondMatrixManager(
58  std::bind(&StdExpansion::CreateStdStaticCondMatrix, this,
59  std::placeholders::_1),
60  std::string("StdExpansionStdStaticCondMatrix"))
61 {
62  switch (m_base.size())
63  {
64  case 3:
66  "NULL Basis attempting to be used.");
68  /* Falls through. */
69  case 2:
71  "NULL Basis attempting to be used.");
73  /* Falls through. */
74  case 1:
76  "NULL Basis attempting to be used.");
78  break;
79  default:
80  break;
81  // ASSERTL0(false, "numbases incorrectly specified");
82  };
83 
84 } // end constructor
85 
87  : std::enable_shared_from_this<StdExpansion>(T), m_base(T.m_base),
88  m_elmt_id(T.m_elmt_id), m_ncoeffs(T.m_ncoeffs),
89  m_stdMatrixManager(T.m_stdMatrixManager),
90  m_stdStaticCondMatrixManager(T.m_stdStaticCondMatrixManager)
91 {
92 }
93 
95 {
96 }
97 
100 {
101  NekDouble val;
102  int ntot = GetTotPoints();
103  Array<OneD, NekDouble> wsp(ntot);
104 
105  if (sol == NullNekDouble1DArray)
106  {
107  Vmath::Vabs(ntot, phys, 1, wsp, 1);
108  }
109  else
110  {
111  Vmath::Vsub(ntot, sol, 1, phys, 1, wsp, 1);
112  Vmath::Vabs(ntot, wsp, 1, wsp, 1);
113  }
114 
115  val = Vmath::Vamax(ntot, wsp, 1);
116  return val;
117 }
118 
120  const Array<OneD, const NekDouble> &sol)
121 {
122  NekDouble val;
123  int ntot = GetTotPoints();
124  Array<OneD, NekDouble> wsp(ntot);
125 
126  if (sol.size() == 0)
127  {
128  Vmath::Vmul(ntot, phys, 1, phys, 1, wsp, 1);
129  }
130  else
131  {
132  Vmath::Vsub(ntot, sol, 1, phys, 1, wsp, 1);
133  Vmath::Vmul(ntot, wsp, 1, wsp, 1, wsp, 1);
134  }
135 
136  val = v_Integral(wsp);
137 
138  // if val too small, sqrt returns nan.
140  {
141  return 0.0;
142  }
143  else
144  {
145  return sqrt(val);
146  }
147 }
148 
150  const Array<OneD, const NekDouble> &sol)
151 {
152  int i;
153  NekDouble val;
154  int ntot = GetTotPoints();
155  int coordim = v_GetCoordim();
156  Array<OneD, NekDouble> wsp(3 * ntot);
157  Array<OneD, NekDouble> wsp_deriv = wsp + ntot;
158  Array<OneD, NekDouble> sum = wsp_deriv + ntot;
159 
160  if (sol == NullNekDouble1DArray)
161  {
162  Vmath::Vcopy(ntot, phys, 1, wsp, 1);
163  Vmath::Vmul(ntot, phys, 1, phys, 1, sum, 1);
164  }
165  else
166  {
167  Vmath::Vsub(ntot, sol, 1, phys, 1, wsp, 1);
168  Vmath::Vmul(ntot, wsp, 1, wsp, 1, sum, 1);
169  }
170 
171  for (i = 0; i < coordim; ++i)
172  {
173  v_PhysDeriv(i, wsp, wsp_deriv);
174  Vmath::Vvtvp(ntot, wsp_deriv, 1, wsp_deriv, 1, sum, 1, sum, 1);
175  }
176 
177  val = sqrt(v_Integral(sum));
178 
179  return val;
180 }
181 
183  const StdMatrixKey &mkey)
184 {
185  DNekBlkMatSharedPtr returnval;
186 
187  DNekMatSharedPtr mat = GetStdMatrix(mkey);
188  int nbdry = NumBndryCoeffs(); // also checks to see if this is a boundary
189  // interior decomposed expansion
190  int nint = m_ncoeffs - nbdry;
196 
197  int i, j;
198 
199  Array<OneD, unsigned int> bmap(nbdry);
200  Array<OneD, unsigned int> imap(nint);
201  GetBoundaryMap(bmap);
202  GetInteriorMap(imap);
203 
204  for (i = 0; i < nbdry; ++i)
205  {
206  for (j = 0; j < nbdry; ++j)
207  {
208  (*A)(i, j) = (*mat)(bmap[i], bmap[j]);
209  }
210 
211  for (j = 0; j < nint; ++j)
212  {
213  (*B)(i, j) = (*mat)(bmap[i], imap[j]);
214  }
215  }
216 
217  for (i = 0; i < nint; ++i)
218  {
219  for (j = 0; j < nbdry; ++j)
220  {
221  (*C)(i, j) = (*mat)(imap[i], bmap[j]);
222  }
223 
224  for (j = 0; j < nint; ++j)
225  {
226  (*D)(i, j) = (*mat)(imap[i], imap[j]);
227  }
228  }
229 
230  // Calculate static condensed system
231  if (nint)
232  {
233  D->Invert();
234  (*B) = (*B) * (*D);
235  (*A) = (*A) - (*B) * (*C);
236  }
237 
238  // set up block matrix system
239  Array<OneD, unsigned int> exp_size(2);
240  exp_size[0] = nbdry;
241  exp_size[1] = nint;
242  returnval =
244 
245  returnval->SetBlock(0, 0, A);
246  returnval->SetBlock(0, 1, B);
247  returnval->SetBlock(1, 0, C);
248  returnval->SetBlock(1, 1, D);
249 
250  return returnval;
251 }
252 
254 {
255  int i;
256  DNekMatSharedPtr returnval;
257 
258  switch (mkey.GetMatrixType())
259  {
260  case eInvMass:
261  {
262  StdMatrixKey masskey(eMass, mkey.GetShapeType(), *this,
264  mkey.GetNodalPointsType());
265  DNekMatSharedPtr mmat = GetStdMatrix(masskey);
266 
268  *mmat); // Populate standard mass matrix.
269  returnval->Invert();
270  }
271  break;
272  case eInvNBasisTrans:
273  {
274  StdMatrixKey tmpkey(eNBasisTrans, mkey.GetShapeType(), *this,
276  mkey.GetNodalPointsType());
277  DNekMatSharedPtr tmpmat = GetStdMatrix(tmpkey);
279  *tmpmat); // Populate matrix.
280  returnval->Invert();
281  }
282  break;
283  case eBwdMat:
284  {
285  int nq = GetTotPoints();
287  Array<OneD, NekDouble> tmpout(nq);
288 
289  returnval =
291  Array<OneD, NekDouble> Bwd_data = returnval->GetPtr();
292 
294  this->DetShapeType(), *this);
295  DNekMatSharedPtr MatBwdTrans = GetStdMatrix(matkey);
296  Array<OneD, NekDouble> BwdTrans_data = MatBwdTrans->GetPtr();
297 
298  for (int i = 0; i < m_ncoeffs; ++i)
299  {
300  Array<OneD, NekDouble> tmpinn = BwdTrans_data + nq * i;
301  Array<OneD, NekDouble> tmpout = Bwd_data + i;
302 
303  Vmath::Vcopy(nq, tmpinn, 1, tmpout, m_ncoeffs);
304  }
305  }
306  break;
307  case eBwdTrans:
308  {
309  int nq = GetTotPoints();
311  Array<OneD, NekDouble> tmpout(nq);
312 
313  returnval =
315 
316  for (int i = 0; i < m_ncoeffs; ++i)
317  {
318  Vmath::Zero(m_ncoeffs, tmpin, 1);
319  tmpin[i] = 1.0;
320 
321  BwdTrans_SumFac(tmpin, tmpout);
322 
323  Vmath::Vcopy(nq, tmpout.get(), 1,
324  returnval->GetRawPtr() + i * nq, 1);
325  }
326  }
327  break;
328  case eIProductWRTBase:
329  {
330  int nq = GetTotPoints();
331  Array<OneD, NekDouble> tmpin(nq);
333 
334  returnval =
336 
337  for (i = 0; i < nq; ++i)
338  {
339  Vmath::Zero(nq, tmpin, 1);
340  tmpin[i] = 1.0;
341 
342  IProductWRTBase_SumFac(tmpin, tmpout);
343 
344  Vmath::Vcopy(m_ncoeffs, tmpout.get(), 1,
345  returnval->GetRawPtr() + i * m_ncoeffs, 1);
346  }
347  }
348  break;
350  {
351  int nq = GetTotPoints();
352  Array<OneD, NekDouble> tmpin(nq);
354 
355  returnval =
357 
358  for (i = 0; i < nq; ++i)
359  {
360  Vmath::Zero(nq, tmpin, 1);
361  tmpin[i] = 1.0;
362 
363  IProductWRTDerivBase_SumFac(0, tmpin, tmpout);
364 
365  Vmath::Vcopy(m_ncoeffs, tmpout.get(), 1,
366  returnval->GetRawPtr() + i * m_ncoeffs, 1);
367  }
368  }
369  break;
371  {
372  int nq = GetTotPoints();
373  Array<OneD, NekDouble> tmpin(nq);
375 
376  returnval =
378 
379  for (i = 0; i < nq; ++i)
380  {
381  Vmath::Zero(nq, tmpin, 1);
382  tmpin[i] = 1.0;
383 
384  IProductWRTDerivBase_SumFac(1, tmpin, tmpout);
385 
386  Vmath::Vcopy(m_ncoeffs, tmpout.get(), 1,
387  returnval->GetRawPtr() + i * m_ncoeffs, 1);
388  }
389  }
390  break;
392  {
393  int nq = GetTotPoints();
394  Array<OneD, NekDouble> tmpin(nq);
396 
397  returnval =
399 
400  for (i = 0; i < nq; ++i)
401  {
402  Vmath::Zero(nq, tmpin, 1);
403  tmpin[i] = 1.0;
404 
405  IProductWRTDerivBase_SumFac(2, tmpin, tmpout);
406 
407  Vmath::Vcopy(m_ncoeffs, tmpout.get(), 1,
408  returnval->GetRawPtr() + i * m_ncoeffs, 1);
409  }
410  }
411  break;
412  case eDerivBase0:
413  {
414  int nq = GetTotPoints();
415  returnval =
417  GenStdMatBwdDeriv(0, returnval);
418  }
419  break;
420  case eDerivBase1:
421  {
422  int nq = GetTotPoints();
423  returnval =
425  GenStdMatBwdDeriv(1, returnval);
426  }
427  break;
428  case eDerivBase2:
429  {
430  int nq = GetTotPoints();
431  returnval =
433  GenStdMatBwdDeriv(2, returnval);
434  }
435  break;
436  case eEquiSpacedToCoeffs:
437  {
438  // check to see if equispaced basis
439  int nummodes = m_base[0]->GetNumModes();
440  bool equispaced = true;
441  for (int i = 1; i < m_base.size(); ++i)
442  {
443  if (m_base[i]->GetNumModes() != nummodes)
444  {
445  equispaced = false;
446  }
447  }
448 
449  ASSERTL0(equispaced,
450  "Currently need to have same num modes in all "
451  "directionmodes to use EquiSpacedToCoeff method");
452 
453  int ntot = GetTotPoints();
454  Array<OneD, NekDouble> qmode(ntot);
456 
457  returnval =
459  for (int i = 0; i < m_ncoeffs; ++i)
460  {
461  // Get mode at quadrature points
462  FillMode(i, qmode);
463 
464  // interpolate to equi spaced
465  PhysInterpToSimplexEquiSpaced(qmode, emode, nummodes);
466 
467  // fill matrix
468  Vmath::Vcopy(m_ncoeffs, &emode[0], 1,
469  returnval->GetRawPtr() + i * m_ncoeffs, 1);
470  }
471  // invert matrix
472  returnval->Invert();
473  }
474  break;
475  case eMass:
476  case eHelmholtz:
477  case eLaplacian:
478  case eLaplacian00:
479  case eLaplacian01:
480  case eLaplacian02:
481  case eLaplacian11:
482  case eLaplacian12:
483  case eLaplacian22:
484  case eWeakDeriv0:
485  case eWeakDeriv1:
486  case eWeakDeriv2:
488  case eMassLevelCurvature:
491  {
493  returnval =
495  DNekMat &Mat = *returnval;
496 
497  for (i = 0; i < m_ncoeffs; ++i)
498  {
499  Vmath::Zero(m_ncoeffs, tmp, 1);
500  tmp[i] = 1.0;
501 
502  GeneralMatrixOp_MatFree(tmp, tmp, mkey);
503 
504  Vmath::Vcopy(m_ncoeffs, &tmp[0], 1,
505  &(Mat.GetPtr())[0] + i * m_ncoeffs, 1);
506  }
507  }
508  break;
509  default:
510  {
512  "This type of matrix, " +
513  static_cast<std::string>(
514  MatrixTypeMap[mkey.GetMatrixType()]) +
515  ", can not be created using a general approach");
516  }
517  break;
518  }
519 
520  return returnval;
521 }
522 
524  Array<OneD, NekDouble> &outarray,
525  const StdMatrixKey &mkey)
526 {
527  switch (mkey.GetMatrixType())
528  {
529  case eMass:
530  MassMatrixOp(inarray, outarray, mkey);
531  break;
532  case eWeakDeriv0:
533  WeakDerivMatrixOp(0, inarray, outarray, mkey);
534  break;
535  case eWeakDeriv1:
536  WeakDerivMatrixOp(1, inarray, outarray, mkey);
537  break;
538  case eWeakDeriv2:
539  WeakDerivMatrixOp(2, inarray, outarray, mkey);
540  break;
542  WeakDirectionalDerivMatrixOp(inarray, outarray, mkey);
543  break;
544  case eMassLevelCurvature:
545  MassLevelCurvatureMatrixOp(inarray, outarray, mkey);
546  break;
548  LinearAdvectionDiffusionReactionMatrixOp(inarray, outarray, mkey,
549  false);
550  break;
552  LinearAdvectionDiffusionReactionMatrixOp(inarray, outarray, mkey);
553  break;
554  case eLaplacian:
555  LaplacianMatrixOp(inarray, outarray, mkey);
556  break;
557  case eLaplacian00:
558  LaplacianMatrixOp(0, 0, inarray, outarray, mkey);
559  break;
560  case eLaplacian01:
561  LaplacianMatrixOp(0, 1, inarray, outarray, mkey);
562  break;
563  case eLaplacian02:
564  LaplacianMatrixOp(0, 2, inarray, outarray, mkey);
565  break;
566  case eLaplacian10:
567  LaplacianMatrixOp(1, 0, inarray, outarray, mkey);
568  break;
569  case eLaplacian11:
570  LaplacianMatrixOp(1, 1, inarray, outarray, mkey);
571  break;
572  case eLaplacian12:
573  LaplacianMatrixOp(1, 2, inarray, outarray, mkey);
574  break;
575  case eLaplacian20:
576  LaplacianMatrixOp(2, 0, inarray, outarray, mkey);
577  break;
578  case eLaplacian21:
579  LaplacianMatrixOp(2, 1, inarray, outarray, mkey);
580  break;
581  case eLaplacian22:
582  LaplacianMatrixOp(2, 2, inarray, outarray, mkey);
583  break;
584  case eHelmholtz:
585  HelmholtzMatrixOp(inarray, outarray, mkey);
586  break;
587  default:
589  "This matrix does not have an operator");
590  break;
591  }
592 }
593 
595  const Array<OneD, const NekDouble> &inarray,
596  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
597 {
598  switch (mkey.GetMatrixType())
599  {
600  case eMass:
601  MassMatrixOp_MatFree(inarray, outarray, mkey);
602  break;
603  case eWeakDeriv0:
604  WeakDerivMatrixOp_MatFree(0, inarray, outarray, mkey);
605  break;
606  case eWeakDeriv1:
607  WeakDerivMatrixOp_MatFree(1, inarray, outarray, mkey);
608  break;
609  case eWeakDeriv2:
610  WeakDerivMatrixOp_MatFree(2, inarray, outarray, mkey);
611  break;
613  WeakDirectionalDerivMatrixOp_MatFree(inarray, outarray, mkey);
614  break;
615  case eMassLevelCurvature:
616  MassLevelCurvatureMatrixOp_MatFree(inarray, outarray, mkey);
617  break;
620  mkey, false);
621  break;
624  mkey);
625  break;
626  case eLaplacian:
627  LaplacianMatrixOp_MatFree(inarray, outarray, mkey);
628  break;
629  case eLaplacian00:
630  LaplacianMatrixOp_MatFree(0, 0, inarray, outarray, mkey);
631  break;
632  case eLaplacian01:
633  LaplacianMatrixOp_MatFree(0, 1, inarray, outarray, mkey);
634  break;
635  case eLaplacian02:
636  LaplacianMatrixOp_MatFree(0, 2, inarray, outarray, mkey);
637  break;
638  case eLaplacian10:
639  LaplacianMatrixOp_MatFree(1, 0, inarray, outarray, mkey);
640  break;
641  case eLaplacian11:
642  LaplacianMatrixOp_MatFree(1, 1, inarray, outarray, mkey);
643  break;
644  case eLaplacian12:
645  LaplacianMatrixOp_MatFree(1, 2, inarray, outarray, mkey);
646  break;
647  case eLaplacian20:
648  LaplacianMatrixOp_MatFree(2, 0, inarray, outarray, mkey);
649  break;
650  case eLaplacian21:
651  LaplacianMatrixOp_MatFree(2, 1, inarray, outarray, mkey);
652  break;
653  case eLaplacian22:
654  LaplacianMatrixOp_MatFree(2, 2, inarray, outarray, mkey);
655  break;
656  case eHelmholtz:
657  HelmholtzMatrixOp_MatFree(inarray, outarray, mkey);
658  break;
659  default:
661  "This matrix does not have an operator");
662  break;
663  }
664 }
665 
667  const Array<OneD, const NekDouble> &inarray,
668  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
669 {
670  int nq = GetTotPoints();
671  Array<OneD, NekDouble> tmp(nq);
672 
673  v_BwdTrans(inarray, tmp);
674 
675  if (mkey.HasVarCoeff(eVarCoeffMass))
676  {
677  Vmath::Vmul(nq, mkey.GetVarCoeff(eVarCoeffMass), 1, tmp, 1, tmp, 1);
678  }
679 
680  v_IProductWRTBase(tmp, outarray);
681 }
682 
684  const int k1, const int k2, const Array<OneD, const NekDouble> &inarray,
685  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
686 {
687  ASSERTL1(k1 >= 0 && k1 < GetCoordim(), "invalid first argument");
688  ASSERTL1(k2 >= 0 && k2 < GetCoordim(), "invalid second argument");
689 
690  int nq = GetTotPoints();
691  Array<OneD, NekDouble> tmp(nq);
692  Array<OneD, NekDouble> dtmp(nq);
693  VarCoeffType varcoefftypes[3][3] = {
697 
698  ConstFactorType constcoefftypes[3][3] = {
702 
703  v_BwdTrans(inarray, tmp);
704  v_PhysDeriv(k2, tmp, dtmp);
705  if (mkey.GetNVarCoeff() && (!mkey.ConstFactorExists(eFactorSVVDiffCoeff)))
706  {
707  if (k1 == k2)
708  {
709  // By default, k1 == k2 has \sigma = 1 (diagonal entries)
710  if (mkey.HasVarCoeff(varcoefftypes[k1][k1]))
711  {
712  Vmath::Vmul(nq, mkey.GetVarCoeff(varcoefftypes[k1][k1]), 1,
713  dtmp, 1, dtmp, 1);
714  }
715  v_IProductWRTDerivBase_SumFac(k1, dtmp, outarray);
716  }
717  else
718  {
719  // By default, k1 != k2 has \sigma = 0 (off-diagonal entries)
720  if (mkey.HasVarCoeff(varcoefftypes[k1][k2]))
721  {
722  Vmath::Vmul(nq, mkey.GetVarCoeff(varcoefftypes[k1][k2]), 1,
723  dtmp, 1, dtmp, 1);
724  v_IProductWRTDerivBase_SumFac(k1, dtmp, outarray);
725  }
726  else
727  {
728  Vmath::Zero(GetNcoeffs(), outarray, 1);
729  }
730  }
731  }
732  else if (mkey.ConstFactorExists(eFactorCoeffD00) &&
734  {
735  if (k1 == k2)
736  {
737  // By default, k1 == k2 has \sigma = 1 (diagonal entries)
738  if (mkey.ConstFactorExists(constcoefftypes[k1][k1]))
739  {
740  Vmath::Smul(nq, mkey.GetConstFactor(constcoefftypes[k1][k1]),
741  dtmp, 1, dtmp, 1);
742  }
743  v_IProductWRTDerivBase(k1, dtmp, outarray);
744  }
745  else
746  {
747  // By default, k1 != k2 has \sigma = 0 (off-diagonal entries)
748  if (mkey.ConstFactorExists(constcoefftypes[k1][k2]))
749  {
750  Vmath::Smul(nq, mkey.GetConstFactor(constcoefftypes[k1][k2]),
751  dtmp, 1, dtmp, 1);
752  v_IProductWRTDerivBase(k1, dtmp, outarray);
753  }
754  else
755  {
756  Vmath::Zero(GetNcoeffs(), outarray, 1);
757  }
758  }
759  }
760  else
761  {
762  // Multiply by svv tensor
764  {
765  Vmath::Vcopy(nq, dtmp, 1, tmp, 1);
766  SVVLaplacianFilter(dtmp, mkey);
767  Vmath::Vadd(nq, tmp, 1, dtmp, 1, dtmp, 1);
768  }
769  v_IProductWRTDerivBase(k1, dtmp, outarray);
770  }
771 }
772 
774  const Array<OneD, const NekDouble> &inarray,
775  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
776 {
777  const int dim = GetCoordim();
778 
779  int i, j;
780 
782  Array<OneD, NekDouble> store2(m_ncoeffs, 0.0);
783 
784  if ((mkey.GetNVarCoeff() == 0 &&
787  {
788  // just call diagonal matrix form of laplcian operator
789  for (i = 0; i < dim; ++i)
790  {
791  LaplacianMatrixOp(i, i, inarray, store, mkey);
792  Vmath::Vadd(m_ncoeffs, store, 1, store2, 1, store2, 1);
793  }
794  }
795  else
796  {
797  const MatrixType mtype[3][3] = {
801  StdMatrixKeySharedPtr mkeyij;
802 
803  for (i = 0; i < dim; i++)
804  {
805  for (j = 0; j < dim; j++)
806  {
808  mkey, mtype[i][j]);
809  LaplacianMatrixOp(i, j, inarray, store, *mkeyij);
810  Vmath::Vadd(m_ncoeffs, store, 1, store2, 1, store2, 1);
811  }
812  }
813  }
814 
815  Vmath::Vcopy(m_ncoeffs, store2.get(), 1, outarray.get(), 1);
816 }
817 
819  const int k1, const Array<OneD, const NekDouble> &inarray,
820  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
821 {
823  int nq = GetTotPoints();
824 
825  v_BwdTrans(inarray, tmp);
826  v_PhysDeriv(k1, tmp, tmp);
827 
829  if (mkey.HasVarCoeff(keys[k1]))
830  {
831  Vmath::Vmul(nq, &(mkey.GetVarCoeff(keys[k1]))[0], 1, &tmp[0], 1,
832  &tmp[0], 1);
833  }
834 
835  v_IProductWRTBase(tmp, outarray);
836 }
837 
839  const Array<OneD, const NekDouble> &inarray,
840  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
841 {
842  int nq = GetTotPoints();
843 
844  Array<OneD, NekDouble> tmp(nq), Dtmp(nq);
845  Array<OneD, NekDouble> Mtmp(nq), Mout(m_ncoeffs);
846 
847  v_BwdTrans(inarray, tmp);
849 
850  v_IProductWRTBase(Dtmp, outarray);
851 
852  // Compte M_{div tv}
853  Vmath::Vmul(nq, &(mkey.GetVarCoeff(eVarCoeffMFDiv))[0], 1, &tmp[0], 1,
854  &Mtmp[0], 1);
855 
856  v_IProductWRTBase(Mtmp, Mout);
857 
858  // Add D_tv + M_{div tv}
859  Vmath::Vadd(m_ncoeffs, &Mout[0], 1, &outarray[0], 1, &outarray[0], 1);
860 }
861 
863  const Array<OneD, const NekDouble> &inarray,
864  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
865 {
866  boost::ignore_unused(inarray, outarray, mkey);
867 }
868 
870  const Array<OneD, const NekDouble> &inarray,
871  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey,
872  bool addDiffusionTerm)
873 {
874 
875  int i;
876  int ndir =
877  mkey.GetNVarCoeff(); // assume num.r consts corresponds to directions
878  ASSERTL0(ndir, "Must define at least one advection velocity");
879 
880  NekDouble lambda = mkey.GetConstFactor(eFactorLambda);
881  int totpts = GetTotPoints();
882  Array<OneD, NekDouble> tmp(3 * totpts);
883  Array<OneD, NekDouble> tmp_deriv = tmp + totpts;
884  Array<OneD, NekDouble> tmp_adv = tmp_deriv + totpts;
885 
886  ASSERTL1(ndir <= GetCoordim(),
887  "Number of constants is larger than coordinate dimensions");
888 
889  v_BwdTrans(inarray, tmp);
890 
891  VarCoeffType varcoefftypes[] = {eVarCoeffVelX, eVarCoeffVelY,
892  eVarCoeffVelZ};
893 
894  // calculate u dx + v dy + ..
895  Vmath::Zero(totpts, tmp_adv, 1);
896  for (i = 0; i < ndir; ++i)
897  {
898  v_PhysDeriv(i, tmp, tmp_deriv);
899  Vmath::Vvtvp(totpts, mkey.GetVarCoeff(varcoefftypes[i]), 1, tmp_deriv,
900  1, tmp_adv, 1, tmp_adv, 1);
901  }
902 
903  if (lambda) // add -lambda*u
904  {
905  Vmath::Svtvp(totpts, -lambda, tmp, 1, tmp_adv, 1, tmp_adv, 1);
906  }
907 
908  if (addDiffusionTerm)
909  {
911  StdMatrixKey mkeylap(eLaplacian, DetShapeType(), *this,
912  mkey.GetConstFactors(), mkey.GetVarCoeffs(),
913  mkey.GetNodalPointsType());
914  LaplacianMatrixOp(inarray, lap, mkeylap);
915 
916  v_IProductWRTBase(tmp_adv, outarray);
917  // Lap v - u.grad v + lambda*u
918  // => (grad u, grad v) + u.grad v - lambda*u
919  Vmath::Vadd(m_ncoeffs, lap, 1, outarray, 1, outarray, 1);
920  }
921  else
922  {
923  v_IProductWRTBase(tmp_adv, outarray);
924  }
925 }
926 
928  const Array<OneD, const NekDouble> &inarray,
929  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
930 {
931  NekDouble lambda = mkey.GetConstFactor(eFactorLambda);
933  StdMatrixKey mkeymass(eMass, DetShapeType(), *this);
934  StdMatrixKey mkeylap(eLaplacian, DetShapeType(), *this,
935  mkey.GetConstFactors(), mkey.GetVarCoeffs(),
936  mkey.GetNodalPointsType());
937 
938  MassMatrixOp(inarray, tmp, mkeymass);
939  LaplacianMatrixOp(inarray, outarray, mkeylap);
940 
941  Blas::Daxpy(m_ncoeffs, lambda, tmp, 1, outarray, 1);
942 }
943 
945  Array<OneD, NekDouble> &outarray)
946 {
947  int nq = GetTotPoints();
948  StdMatrixKey bwdtransmatkey(eBwdTrans, DetShapeType(), *this);
949  DNekMatSharedPtr bwdtransmat = GetStdMatrix(bwdtransmatkey);
950 
951  Blas::Dgemv('N', nq, m_ncoeffs, 1.0, bwdtransmat->GetPtr().get(), nq,
952  inarray.get(), 1, 0.0, outarray.get(), 1);
953 }
954 
955 // VIRTUAL INLINE FUNCTIONS FROM HEADER FILE
957  const Array<OneD, const NekDouble> &Lcoord,
958  const Array<OneD, const NekDouble> &physvals)
959 {
960  return v_StdPhysEvaluate(Lcoord, physvals);
961 }
962 
964  const std::vector<unsigned int> &nummodes, int &modes_offset)
965 {
966  boost::ignore_unused(nummodes, modes_offset);
967  NEKERROR(ErrorUtil::efatal, "This function is not defined for this class");
968  return 0;
969 }
970 
973 {
974  boost::ignore_unused(Fx, outarray);
975  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
976 }
977 
981 {
982  boost::ignore_unused(Fx, Fy, outarray);
983  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
984 }
985 
990 {
991  boost::ignore_unused(Fx, Fy, Fz, outarray);
992  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
993 }
994 
996  const Array<OneD, const Array<OneD, NekDouble>> &Fvec,
997  Array<OneD, NekDouble> &outarray)
998 {
999  boost::ignore_unused(Fvec, outarray);
1000  NEKERROR(ErrorUtil::efatal, "This function is not valid for this class");
1001 }
1002 
1004  const LocalRegions::MatrixKey &mkey)
1005 {
1006  boost::ignore_unused(mkey);
1007  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
1009 }
1010 
1012  const LocalRegions::MatrixKey &mkey)
1013 {
1014  boost::ignore_unused(mkey);
1015  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
1016 }
1017 
1020  Array<OneD, NekDouble> &outarray)
1021 {
1022  boost::ignore_unused(dir, inarray, outarray);
1023  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
1024 }
1025 
1028 {
1029  boost::ignore_unused(coeffs, dir);
1030  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
1031 }
1032 
1034  const Array<OneD, const NekDouble> &Lcoord,
1035  const Array<OneD, const NekDouble> &physvals)
1036 
1037 {
1038  boost::ignore_unused(Lcoord, physvals);
1039  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
1040  return 0;
1041 }
1042 
1045 {
1046  boost::ignore_unused(xi, eta);
1047  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
1048 }
1049 
1052 {
1053  boost::ignore_unused(eta, xi);
1054  NEKERROR(ErrorUtil::efatal, "This function is not defined for this shape");
1055 }
1056 
1058 {
1059  ASSERTL0(false, "This function is needs defining for this shape");
1060  return 0;
1061 }
1062 
1064 {
1065  ASSERTL0(false, "This function is needs defining for this shape");
1066  return 0;
1067 }
1068 
1070 {
1071  ASSERTL0(false, "This function is needs defining for this shape");
1072  return 0;
1073 }
1074 
1076  const int k) const
1077 {
1078  boost::ignore_unused(i, k);
1079  ASSERTL0(false, "This function is not valid or not defined");
1081 }
1082 
1084 {
1085  boost::ignore_unused(i);
1086  ASSERTL0(false, "This function is not valid or not defined");
1087  return 0;
1088 }
1089 
1090 int StdExpansion::v_GetTraceNcoeffs(const int i) const
1091 {
1092  boost::ignore_unused(i);
1093  ASSERTL0(false, "This function is not valid or not defined");
1094  return 0;
1095 }
1096 
1098 {
1099  boost::ignore_unused(i);
1100  ASSERTL0(false, "This function is not valid or not defined");
1101  return 0;
1102 }
1103 
1105 {
1106  ASSERTL0(false, "This function is not valid or not defined");
1107  return 0;
1108 }
1109 
1111  const int j) const
1112 {
1113  boost::ignore_unused(i, j);
1114  ASSERTL0(false, "This function is not valid or not defined");
1116 }
1117 
1119 {
1120  ASSERTL0(false, "This function is not valid or not defined");
1121 
1123 }
1124 
1126 {
1127  ASSERTL0(false, "This expansion does not have a shape type defined");
1129 }
1130 
1131 std::shared_ptr<StdExpansion> StdExpansion::v_GetStdExp(void) const
1132 {
1133  ASSERTL0(false, "This method is not defined for this expansion");
1134  StdExpansionSharedPtr returnval;
1135  return returnval;
1136 }
1137 
1138 std::shared_ptr<StdExpansion> StdExpansion::v_GetLinStdExp(void) const
1139 {
1140  ASSERTL0(false, "This method is not defined for this expansion");
1141  StdExpansionSharedPtr returnval;
1142  return returnval;
1143 }
1144 
1146 {
1147  ASSERTL0(false, "This function is not valid or not defined");
1148  return 0;
1149 }
1150 
1152 {
1153  ASSERTL0(false, "This function has not been defined for this expansion");
1154  return false;
1155 }
1156 
1158 {
1159  return false;
1160 }
1161 
1163  const int dir, const Array<OneD, const NekDouble> &inarray,
1164  Array<OneD, NekDouble> &outarray)
1165 {
1166  boost::ignore_unused(dir, inarray, outarray);
1167  NEKERROR(ErrorUtil::efatal, "This method has not been defined");
1168 }
1169 
1170 /**
1171  *
1172  */
1174  const Array<OneD, const NekDouble> &direction,
1175  const Array<OneD, const NekDouble> &inarray,
1176  Array<OneD, NekDouble> &outarray)
1177 {
1178  boost::ignore_unused(direction, inarray, outarray);
1179  NEKERROR(ErrorUtil::efatal, "This method has not been defined");
1180 }
1181 
1182 /**
1183  *
1184  */
1186  const Array<OneD, const NekDouble> &inarray,
1187  Array<OneD, NekDouble> &outarray)
1188 {
1189  boost::ignore_unused(inarray, outarray);
1190  NEKERROR(ErrorUtil::efatal, "This method has not been defined");
1191 }
1192 
1193 /**
1194  * @brief Integrates the specified function over the domain.
1195  * @see StdRegions#StdExpansion#Integral.
1196  */
1198 {
1199  boost::ignore_unused(inarray);
1200  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1201  "local expansions");
1202  return 0;
1203 }
1204 
1205 /**
1206  * @brief Calculate the derivative of the physical points
1207  * @see StdRegions#StdExpansion#PhysDeriv
1208  */
1210  Array<OneD, NekDouble> &out_d1,
1211  Array<OneD, NekDouble> &out_d2,
1212  Array<OneD, NekDouble> &out_d3)
1213 {
1214  boost::ignore_unused(inarray, out_d1, out_d2, out_d3);
1215  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1216  "local expansions");
1217 }
1218 
1220  Array<OneD, NekDouble> &out_ds)
1221 {
1222  boost::ignore_unused(inarray, out_ds);
1223  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1224  "local expansions");
1225 }
1227  Array<OneD, NekDouble> &out_dn)
1228 {
1229  boost::ignore_unused(inarray, out_dn);
1230  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1231  "local expansions");
1232 }
1233 
1234 /**
1235  * @brief Calculate the derivative of the physical points in a
1236  * given direction
1237  * @see StdRegions#StdExpansion#PhysDeriv
1238  */
1239 void StdExpansion::v_PhysDeriv(const int dir,
1240  const Array<OneD, const NekDouble> &inarray,
1241  Array<OneD, NekDouble> &out_d0)
1242 
1243 {
1244  boost::ignore_unused(dir, inarray, out_d0);
1245  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1246  "specific element types");
1247 }
1248 
1249 /**
1250  * @brief Physical derivative along a direction vector.
1251  * @see StdRegions#StdExpansion#PhysDirectionalDeriv
1252  */
1254  const Array<OneD, const NekDouble> &inarray,
1255  const Array<OneD, const NekDouble> &direction,
1256  Array<OneD, NekDouble> &outarray)
1257 {
1258  boost::ignore_unused(inarray, direction, outarray);
1259  NEKERROR(ErrorUtil::efatal, "This function is only valid for "
1260  "specific element types");
1261 }
1262 
1264  Array<OneD, NekDouble> &out_d1,
1265  Array<OneD, NekDouble> &out_d2,
1266  Array<OneD, NekDouble> &out_d3)
1267 {
1268  boost::ignore_unused(inarray, out_d1, out_d2, out_d3);
1269  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1270 }
1271 
1272 void StdExpansion::v_StdPhysDeriv(const int dir,
1273  const Array<OneD, const NekDouble> &inarray,
1274  Array<OneD, NekDouble> &outarray)
1275 {
1276  boost::ignore_unused(dir, inarray, outarray);
1277  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1278 }
1279 
1281  const Array<OneD, const NekDouble> &coords,
1282  const Array<OneD, const NekDouble> &physvals)
1283 {
1284  boost::ignore_unused(coords, physvals);
1285  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1286  return 0;
1287 }
1288 
1291  const Array<OneD, const NekDouble> &physvals)
1292 {
1293  boost::ignore_unused(I, physvals);
1294  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1295  return 0;
1296 }
1297 
1299  const Array<OneD, const NekDouble> &coords, int mode)
1300 {
1301  boost::ignore_unused(coords, mode);
1302  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1303  return 0;
1304 }
1305 
1306 void StdExpansion::v_FillMode(const int mode, Array<OneD, NekDouble> &outarray)
1307 {
1308  boost::ignore_unused(mode, outarray);
1309  NEKERROR(ErrorUtil::efatal, "This function has not "
1310  "been defined for this shape");
1311 }
1312 
1314 {
1315  boost::ignore_unused(mkey);
1316  NEKERROR(ErrorUtil::efatal, "This function has not "
1317  "been defined for this element");
1318  DNekMatSharedPtr returnval;
1319  return returnval;
1320 }
1321 
1323 {
1324  boost::ignore_unused(mkey);
1325  NEKERROR(ErrorUtil::efatal, "This function has not "
1326  "been defined for this element");
1327  DNekMatSharedPtr returnval;
1328  return returnval;
1329 }
1330 
1332  Array<OneD, NekDouble> &coords_1,
1333  Array<OneD, NekDouble> &coords_2)
1334 {
1335  boost::ignore_unused(coords_0, coords_1, coords_2);
1336  NEKERROR(ErrorUtil::efatal, "Write coordinate definition method");
1337 }
1338 
1340  Array<OneD, NekDouble> &coord)
1341 {
1342  boost::ignore_unused(Lcoord, coord);
1343  NEKERROR(ErrorUtil::efatal, "Write coordinate definition method");
1344 }
1345 
1347 {
1348  NEKERROR(ErrorUtil::efatal, "Write method");
1349  return -1;
1350 }
1351 
1353 {
1354  boost::ignore_unused(outarray);
1355  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1356 }
1357 
1359 {
1360  boost::ignore_unused(outarray);
1361  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1362 }
1363 
1364 int StdExpansion::v_GetVertexMap(const int localVertexId, bool useCoeffPacking)
1365 {
1366  boost::ignore_unused(localVertexId, useCoeffPacking);
1367  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1368  return 0;
1369 }
1370 
1372  Array<OneD, unsigned int> &maparray,
1373  Array<OneD, int> &signarray,
1374  Orientation traceOrient, int P, int Q)
1375 {
1376  boost::ignore_unused(tid, maparray, signarray, traceOrient, P, Q);
1377  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1378 }
1379 
1380 void StdExpansion::v_GetTraceCoeffMap(const unsigned int traceid,
1381  Array<OneD, unsigned int> &maparray)
1382 {
1383  boost::ignore_unused(traceid, maparray);
1384  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1385 }
1386 
1387 void StdExpansion::v_GetElmtTraceToTraceMap(const unsigned int tid,
1388  Array<OneD, unsigned int> &maparray,
1389  Array<OneD, int> &signarray,
1390  Orientation traceOrient, int P,
1391  int Q)
1392 {
1393  boost::ignore_unused(tid, maparray, signarray, traceOrient, P, Q);
1394  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1395 }
1396 
1398  const int tid, Array<OneD, unsigned int> &maparray,
1399  Array<OneD, int> &signarray, const Orientation traceOrient)
1400 {
1401  boost::ignore_unused(tid, maparray, signarray, traceOrient);
1402  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1403 }
1404 
1405 void StdExpansion::v_GetTraceNumModes(const int tid, int &numModes0,
1406  int &numModes1, Orientation traceOrient)
1407 {
1408  boost::ignore_unused(tid, traceOrient, numModes0, numModes1);
1409  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1410 }
1411 
1413  const int vertex, const Array<OneD, const NekDouble> &inarray,
1414  NekDouble &outarray)
1415 {
1416  boost::ignore_unused(vertex, inarray, outarray);
1417  NEKERROR(ErrorUtil::efatal, "Method does not exist for "
1418  "this shape or library");
1419 }
1420 
1422  const Array<OneD, const NekDouble> &inarray,
1423  Array<OneD, NekDouble> &outarray)
1424 {
1425  v_MultiplyByStdQuadratureMetric(inarray, outarray);
1426 }
1427 
1429  const Array<OneD, const NekDouble> &inarray,
1430  Array<OneD, NekDouble> &outarray)
1431 {
1432  boost::ignore_unused(inarray, outarray);
1434  "Method does not exist for this shape or library");
1435 }
1436 
1438  const Array<OneD, const NekDouble> &inarray,
1439  Array<OneD, NekDouble> &outarray)
1440 {
1441  boost::ignore_unused(inarray, outarray);
1442  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1443 }
1444 
1446  const Array<OneD, const NekDouble> &inarray,
1447  Array<OneD, NekDouble> &outarray, bool multiplybyweights)
1448 {
1449  boost::ignore_unused(inarray, outarray, multiplybyweights);
1450  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1451 }
1452 
1453 /**
1454  *
1455  */
1457  const Array<OneD, const NekDouble> &direction,
1458  const Array<OneD, const NekDouble> &inarray,
1459  Array<OneD, NekDouble> &outarray)
1460 {
1461  boost::ignore_unused(direction, inarray, outarray);
1462  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1463 }
1464 
1466  const int dir, const Array<OneD, const NekDouble> &inarray,
1467  Array<OneD, NekDouble> &outarray)
1468 {
1469  boost::ignore_unused(dir, inarray, outarray);
1470  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
1471 }
1472 
1474  Array<OneD, NekDouble> &outarray,
1475  const StdMatrixKey &mkey)
1476 {
1477  // If this function is not reimplemented on shape level, the function
1478  // below will be called
1479  MassMatrixOp_MatFree(inarray, outarray, mkey);
1480 }
1481 
1483  const Array<OneD, const NekDouble> &inarray,
1484  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1485 {
1486  // If this function is not reimplemented on shape level, the function
1487  // below will be called
1488  LaplacianMatrixOp_MatFree(inarray, outarray, mkey);
1489 }
1490 
1492  const StdMatrixKey &mkey)
1493 {
1494  boost::ignore_unused(array, mkey);
1495  ASSERTL0(false, "This function is not defined in StdExpansion.");
1496 }
1497 
1499  const NekDouble alpha,
1500  const NekDouble exponent,
1501  const NekDouble cutoff)
1502 {
1503  boost::ignore_unused(array, alpha, exponent, cutoff);
1504  ASSERTL0(false, "This function is not defined in StdExpansion.");
1505 }
1506 
1508  int numMin, const Array<OneD, const NekDouble> &inarray,
1509  Array<OneD, NekDouble> &outarray)
1510 {
1511  boost::ignore_unused(numMin, inarray, outarray);
1512  ASSERTL0(false, "This function is not defined in StdExpansion.");
1513 }
1514 
1516  const int k1, const int k2, const Array<OneD, const NekDouble> &inarray,
1517  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1518 {
1519  // If this function is not reimplemented on shape level, the function
1520  // below will be called
1521  LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
1522 }
1523 
1525  const int i, const Array<OneD, const NekDouble> &inarray,
1526  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1527 {
1528  // If this function is not reimplemented on shape level, the function
1529  // below will be called
1530  WeakDerivMatrixOp_MatFree(i, inarray, outarray, mkey);
1531 }
1532 
1534  const Array<OneD, const NekDouble> &inarray,
1535  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1536 {
1537  // If this function is not reimplemented on shape level, the function
1538  // below will be called
1539  WeakDirectionalDerivMatrixOp_MatFree(inarray, outarray, mkey);
1540 }
1541 
1543  const Array<OneD, const NekDouble> &inarray,
1544  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1545 {
1546  // If this function is not reimplemented on shape level, the function
1547  // below will be called
1548  MassLevelCurvatureMatrixOp_MatFree(inarray, outarray, mkey);
1549 }
1550 
1552  const Array<OneD, const NekDouble> &inarray,
1553  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey,
1554  bool addDiffusionTerm)
1555 {
1556  // If this function is not reimplemented on shape level, the function
1557  // below will be called
1558  LinearAdvectionDiffusionReactionMatrixOp_MatFree(inarray, outarray, mkey,
1559  addDiffusionTerm);
1560 }
1561 
1563  const Array<OneD, const NekDouble> &inarray,
1564  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1565 {
1566  // If this function is not reimplemented on shape level, the function
1567  // below will be called
1568  HelmholtzMatrixOp_MatFree(inarray, outarray, mkey);
1569 }
1570 
1572  const Array<OneD, const NekDouble> &inarray,
1573  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1574 {
1575  // If this function is not reimplemented on shape level, the function
1576  // below will be called
1577  LaplacianMatrixOp_MatFree_GenericImpl(inarray, outarray, mkey);
1578 }
1579 
1581  const Array<OneD, const NekDouble> &inarray,
1583 {
1584  boost::ignore_unused(inarray, outarray, wsp);
1585  ASSERTL0(false, "Not implemented.");
1586 }
1587 
1589  const Array<OneD, const NekDouble> &inarray,
1590  Array<OneD, NekDouble> &outarray, const StdMatrixKey &mkey)
1591 {
1592  // If this function is not reimplemented on shape level, the function
1593  // below will be called
1594  HelmholtzMatrixOp_MatFree_GenericImpl(inarray, outarray, mkey);
1595 }
1596 
1598  const DNekScalMatSharedPtr &m_transformationmatrix)
1599 {
1600  boost::ignore_unused(m_transformationmatrix);
1601  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
1602  return NullDNekMatSharedPtr;
1603 }
1604 
1606  const Array<OneD, const NekDouble> &inarray,
1607  Array<OneD, NekDouble> &outarray, int npset)
1608 {
1610  DNekMatSharedPtr intmat;
1611 
1612  int nqtot = GetTotPoints();
1613  int np = 0;
1614  if (npset == -1) // use values from basis num points()
1615  {
1616  int nqbase;
1617  for (int i = 0; i < m_base.size(); ++i)
1618  {
1619  nqbase = m_base[i]->GetNumPoints();
1620  np = std::max(np, nqbase);
1621  }
1622 
1623  StdMatrixKey Ikey(ePhysInterpToEquiSpaced, shape, *this);
1624  intmat = GetStdMatrix(Ikey);
1625  }
1626  else
1627  {
1628  np = npset;
1629 
1630  ConstFactorMap cmap;
1631  cmap[eFactorConst] = np;
1632  StdMatrixKey Ikey(ePhysInterpToEquiSpaced, shape, *this, cmap);
1633  intmat = GetStdMatrix(Ikey);
1634  }
1635 
1636  NekVector<NekDouble> in(nqtot, inarray, eWrapper);
1638  LibUtilities::GetNumberOfCoefficients(shape, np, np, np), outarray,
1639  eWrapper);
1640  out = (*intmat) * in;
1641 }
1642 
1644  bool standard)
1645 {
1646  boost::ignore_unused(conn, standard);
1648  "GetSimplexEquiSpacedConnectivity not"
1649  " implemented for " +
1650  static_cast<std::string>(
1652 }
1653 
1655  const Array<OneD, const NekDouble> &inarray,
1656  Array<OneD, NekDouble> &outarray)
1657 {
1659 
1660  // inarray has to be consistent with NumModes definition
1661  // There is also a check in GetStdMatrix to see if all
1662  // modes are of the same size
1663  ConstFactorMap cmap;
1664 
1665  cmap[eFactorConst] = m_base[0]->GetNumModes();
1666  StdMatrixKey Ikey(eEquiSpacedToCoeffs, shape, *this, cmap);
1667  DNekMatSharedPtr intmat = GetStdMatrix(Ikey);
1668 
1669  NekVector<NekDouble> in(m_ncoeffs, inarray, eWrapper);
1670  NekVector<NekDouble> out(m_ncoeffs, outarray, eWrapper);
1671  out = (*intmat) * in;
1672 }
1673 
1674 } // namespace StdRegions
1675 } // 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
Describes the specification for a Basis.
Definition: Basis.h:50
Defines a specification for a set of points.
Definition: Points.h:59
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
The base class for all shapes.
Definition: StdExpansion.h:71
void GetBoundaryMap(Array< OneD, unsigned int > &outarray)
Definition: StdExpansion.h:677
virtual void v_LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual int v_GetVertexMap(int localVertexId, bool useCoeffPacking=false)
virtual void v_LocCollapsedToLocCoord(const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi)
virtual ~StdExpansion()
Destructor.
virtual void v_PhysDeriv_n(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dn)
virtual void v_PhysDirectionalDeriv(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &direction, Array< OneD, NekDouble > &outarray)
Physical derivative along a direction vector.
StdExpansion()
Default Constructor.
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
Definition: StdExpansion.h:130
virtual void v_GetTraceToElementMap(const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:140
virtual void v_GetCoord(const Array< OneD, const NekDouble > &Lcoord, Array< OneD, NekDouble > &coord)
void GeneralMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void FillMode(const int mode, Array< OneD, NekDouble > &outarray)
This function fills the array outarray with the mode-th mode of the expansion.
Definition: StdExpansion.h:499
virtual NekDouble v_StdPhysEvaluate(const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
void PhysInterpToSimplexEquiSpaced(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, int npset=-1)
This function performs an interpolation from the physical space points provided at input into an arra...
void WeakDirectionalDerivMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:809
void WeakDerivMatrixOp_MatFree(const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
DNekBlkMatSharedPtr CreateStdStaticCondMatrix(const StdMatrixKey &mkey)
Create the static condensation of a matrix when using a boundary interior decomposition.
void LaplacianMatrixOp_MatFree_GenericImpl(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_IProductWRTDirectionalDerivBase_SumFac(const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual int v_GetTraceNumPoints(const int i) const
void HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void BwdTrans_MatOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void MassLevelCurvatureMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:816
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_0, Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2)
void LinearAdvectionDiffusionReactionMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
virtual std::shared_ptr< StdExpansion > v_GetStdExp(void) const
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:611
virtual void v_GetElmtTraceToTraceMap(const unsigned int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
virtual bool v_IsBoundaryInteriorExpansion()
virtual void v_NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
virtual int v_GetTraceNcoeffs(const int i) const
virtual void v_GetSimplexEquiSpacedConnectivity(Array< OneD, int > &conn, bool standard=true)
virtual void v_SVVLaplacianFilter(Array< OneD, NekDouble > &array, const StdMatrixKey &mkey)
virtual std::shared_ptr< StdExpansion > v_GetLinStdExp(void) const
void WeakDerivMatrixOp(const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:801
virtual void v_GetInteriorMap(Array< OneD, unsigned int > &outarray)
virtual void v_SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
NekDouble Linf(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
Function to evaluate the discrete error where is given by the array sol.
void EquiSpacedToCoeffs(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs a projection/interpolation from the equispaced points sometimes used in post-p...
void MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:760
virtual int v_GetTotalTraceIntNcoeffs() const
virtual void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_WeakDirectionalDerivMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void BwdTrans_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
virtual DNekMatSharedPtr v_BuildInverseTransformationMatrix(const DNekScalMatSharedPtr &m_transformationmatrix)
NekDouble L2(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
Function to evaluate the discrete error, where is given by the array sol.
virtual void v_GetVertexPhysVals(const int vertex, const Array< OneD, const NekDouble > &inarray, NekDouble &outarray)
virtual int v_GetShapeDimension() const
void IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual int v_GetTraceIntNcoeffs(const int i) const
DNekMatSharedPtr CreateGeneralMatrix(const StdMatrixKey &mkey)
this function generates the mass matrix
virtual int v_NumBndryCoeffs() const
virtual void v_MassLevelCurvatureMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_GetTraceCoeffMap(const unsigned int traceid, Array< OneD, unsigned int > &maparray)
virtual void v_MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual LibUtilities::ShapeType v_DetShapeType() const
void LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:767
virtual void v_FillMode(const int mode, Array< OneD, NekDouble > &outarray)
void LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_MultiplyByStdQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_ExponentialFilter(Array< OneD, NekDouble > &array, const NekDouble alpha, const NekDouble exponent, const NekDouble cutoff)
void MassLevelCurvatureMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:375
virtual void v_LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual DNekMatSharedPtr v_GenMatrix(const StdMatrixKey &mkey)
void GetInteriorMap(Array< OneD, unsigned int > &outarray)
Definition: StdExpansion.h:682
virtual DNekMatSharedPtr v_CreateStdMatrix(const StdMatrixKey &mkey)
virtual void v_BwdTrans_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual int v_NumDGBndryCoeffs() const
virtual void v_GetTraceNumModes(const int fid, int &numModes0, int &numModes1, Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2)
virtual void v_PhysDeriv_s(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_ds)
virtual void v_FwdTransBndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_GetBoundaryMap(Array< OneD, unsigned int > &outarray)
void GenStdMatBwdDeriv(const int dir, DNekMatSharedPtr &mat)
virtual void v_MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual const LibUtilities::PointsKey v_GetNodalPointsKey() const
virtual void v_WeakDerivMatrixOp(const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void LinearAdvectionDiffusionReactionMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
Definition: StdExpansion.h:823
void IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
virtual void v_LocCoordToLocCollapsed(const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
virtual NekDouble v_PhysEvaluateBasis(const Array< OneD, const NekDouble > &coords, int mode)
virtual void v_GetTraceInteriorToElementMap(const int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eForwards)
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)=0
Calculates the inner product of a given function f with the different modes of the expansion.
virtual void v_LinearAdvectionDiffusionReactionMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
void HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:838
virtual void v_ReduceOrderCoeffs(int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_DropLocStaticCondMatrix(const LocalRegions::MatrixKey &mkey)
void WeakDirectionalDerivMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_LaplacianMatrixOp_MatFree_Kernel(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
NekDouble H1(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
Function to evaluate the discrete error, where is given by the array sol.
virtual int v_CalcNumberOfCoefficients(const std::vector< unsigned int > &nummodes, int &modes_offset)
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Integrates the specified function over the domain.
void SVVLaplacianFilter(Array< OneD, NekDouble > &array, const StdMatrixKey &mkey)
Definition: StdExpansion.h:781
void HelmholtzMatrixOp_MatFree_GenericImpl(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)=0
virtual LibUtilities::PointsKey v_GetTracePointsKey(const int i, const int j) const
Array< OneD, LibUtilities::BasisSharedPtr > m_base
virtual void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
virtual void v_IProductWRTDirectionalDerivBase(const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2, Array< OneD, NekDouble > &out_d3)
Calculate the derivative of the physical points.
virtual const LibUtilities::BasisKey v_GetTraceBasisKey(const int i, const int k) const
virtual void v_StdPhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2, Array< OneD, NekDouble > &out_d3)
NekDouble StdPhysEvaluate(const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
void MassMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void GeneralMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
virtual DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix(const LocalRegions::MatrixKey &mkey)
const ConstFactorMap & GetConstFactors() const
Definition: StdMatrixKey.h:140
LibUtilities::ShapeType GetShapeType() const
Definition: StdMatrixKey.h:90
const VarCoeffMap & GetVarCoeffs() const
Definition: StdMatrixKey.h:167
MatrixType GetMatrixType() const
Definition: StdMatrixKey.h:85
const Array< OneD, const NekDouble > & GetVarCoeff(const StdRegions::VarCoeffType &coeff) const
Definition: StdMatrixKey.h:150
bool HasVarCoeff(const StdRegions::VarCoeffType &coeff) const
Definition: StdMatrixKey.h:172
LibUtilities::PointsType GetNodalPointsType() const
Definition: StdMatrixKey.h:95
NekDouble GetConstFactor(const ConstFactorType &factor) const
Definition: StdMatrixKey.h:126
bool ConstFactorExists(const ConstFactorType &factor) const
Definition: StdMatrixKey.h:135
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].
Definition: Blas.hpp:246
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.
Definition: Blas.hpp:154
const char *const ShapeTypeMap[SIZE_ShapeType]
Definition: ShapeType.hpp:77
BasisManagerT & BasisManager(void)
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
int GetNumberOfCoefficients(ShapeType shape, std::vector< unsigned int > &modes, int offset=0)
Definition: ShapeType.hpp:308
static const PointsKey NullPointsKey(0, eNoPointsType)
static const NekDouble kNekSqrtTol
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
const char *const MatrixTypeMap[]
Definition: StdRegions.hpp:140
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:282
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:283
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:241
std::shared_ptr< StdMatrixKey > StdMatrixKeySharedPtr
Definition: StdMatrixKey.h:197
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
static DNekScalBlkMatSharedPtr NullDNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:85
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:79
static DNekMatSharedPtr NullDNekMatSharedPtr
Definition: NekTypeDefs.hpp:83
static Array< OneD, NekDouble > NullNekDouble1DArray
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 Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
Definition: Vmath.cpp:553
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 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.
Definition: Vmath.cpp:359
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:248
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492
T Vamax(int n, const T *x, const int incx)
Return the maximum absolute element in x called vamax to avoid conflict with max.
Definition: Vmath.cpp:999
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255
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.
Definition: Vmath.cpp:419
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:291