Nektar++
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
39
40#ifdef max
41#undef max
42#endif
43
44namespace Nektar::StdRegions
45{
46
48 [[maybe_unused]] int numcoeffs,
49 [[maybe_unused]] const LibUtilities::BasisKey &Ba,
50 [[maybe_unused]] const LibUtilities::BasisKey &Bb,
51 [[maybe_unused]] const LibUtilities::BasisKey &Bc)
52{
53}
54
58{
59 const int nquad0 = m_base[0]->GetNumPoints();
60 const int nquad1 = m_base[1]->GetNumPoints();
61 const int nquad2 = m_base[2]->GetNumPoints();
62
63 Array<OneD, NekDouble> wsp(nquad0 * nquad1 * nquad2);
64
65 // copy inarray to wsp in case inarray is used as outarray
66 Vmath::Vcopy(nquad0 * nquad1 * nquad2, &inarray[0], 1, &wsp[0], 1);
67
68 if (out_dx.size() > 0)
69 {
70 NekDouble *D0 = &((m_base[0]->GetD())->GetPtr())[0];
71
72 Blas::Dgemm('N', 'N', nquad0, nquad1 * nquad2, nquad0, 1.0, D0, nquad0,
73 &wsp[0], nquad0, 0.0, &out_dx[0], nquad0);
74 }
75
76 if (out_dy.size() > 0)
77 {
78 NekDouble *D1 = &((m_base[1]->GetD())->GetPtr())[0];
79 for (int j = 0; j < nquad2; ++j)
80 {
81 Blas::Dgemm('N', 'T', nquad0, nquad1, nquad1, 1.0,
82 &wsp[j * nquad0 * nquad1], nquad0, D1, nquad1, 0.0,
83 &out_dy[j * nquad0 * nquad1], nquad0);
84 }
85 }
86
87 if (out_dz.size() > 0)
88 {
89 NekDouble *D2 = &((m_base[2]->GetD())->GetPtr())[0];
90
91 Blas::Dgemm('N', 'T', nquad0 * nquad1, nquad2, nquad2, 1.0, &wsp[0],
92 nquad0 * nquad1, D2, nquad2, 0.0, &out_dz[0],
93 nquad0 * nquad1);
94 }
95}
96
100 const Array<OneD, const NekDouble> &base2,
101 const Array<OneD, const NekDouble> &inarray,
103 bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
104{
105 v_BwdTrans_SumFacKernel(base0, base1, base2, inarray, outarray, wsp,
106 doCheckCollDir0, doCheckCollDir1, doCheckCollDir2);
107}
108
110 const Array<OneD, const NekDouble> &base0,
111 const Array<OneD, const NekDouble> &base1,
112 const Array<OneD, const NekDouble> &base2,
113 const Array<OneD, const NekDouble> &inarray,
115 bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
116{
117 v_IProductWRTBase_SumFacKernel(base0, base1, base2, inarray, outarray, wsp,
118 doCheckCollDir0, doCheckCollDir1,
119 doCheckCollDir2);
120}
121
123{
124 ASSERTL1((dir == 0) || (dir == 1) || (dir == 2), "Invalid direction.");
125
126 const int nq0 = m_base[0]->GetNumPoints();
127 const int nq1 = m_base[1]->GetNumPoints();
128 const int nq2 = m_base[2]->GetNumPoints();
129 const int nq = nq0 * nq1 * nq2;
130 const int nm0 = m_base[0]->GetNumModes();
131 const int nm1 = m_base[1]->GetNumModes();
132
133 Array<OneD, NekDouble> alloc(4 * nq + m_ncoeffs + nm0 * nq2 * (nq1 + nm1),
134 0.0);
135 Array<OneD, NekDouble> tmp1(alloc); // Quad metric
136 Array<OneD, NekDouble> tmp2(alloc + nq); // Dir1 metric
137 Array<OneD, NekDouble> tmp3(alloc + 2 * nq); // Dir2 metric
138 Array<OneD, NekDouble> tmp4(alloc + 3 * nq); // Dir3 metric
139 Array<OneD, NekDouble> tmp5(alloc + 4 * nq); // iprod tmp
140 Array<OneD, NekDouble> wsp(tmp5 + m_ncoeffs); // Wsp
141 switch (dir)
142 {
143 case 0:
144 for (int i = 0; i < nq; i++)
145 {
146 tmp2[i] = 1.0;
148 m_base[0]->GetDbdata(), m_base[1]->GetBdata(),
149 m_base[2]->GetBdata(), tmp2, tmp5, wsp, false, true, true);
150
151 tmp2[i] = 0.0;
152
153 for (int j = 0; j < m_ncoeffs; j++)
154 {
155 (*mat)(j, i) = tmp5[j];
156 }
157 }
158 break;
159 case 1:
160 for (int i = 0; i < nq; i++)
161 {
162 tmp2[i] = 1.0;
164 m_base[0]->GetBdata(), m_base[1]->GetDbdata(),
165 m_base[2]->GetBdata(), tmp2, tmp5, wsp, true, false, true);
166
167 tmp2[i] = 0.0;
168
169 for (int j = 0; j < m_ncoeffs; j++)
170 {
171 (*mat)(j, i) = tmp5[j];
172 }
173 }
174 break;
175 case 2:
176 for (int i = 0; i < nq; i++)
177 {
178 tmp2[i] = 1.0;
180 m_base[0]->GetBdata(), m_base[1]->GetBdata(),
181 m_base[2]->GetDbdata(), tmp2, tmp5, wsp, true, true, false);
182 tmp2[i] = 0.0;
183
184 for (int j = 0; j < m_ncoeffs; j++)
185 {
186 (*mat)(j, i) = tmp5[j];
187 }
188 }
189 break;
190 default:
191 NEKERROR(ErrorUtil::efatal, "Not a 2D expansion.");
192 break;
193 }
194}
195
197 const Array<OneD, const NekDouble> &coords,
198 const Array<OneD, const NekDouble> &physvals)
199{
201
202 WARNINGL2(coords[0] >= -1 - NekConstants::kNekZeroTol, "coord[0] < -1");
203 WARNINGL2(coords[0] <= 1 + NekConstants::kNekZeroTol, "coord[0] > 1");
204 WARNINGL2(coords[1] >= -1 - NekConstants::kNekZeroTol, "coord[1] < -1");
205 WARNINGL2(coords[1] <= 1 + NekConstants::kNekZeroTol, "coord[1] > 1");
206 WARNINGL2(coords[2] >= -1 - NekConstants::kNekZeroTol, "coord[2] < -1");
207 WARNINGL2(coords[2] <= 1 + NekConstants::kNekZeroTol, "coord[2] > 1");
208
209 // Obtain local collapsed coordinate from Cartesian coordinate.
210 LocCoordToLocCollapsed(coords, eta);
211
212 const int nq0 = m_base[0]->GetNumPoints();
213 const int nq1 = m_base[1]->GetNumPoints();
214 const int nq2 = m_base[2]->GetNumPoints();
215
216 Array<OneD, NekDouble> wsp1(nq1 * nq2), wsp2(nq2);
217
218 // Construct the 2D square...
219 const NekDouble *ptr = &physvals[0];
220 for (int i = 0; i < nq1 * nq2; ++i, ptr += nq0)
221 {
222 wsp1[i] = StdExpansion::BaryEvaluate<0>(eta[0], ptr);
223 }
224
225 for (int i = 0; i < nq2; ++i)
226 {
227 wsp2[i] = StdExpansion::BaryEvaluate<1>(eta[1], &wsp1[i * nq1]);
228 }
229
230 return StdExpansion::BaryEvaluate<2>(eta[2], &wsp2[0]);
231}
232
235 const Array<OneD, const NekDouble> &physvals)
236{
237 NekDouble value;
238
239 int Qx = m_base[0]->GetNumPoints();
240 int Qy = m_base[1]->GetNumPoints();
241 int Qz = m_base[2]->GetNumPoints();
242
243 Array<OneD, NekDouble> sumFactorization_qr =
244 Array<OneD, NekDouble>(Qy * Qz);
245 Array<OneD, NekDouble> sumFactorization_r = Array<OneD, NekDouble>(Qz);
246
247 // Lagrangian interpolation matrix
248 NekDouble *interpolatingNodes = nullptr;
249
250 // Interpolate first coordinate direction
251 interpolatingNodes = &I[0]->GetPtr()[0];
252
253 Blas::Dgemv('T', Qx, Qy * Qz, 1.0, &physvals[0], Qx, &interpolatingNodes[0],
254 1, 0.0, &sumFactorization_qr[0], 1);
255
256 // Interpolate in second coordinate direction
257 interpolatingNodes = &I[1]->GetPtr()[0];
258
259 Blas::Dgemv('T', Qy, Qz, 1.0, &sumFactorization_qr[0], Qy,
260 &interpolatingNodes[0], 1, 0.0, &sumFactorization_r[0], 1);
261
262 // Interpolate in third coordinate direction
263 interpolatingNodes = &I[2]->GetPtr()[0];
264 value = Blas::Ddot(Qz, interpolatingNodes, 1, &sumFactorization_r[0], 1);
265
266 return value;
267}
268
269/**
270 * @param inarray Input coefficients.
271 * @param output Output coefficients.
272 * @param mkey Matrix key
273 */
275 const Array<OneD, const NekDouble> &inarray,
277{
278 if (mkey.GetNVarCoeff() == 0 &&
281 {
282 // This implementation is only valid when there are no
283 // coefficients associated to the Laplacian operator
284 int nqtot = GetTotPoints();
285
286 const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
287 const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
288 const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
289
290 // Allocate temporary storage
291 Array<OneD, NekDouble> wsp0(7 * nqtot);
292 Array<OneD, NekDouble> wsp1(wsp0 + nqtot);
293
294 if (!(m_base[0]->Collocation() && m_base[1]->Collocation() &&
295 m_base[2]->Collocation()))
296 {
297 // LAPLACIAN MATRIX OPERATION
298 // wsp0 = u = B * u_hat
299 // wsp1 = du_dxi1 = D_xi1 * wsp0 = D_xi1 * u
300 // wsp2 = du_dxi2 = D_xi2 * wsp0 = D_xi2 * u
301 BwdTrans_SumFacKernel(base0, base1, base2, inarray, wsp0, wsp1,
302 true, true, true);
303 LaplacianMatrixOp_MatFree_Kernel(wsp0, outarray, wsp1);
304 }
305 else
306 {
307 LaplacianMatrixOp_MatFree_Kernel(inarray, outarray, wsp1);
308 }
309 }
310 else
311 {
313 mkey);
314 }
315}
316
318 const Array<OneD, const NekDouble> &inarray,
320{
321 if (mkey.GetNVarCoeff() == 0 &&
323 {
324 using std::max;
325
326 int nquad0 = m_base[0]->GetNumPoints();
327 int nquad1 = m_base[1]->GetNumPoints();
328 int nquad2 = m_base[2]->GetNumPoints();
329 int nmodes0 = m_base[0]->GetNumModes();
330 int nmodes1 = m_base[1]->GetNumModes();
331 int nmodes2 = m_base[2]->GetNumModes();
332 int wspsize = max(nquad0 * nmodes2 * (nmodes1 + nquad1),
333 nquad0 * nquad1 * (nquad2 + nmodes0) +
334 nmodes0 * nmodes1 * nquad2);
335
337
338 const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
339 const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
340 const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
341 Array<OneD, NekDouble> wsp0(8 * wspsize);
342 Array<OneD, NekDouble> wsp1(wsp0 + 1 * wspsize);
343 Array<OneD, NekDouble> wsp2(wsp0 + 2 * wspsize);
344
345 if (!(m_base[0]->Collocation() && m_base[1]->Collocation() &&
346 m_base[2]->Collocation()))
347 {
348 // MASS MATRIX OPERATION
349 // The following is being calculated:
350 // wsp0 = B * u_hat = u
351 // wsp1 = W * wsp0
352 // outarray = B^T * wsp1 = B^T * W * B * u_hat = M * u_hat
353 BwdTrans_SumFacKernel(base0, base1, base2, inarray, wsp0, wsp2,
354 true, true, true);
355 MultiplyByQuadratureMetric(wsp0, wsp1);
356 IProductWRTBase_SumFacKernel(base0, base1, base2, wsp1, outarray,
357 wsp2, true, true, true);
358 LaplacianMatrixOp_MatFree_Kernel(wsp0, wsp1, wsp2);
359 }
360 else
361 {
362 // specialised implementation for the classical spectral
363 // element method
364 MultiplyByQuadratureMetric(inarray, outarray);
365 LaplacianMatrixOp_MatFree_Kernel(inarray, wsp1, wsp2);
366 }
367
368 // outarray = lambda * outarray + wsp1
369 // = (lambda * M + L ) * u_hat
370 Vmath::Svtvp(m_ncoeffs, lambda, &outarray[0], 1, &wsp1[0], 1,
371 &outarray[0], 1);
372 }
373 else
374 {
376 mkey);
377 }
378}
379
381 const Array<OneD, const NekDouble> &inarray)
382{
383 const int nqtot = GetTotPoints();
386 return Vmath::Vsum(nqtot, tmp, 1);
387}
388
390{
391 NEKERROR(ErrorUtil::efatal, "This function is not valid or not defined");
392 return 0;
393}
394
395int StdExpansion3D::v_GetEdgeNcoeffs([[maybe_unused]] const int i) const
396{
397 NEKERROR(ErrorUtil::efatal, "This function is not valid or not defined");
398 return 0;
399}
400
402 [[maybe_unused]] const int tid,
403 [[maybe_unused]] Array<OneD, unsigned int> &maparray,
404 [[maybe_unused]] Array<OneD, int> &signarray,
405 [[maybe_unused]] Orientation traceOrient)
406{
407 NEKERROR(ErrorUtil::efatal, "Method does not exist for this shape");
408}
409
412 Array<OneD, int> &signarray,
413 Orientation traceOrient, int P,
414 int Q)
415{
416 Array<OneD, unsigned int> map1, map2;
417 GetTraceCoeffMap(tid, map1);
418 GetElmtTraceToTraceMap(tid, map2, signarray, traceOrient, P, Q);
419
420 if (maparray.size() != map2.size())
421 {
422 maparray = Array<OneD, unsigned int>(map2.size());
423 }
424
425 for (int i = 0; i < map2.size(); ++i)
426 {
427 maparray[i] = map1[map2[i]];
428 }
429}
430
432 [[maybe_unused]] const int facedir,
433 const LibUtilities::BasisType faceDirBasisType, const int numpoints,
434 const int nummodes)
435{
436
437 switch (faceDirBasisType)
438 {
440 {
441 const LibUtilities::PointsKey pkey(
444 pkey);
445 }
448 {
449 const LibUtilities::PointsKey pkey(
452 pkey);
453 }
455 {
456 const LibUtilities::PointsKey pkey(
459 pkey);
460 }
462 {
463 const LibUtilities::PointsKey pkey(
466 pkey);
467 }
470 {
471 const LibUtilities::PointsKey pkey(
474 pkey);
475 }
476 default:
477 {
478 NEKERROR(ErrorUtil::efatal, "expansion type unknown");
479 break;
480 }
481 }
482
483 // Keep things happy by returning a value.
485}
486
488 const int facedir, const LibUtilities::BasisType faceDirBasisType,
489 const int numpoints, const int nummodes, bool UseGLL)
490{
491 switch (faceDirBasisType)
492 {
494 {
495 const LibUtilities::PointsKey pkey(
498 pkey);
499 }
503 {
504 switch (facedir)
505 {
506 case 0:
507 {
508 const LibUtilities::PointsKey pkey(
511 nummodes, pkey);
512 }
513 case 1:
514 {
516
517 if (UseGLL)
518 {
521 }
522 else
523 {
525 numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
526 }
528 nummodes, pkey);
529 }
530 default:
531 {
532
533 NEKERROR(ErrorUtil::efatal, "invalid value to flag");
534 break;
535 }
536 }
537 break;
538 }
539
541 {
542 switch (facedir)
543 {
544 case 0:
545 {
546 const LibUtilities::PointsKey pkey(
549 nummodes, pkey);
550 }
551 case 1:
552 {
553 const LibUtilities::PointsKey pkey(
554 numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
556 nummodes, pkey);
557 }
558 default:
559 {
560 NEKERROR(ErrorUtil::efatal, "invalid value to flag");
561 break;
562 }
563 }
564 break;
565 }
566
571 {
572 switch (facedir)
573 {
574 case 0:
575 {
576 const LibUtilities::PointsKey pkey(
579 nummodes, pkey);
580 }
581 case 1:
582 {
583 const LibUtilities::PointsKey pkey(
584 numpoints, LibUtilities::eGaussRadauMAlpha1Beta0);
586 nummodes, pkey);
587 }
588 default:
589 {
590 NEKERROR(ErrorUtil::efatal, "invalid value to flag");
591 break;
592 }
593 }
594 break;
595 }
596 default:
597 {
598 NEKERROR(ErrorUtil::efatal, "expansion type unknown");
599 break;
600 }
601 }
602
603 // Keep things happy by returning a value.
605}
606
607void StdExpansion3D::v_PhysInterp(std::shared_ptr<StdExpansion> fromExp,
608 const Array<OneD, const NekDouble> &fromData,
610{
611
612 LibUtilities::Interp3D(fromExp->GetBasis(0)->GetPointsKey(),
613 fromExp->GetBasis(1)->GetPointsKey(),
614 fromExp->GetBasis(2)->GetPointsKey(), fromData,
615 m_base[0]->GetPointsKey(), m_base[1]->GetPointsKey(),
616 m_base[2]->GetPointsKey(), toData);
617}
618
619} // namespace Nektar::StdRegions
#define WARNINGL2(condition, msg)
Definition: ErrorUtil.hpp:266
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
Describes the specification for a Basis.
Definition: Basis.h:45
Defines a specification for a set of points.
Definition: Points.h:50
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
void v_HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
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
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.
NekDouble v_PhysEvaluateInterp(const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals) override
NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray) override
Integrates the specified function over the domain.
void v_PhysInterp(std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData) override
void v_GenStdMatBwdDeriv(const int dir, DNekMatSharedPtr &mat) override
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.
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:134
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:708
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:702
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:732
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:124
bool ConstFactorExists(const ConstFactorType &factor) const
Definition: StdMatrixKey.h:133
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 = alpha A x plus beta y where A[m x n].
Definition: Blas.hpp:211
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:163
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:383
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
void Interp3D(const BasisKey &fbasis0, const BasisKey &fbasis1, const BasisKey &fbasis2, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, const BasisKey &tbasis2, Array< OneD, NekDouble > &to)
this function interpolates a 3D function evaluated at the quadrature points of the 3D basis,...
Definition: Interp.cpp:162
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51
@ eModified_B
Principle Modified Functions .
Definition: BasisType.h:49
@ P
Monomial polynomials .
Definition: BasisType.h:62
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:42
@ eModified_C
Principle Modified Functions .
Definition: BasisType.h:50
@ eGLL_Lagrange
Lagrange for SEM basis .
Definition: BasisType.h:56
@ eOrtho_C
Principle Orthogonal Functions .
Definition: BasisType.h:46
@ eModifiedPyr_C
Principle Modified Functions.
Definition: BasisType.h:53
@ eOrtho_B
Principle Orthogonal Functions .
Definition: BasisType.h:44
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48
@ eOrthoPyr_C
Principle Orthogonal Functions .
Definition: BasisType.h:51
static const NekDouble kNekZeroTol
LibUtilities::BasisKey EvaluateTriFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes, bool UseGLL)
LibUtilities::BasisKey EvaluateQuadFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes)
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.hpp:396
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.hpp:608
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825