Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
StdExpansion3D.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: StdExpansion3D.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: Daughter of StdExpansion. This class contains routine
32 // which are common to 3D expansion. Typically this inolves physiocal
33 // space operations.
34 //
35 ///////////////////////////////////////////////////////////////////////////////
36 
37 #include <boost/core/ignore_unused.hpp>
38 
40 
41 #ifdef max
42 #undef max
43 #endif
44 
45 namespace Nektar
46 {
47 namespace StdRegions
48 {
50 {
51 }
52 
54  const LibUtilities::BasisKey &Bb,
55  const LibUtilities::BasisKey &Bc)
56  : StdExpansion(numcoeffs, 3, Ba, Bb, Bc)
57 {
58 }
59 
61 {
62 }
63 
65 {
66 }
67 
71 {
72  const int nquad0 = m_base[0]->GetNumPoints();
73  const int nquad1 = m_base[1]->GetNumPoints();
74  const int nquad2 = m_base[2]->GetNumPoints();
75 
76  Array<OneD, NekDouble> wsp(nquad0 * nquad1 * nquad2);
77 
78  // copy inarray to wsp in case inarray is used as outarray
79  Vmath::Vcopy(nquad0 * nquad1 * nquad2, &inarray[0], 1, &wsp[0], 1);
80 
81  if (out_dx.size() > 0)
82  {
83  NekDouble *D0 = &((m_base[0]->GetD())->GetPtr())[0];
84 
85  Blas::Dgemm('N', 'N', nquad0, nquad1 * nquad2, nquad0, 1.0, D0, nquad0,
86  &wsp[0], nquad0, 0.0, &out_dx[0], nquad0);
87  }
88 
89  if (out_dy.size() > 0)
90  {
91  NekDouble *D1 = &((m_base[1]->GetD())->GetPtr())[0];
92  for (int j = 0; j < nquad2; ++j)
93  {
94  Blas::Dgemm('N', 'T', nquad0, nquad1, nquad1, 1.0,
95  &wsp[j * nquad0 * nquad1], nquad0, D1, nquad1, 0.0,
96  &out_dy[j * nquad0 * nquad1], nquad0);
97  }
98  }
99 
100  if (out_dz.size() > 0)
101  {
102  NekDouble *D2 = &((m_base[2]->GetD())->GetPtr())[0];
103 
104  Blas::Dgemm('N', 'T', nquad0 * nquad1, nquad2, nquad2, 1.0, &wsp[0],
105  nquad0 * nquad1, D2, nquad2, 0.0, &out_dz[0],
106  nquad0 * nquad1);
107  }
108 }
109 
111  const Array<OneD, const NekDouble> &base0,
112  const Array<OneD, const NekDouble> &base1,
113  const Array<OneD, const NekDouble> &base2,
114  const Array<OneD, const NekDouble> &inarray,
116  bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
117 {
118  v_BwdTrans_SumFacKernel(base0, base1, base2, inarray, outarray, wsp,
119  doCheckCollDir0, doCheckCollDir1, doCheckCollDir2);
120 }
121 
123  const Array<OneD, const NekDouble> &base0,
124  const Array<OneD, const NekDouble> &base1,
125  const Array<OneD, const NekDouble> &base2,
126  const Array<OneD, const NekDouble> &inarray,
128  bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
129 {
130  v_IProductWRTBase_SumFacKernel(base0, base1, base2, inarray, outarray, wsp,
131  doCheckCollDir0, doCheckCollDir1,
132  doCheckCollDir2);
133 }
134 
136 {
137  ASSERTL1((dir == 0) || (dir == 1) || (dir == 2), "Invalid direction.");
138 
139  const int nq0 = m_base[0]->GetNumPoints();
140  const int nq1 = m_base[1]->GetNumPoints();
141  const int nq2 = m_base[2]->GetNumPoints();
142  const int nq = nq0 * nq1 * nq2;
143  const int nm0 = m_base[0]->GetNumModes();
144  const int nm1 = m_base[1]->GetNumModes();
145 
146  Array<OneD, NekDouble> alloc(4 * nq + m_ncoeffs + nm0 * nq2 * (nq1 + nm1),
147  0.0);
148  Array<OneD, NekDouble> tmp1(alloc); // Quad metric
149  Array<OneD, NekDouble> tmp2(alloc + nq); // Dir1 metric
150  Array<OneD, NekDouble> tmp3(alloc + 2 * nq); // Dir2 metric
151  Array<OneD, NekDouble> tmp4(alloc + 3 * nq); // Dir3 metric
152  Array<OneD, NekDouble> tmp5(alloc + 4 * nq); // iprod tmp
153  Array<OneD, NekDouble> wsp(tmp5 + m_ncoeffs); // Wsp
154  switch (dir)
155  {
156  case 0:
157  for (int i = 0; i < nq; i++)
158  {
159  tmp2[i] = 1.0;
161  m_base[0]->GetDbdata(), m_base[1]->GetBdata(),
162  m_base[2]->GetBdata(), tmp2, tmp5, wsp, false, true, true);
163 
164  tmp2[i] = 0.0;
165 
166  for (int j = 0; j < m_ncoeffs; j++)
167  {
168  (*mat)(j, i) = tmp5[j];
169  }
170  }
171  break;
172  case 1:
173  for (int i = 0; i < nq; i++)
174  {
175  tmp2[i] = 1.0;
177  m_base[0]->GetBdata(), m_base[1]->GetDbdata(),
178  m_base[2]->GetBdata(), tmp2, tmp5, wsp, true, false, true);
179 
180  tmp2[i] = 0.0;
181 
182  for (int j = 0; j < m_ncoeffs; j++)
183  {
184  (*mat)(j, i) = tmp5[j];
185  }
186  }
187  break;
188  case 2:
189  for (int i = 0; i < nq; i++)
190  {
191  tmp2[i] = 1.0;
193  m_base[0]->GetBdata(), m_base[1]->GetBdata(),
194  m_base[2]->GetDbdata(), tmp2, tmp5, wsp, true, true, false);
195  tmp2[i] = 0.0;
196 
197  for (int j = 0; j < m_ncoeffs; j++)
198  {
199  (*mat)(j, i) = tmp5[j];
200  }
201  }
202  break;
203  default:
204  NEKERROR(ErrorUtil::efatal, "Not a 2D expansion.");
205  break;
206  }
207 }
208 
210  const Array<OneD, const NekDouble> &coords,
211  const Array<OneD, const NekDouble> &physvals)
212 {
213  Array<OneD, NekDouble> eta(3);
214 
215  WARNINGL2(coords[0] >= -1 - NekConstants::kNekZeroTol, "coord[0] < -1");
216  WARNINGL2(coords[0] <= 1 + NekConstants::kNekZeroTol, "coord[0] > 1");
217  WARNINGL2(coords[1] >= -1 - NekConstants::kNekZeroTol, "coord[1] < -1");
218  WARNINGL2(coords[1] <= 1 + NekConstants::kNekZeroTol, "coord[1] > 1");
219  WARNINGL2(coords[2] >= -1 - NekConstants::kNekZeroTol, "coord[2] < -1");
220  WARNINGL2(coords[2] <= 1 + NekConstants::kNekZeroTol, "coord[2] > 1");
221 
222  // Obtain local collapsed corodinate from Cartesian coordinate.
223  LocCoordToLocCollapsed(coords, eta);
224 
225  const int nq0 = m_base[0]->GetNumPoints();
226  const int nq1 = m_base[1]->GetNumPoints();
227  const int nq2 = m_base[2]->GetNumPoints();
228 
229  Array<OneD, NekDouble> wsp1(nq1 * nq2), wsp2(nq2);
230 
231  // Construct the 2D square...
232  const NekDouble *ptr = &physvals[0];
233  for (int i = 0; i < nq1 * nq2; ++i, ptr += nq0)
234  {
235  wsp1[i] = StdExpansion::BaryEvaluate<0>(eta[0], ptr);
236  }
237 
238  for (int i = 0; i < nq2; ++i)
239  {
240  wsp2[i] = StdExpansion::BaryEvaluate<1>(eta[1], &wsp1[i * nq1]);
241  }
242 
243  return StdExpansion::BaryEvaluate<2>(eta[2], &wsp2[0]);
244 }
245 
248  const Array<OneD, const NekDouble> &physvals)
249 {
250  NekDouble value;
251 
252  int Qx = m_base[0]->GetNumPoints();
253  int Qy = m_base[1]->GetNumPoints();
254  int Qz = m_base[2]->GetNumPoints();
255 
256  Array<OneD, NekDouble> sumFactorization_qr =
257  Array<OneD, NekDouble>(Qy * Qz);
258  Array<OneD, NekDouble> sumFactorization_r = Array<OneD, NekDouble>(Qz);
259 
260  // Lagrangian interpolation matrix
261  NekDouble *interpolatingNodes = 0;
262 
263  // Interpolate first coordinate direction
264  interpolatingNodes = &I[0]->GetPtr()[0];
265 
266  Blas::Dgemv('T', Qx, Qy * Qz, 1.0, &physvals[0], Qx, &interpolatingNodes[0],
267  1, 0.0, &sumFactorization_qr[0], 1);
268 
269  // Interpolate in second coordinate direction
270  interpolatingNodes = &I[1]->GetPtr()[0];
271 
272  Blas::Dgemv('T', Qy, Qz, 1.0, &sumFactorization_qr[0], Qy,
273  &interpolatingNodes[0], 1, 0.0, &sumFactorization_r[0], 1);
274 
275  // Interpolate in third coordinate direction
276  interpolatingNodes = &I[2]->GetPtr()[0];
277  value = Blas::Ddot(Qz, interpolatingNodes, 1, &sumFactorization_r[0], 1);
278 
279  return value;
280 }
281 
283  const Array<OneD, NekDouble> &coord,
284  const Array<OneD, const NekDouble> &inarray,
285  std::array<NekDouble, 3> &firstOrderDerivs)
286 {
287  boost::ignore_unused(coord, inarray, firstOrderDerivs);
288  return 0;
289 }
290 
291 /**
292  * @param inarray Input coefficients.
293  * @param output Output coefficients.
294  * @param mkey Matrix key
295  */
297  const Array<OneD, const NekDouble> &inarray,
298  Array<OneD, NekDouble> &outarray, const StdRegions::StdMatrixKey &mkey)
299 {
300  if (mkey.GetNVarCoeff() == 0 &&
303  {
304  // This implementation is only valid when there are no
305  // coefficients associated to the Laplacian operator
306  int nqtot = GetTotPoints();
307 
308  const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
309  const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
310  const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
311 
312  // Allocate temporary storage
313  Array<OneD, NekDouble> wsp0(7 * nqtot);
314  Array<OneD, NekDouble> wsp1(wsp0 + nqtot);
315 
316  if (!(m_base[0]->Collocation() && m_base[1]->Collocation() &&
317  m_base[2]->Collocation()))
318  {
319  // LAPLACIAN MATRIX OPERATION
320  // wsp0 = u = B * u_hat
321  // wsp1 = du_dxi1 = D_xi1 * wsp0 = D_xi1 * u
322  // wsp2 = du_dxi2 = D_xi2 * wsp0 = D_xi2 * u
323  BwdTrans_SumFacKernel(base0, base1, base2, inarray, wsp0, wsp1,
324  true, true, true);
325  LaplacianMatrixOp_MatFree_Kernel(wsp0, outarray, wsp1);
326  }
327  else
328  {
329  LaplacianMatrixOp_MatFree_Kernel(inarray, outarray, wsp1);
330  }
331  }
332  else
333  {
335  mkey);
336  }
337 }
338 
340  const Array<OneD, const NekDouble> &inarray,
341  Array<OneD, NekDouble> &outarray, const StdRegions::StdMatrixKey &mkey)
342 {
343  if (mkey.GetNVarCoeff() == 0 &&
345  {
346  using std::max;
347 
348  int nquad0 = m_base[0]->GetNumPoints();
349  int nquad1 = m_base[1]->GetNumPoints();
350  int nquad2 = m_base[2]->GetNumPoints();
351  int nmodes0 = m_base[0]->GetNumModes();
352  int nmodes1 = m_base[1]->GetNumModes();
353  int nmodes2 = m_base[2]->GetNumModes();
354  int wspsize = max(nquad0 * nmodes2 * (nmodes1 + nquad1),
355  nquad0 * nquad1 * (nquad2 + nmodes0) +
356  nmodes0 * nmodes1 * nquad2);
357 
359 
360  const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
361  const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
362  const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
363  Array<OneD, NekDouble> wsp0(8 * wspsize);
364  Array<OneD, NekDouble> wsp1(wsp0 + 1 * wspsize);
365  Array<OneD, NekDouble> wsp2(wsp0 + 2 * wspsize);
366 
367  if (!(m_base[0]->Collocation() && m_base[1]->Collocation() &&
368  m_base[2]->Collocation()))
369  {
370  // MASS MATRIX OPERATION
371  // The following is being calculated:
372  // wsp0 = B * u_hat = u
373  // wsp1 = W * wsp0
374  // outarray = B^T * wsp1 = B^T * W * B * u_hat = M * u_hat
375  BwdTrans_SumFacKernel(base0, base1, base2, inarray, wsp0, wsp2,
376  true, true, true);
377  MultiplyByQuadratureMetric(wsp0, wsp1);
378  IProductWRTBase_SumFacKernel(base0, base1, base2, wsp1, outarray,
379  wsp2, true, true, true);
380  LaplacianMatrixOp_MatFree_Kernel(wsp0, wsp1, wsp2);
381  }
382  else
383  {
384  // specialised implementation for the classical spectral
385  // element method
386  MultiplyByQuadratureMetric(inarray, outarray);
387  LaplacianMatrixOp_MatFree_Kernel(inarray, wsp1, wsp2);
388  }
389 
390  // outarray = lambda * outarray + wsp1
391  // = (lambda * M + L ) * u_hat
392  Vmath::Svtvp(m_ncoeffs, lambda, &outarray[0], 1, &wsp1[0], 1,
393  &outarray[0], 1);
394  }
395  else
396  {
398  mkey);
399  }
400 }
401 
403  const Array<OneD, const NekDouble> &inarray)
404 {
405  const int nqtot = GetTotPoints();
407  v_MultiplyByStdQuadratureMetric(inarray, tmp);
408  return Vmath::Vsum(nqtot, tmp, 1);
409 }
410 
412 {
413  NEKERROR(ErrorUtil::efatal, "This function is not valid or not defined");
414  return 0;
415 }
416 
417 int StdExpansion3D::v_GetEdgeNcoeffs(const int i) const
418 {
419  boost::ignore_unused(i);
420  NEKERROR(ErrorUtil::efatal, "This function is not valid or not defined");
421  return 0;
422 }
423 
425  const int tid, Array<OneD, unsigned int> &maparray,
426  Array<OneD, int> &signarray, Orientation traceOrient)
427 {
428  boost::ignore_unused(tid, maparray, signarray, traceOrient);
429  NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
430 }
431 
433  Array<OneD, unsigned int> &maparray,
434  Array<OneD, int> &signarray,
435  Orientation traceOrient, int P,
436  int Q)
437 {
438  Array<OneD, unsigned int> map1, map2;
439  GetTraceCoeffMap(tid, map1);
440  GetElmtTraceToTraceMap(tid, map2, signarray, traceOrient, P, Q);
441 
442  if (maparray.size() != map2.size())
443  {
444  maparray = Array<OneD, unsigned int>(map2.size());
445  }
446 
447  for (int i = 0; i < map2.size(); ++i)
448  {
449  maparray[i] = map1[map2[i]];
450  }
451 }
452 
454  const int facedir, const LibUtilities::BasisType faceDirBasisType,
455  const int numpoints, const int nummodes)
456 {
457  boost::ignore_unused(facedir);
458 
459  switch (faceDirBasisType)
460  {
462  {
463  const LibUtilities::PointsKey pkey(
466  pkey);
467  }
470  {
471  const LibUtilities::PointsKey pkey(
472  numpoints + 1, LibUtilities::eGaussLobattoLegendre);
474  pkey);
475  }
477  {
478  const LibUtilities::PointsKey pkey(
481  pkey);
482  }
484  {
485  const LibUtilities::PointsKey pkey(
488  pkey);
489  }
492  {
493  const LibUtilities::PointsKey pkey(
494  numpoints + 1, LibUtilities::eGaussLobattoLegendre);
496  pkey);
497  }
498  default:
499  {
500  NEKERROR(ErrorUtil::efatal, "expansion type unknown");
501  break;
502  }
503  }
504 
505  // Keep things happy by returning a value.
507 }
508 
510  const int facedir, const LibUtilities::BasisType faceDirBasisType,
511  const int numpoints, const int nummodes)
512 {
513  switch (faceDirBasisType)
514  {
516  {
517  const LibUtilities::PointsKey pkey(
520  pkey);
521  }
525  {
526  switch (facedir)
527  {
528  case 0:
529  {
530  const LibUtilities::PointsKey pkey(
531  numpoints + 1, LibUtilities::eGaussLobattoLegendre);
533  nummodes, pkey);
534  }
535  case 1:
536  {
537  // const LibUtilities::PointsKey pkey(
538  // numpoints+1,
539  // LibUtilities::eGaussLobattoLegendre);
540  const LibUtilities::PointsKey pkey(
541  numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
543  nummodes, pkey);
544  }
545  default:
546  {
547 
548  NEKERROR(ErrorUtil::efatal, "invalid value to flag");
549  break;
550  }
551  }
552  break;
553  }
554 
556  {
557  switch (facedir)
558  {
559  case 0:
560  {
561  const LibUtilities::PointsKey pkey(
564  nummodes, pkey);
565  }
566  case 1:
567  {
568  const LibUtilities::PointsKey pkey(
569  numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
571  nummodes, pkey);
572  }
573  default:
574  {
575  NEKERROR(ErrorUtil::efatal, "invalid value to flag");
576  break;
577  }
578  }
579  break;
580  }
581 
586  {
587  switch (facedir)
588  {
589  case 0:
590  {
591  const LibUtilities::PointsKey pkey(
594  nummodes, pkey);
595  }
596  case 1:
597  {
598  const LibUtilities::PointsKey pkey(
599  numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
601  nummodes, pkey);
602  }
603  default:
604  {
605  NEKERROR(ErrorUtil::efatal, "invalid value to flag");
606  break;
607  }
608  }
609  break;
610  }
611  default:
612  {
613  NEKERROR(ErrorUtil::efatal, "expansion type unknown");
614  break;
615  }
616  }
617 
618  // Keep things happy by returning a value.
620 }
621 } // namespace StdRegions
622 } // namespace Nektar
#define WARNINGL2(condition, msg)
Definition: ErrorUtil.hpp:273
#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
Describes the specification for a Basis.
Definition: Basis.h:50
Defines a specification for a set of points.
Definition: Points.h:59
virtual int v_GetNedges(void) const
void BwdTrans_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
virtual void v_GetEdgeInteriorToElementMap(const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards)
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
virtual void v_BwdTrans_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)=0
virtual void v_HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
virtual void v_LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
virtual int v_GetEdgeNcoeffs(const int i) const
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
This function evaluates the expansion at a single (arbitrary) point of the domain.
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray) override
Integrates the specified function over the domain.
virtual void v_GenStdMatBwdDeriv(const int dir, DNekMatSharedPtr &mat) override
virtual void v_GetTraceToElementMap(const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient, int P, int Q) override
virtual void v_IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)=0
void PhysTensorDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray_d1, Array< OneD, NekDouble > &outarray_d2, Array< OneD, NekDouble > &outarray_d3)
Calculate the 3D derivative in the local tensor/collapsed coordinate at the physical points.
The base class for all shapes.
Definition: StdExpansion.h:71
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:140
void GetElmtTraceToTraceMap(const unsigned int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
Definition: StdExpansion.h:705
void LaplacianMatrixOp_MatFree_GenericImpl(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
void LocCoordToLocCollapsed(const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
Convert local cartesian coordinate xi into local collapsed coordinates eta.
void GetTraceCoeffMap(const unsigned int traceid, Array< OneD, unsigned int > &maparray)
Definition: StdExpansion.h:699
virtual void v_MultiplyByStdQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:729
void HelmholtzMatrixOp_MatFree_GenericImpl(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Array< OneD, LibUtilities::BasisSharedPtr > m_base
void LaplacianMatrixOp_MatFree_Kernel(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
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 double Ddot(const int &n, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: output = .
Definition: Blas.hpp:182
static void Dgemm(const char &transa, const char &transb, const int &m, const int &n, const int &k, const double &alpha, const double *a, const int &lda, const double *b, const int &ldb, const double &beta, double *c, const int &ldc)
BLAS level 3: Matrix-matrix multiply C = A x B where op(A)[m x k], op(B)[k x n], C[m x n] DGEMM perfo...
Definition: Blas.hpp:368
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:53
@ eModified_B
Principle Modified Functions .
Definition: BasisType.h:51
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:44
@ eModified_C
Principle Modified Functions .
Definition: BasisType.h:52
@ eGLL_Lagrange
Lagrange for SEM basis .
Definition: BasisType.h:58
@ eOrtho_C
Principle Orthogonal Functions .
Definition: BasisType.h:48
@ eModifiedPyr_C
Principle Modified Functions.
Definition: BasisType.h:55
@ eOrtho_B
Principle Orthogonal Functions .
Definition: BasisType.h:46
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:50
@ eOrthoPyr_C
Principle Orthogonal Functions .
Definition: BasisType.h:53
static const NekDouble kNekZeroTol
LibUtilities::BasisKey EvaluateTriFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes)
LibUtilities::BasisKey EvaluateQuadFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes)
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:2
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble
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
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.cpp:895
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255