Nektar++
StdPointExp.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: StdPointExp.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 a Point expansion
32//
33///////////////////////////////////////////////////////////////////////////////
34
36
37using namespace std;
38
39namespace Nektar::StdRegions
40{
41
43 : StdExpansion(Ba.GetNumModes(), 1, Ba),
44 StdExpansion0D(Ba.GetNumModes(), Ba)
45{
46}
47
49{
51}
52
54 [[maybe_unused]] Array<OneD, NekDouble> &coords_1,
55 [[maybe_unused]] Array<OneD, NekDouble> &coords_2)
56{
57 Blas::Dcopy(GetNumPoints(0), (m_base[0]->GetZ()).data(), 1, &coords_0[0],
58 1);
59}
60
62 Array<OneD, NekDouble> &outarray)
63{
64 int nquad = m_base[0]->GetNumPoints();
65
66 if (m_base[0]->Collocation())
67 {
68 Vmath::Vcopy(nquad, inarray, 1, outarray, 1);
69 }
70 else
71 {
72 Blas::Dgemv('N', nquad, m_base[0]->GetNumModes(), 1.0,
73 (m_base[0]->GetBdata()).data(), nquad, &inarray[0], 1, 0.0,
74 &outarray[0], 1);
75 }
76}
77
79 Array<OneD, NekDouble> &outarray)
80{
81 if (m_base[0]->Collocation())
82 {
83 Vmath::Vcopy(m_ncoeffs, inarray, 1, outarray, 1);
84 }
85 else
86 {
87 v_IProductWRTBase(inarray, outarray);
88
89 // get Mass matrix inverse
90 StdMatrixKey masskey(eInvMass, v_DetShapeType(), *this);
91 DNekMatSharedPtr matsys = GetStdMatrix(masskey);
92
95
96 out = (*matsys) * in;
97 }
98}
99
101 Array<OneD, NekDouble> &outarray)
102{
103 v_BwdTrans(inarray, outarray);
104}
105
106// Inner product
108 const Array<OneD, const NekDouble> &inarray,
109 Array<OneD, NekDouble> &outarray,
110 int coll_check)
111{
112 int nquad = m_base[0]->GetNumPoints();
113 Array<OneD, NekDouble> tmp(nquad);
116
117 Vmath::Vmul(nquad, inarray, 1, w, 1, tmp, 1);
118
119 if (coll_check && m_base[0]->Collocation())
120 {
121 Vmath::Vcopy(nquad, tmp, 1, outarray, 1);
122 }
123 else
124 {
125 Blas::Dgemv('T', nquad, m_ncoeffs, 1.0, base.data(), nquad, &tmp[0], 1,
126 0.0, outarray.data(), 1);
127 }
128}
129
131 Array<OneD, NekDouble> &outarray)
132{
133 v_IProductWRTBase(m_base[0]->GetBdata(), inarray, outarray, 1);
134}
135
137 [[maybe_unused]] const int dir, const Array<OneD, const NekDouble> &inarray,
138 Array<OneD, NekDouble> &outarray)
139{
140 ASSERTL1(dir >= 0 && dir < 1, "input dir is out of range");
141 v_IProductWRTBase(m_base[0]->GetDbdata(), inarray, outarray, 1);
142}
143
145 const Array<OneD, const NekDouble> &inarray,
146 Array<OneD, NekDouble> &outarray, [[maybe_unused]] bool multiplybyweights)
147{
148 v_IProductWRTBase(m_base[0]->GetBdata(), inarray, outarray, 1);
149}
150
152{
154 MatrixType mattype;
155
156 switch (mattype = mkey.GetMatrixType())
157 {
158 case eFwdTrans:
159 {
160 Mat =
162 StdMatrixKey iprodkey(eIProductWRTBase, v_DetShapeType(), *this);
163 DNekMat &Iprod = *GetStdMatrix(iprodkey);
164 StdMatrixKey imasskey(eInvMass, v_DetShapeType(), *this);
165 DNekMat &Imass = *GetStdMatrix(imasskey);
166
167 (*Mat) = Imass * Iprod;
168 }
169 break;
170 default:
171 {
173 }
174 break;
175 }
176
177 return Mat;
178}
179
181{
182 return v_GenMatrix(mkey);
183}
184
185} // namespace Nektar::StdRegions
#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
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:65
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:612
DNekMatSharedPtr CreateGeneralMatrix(const StdMatrixKey &mkey)
this function generates the mass matrix
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
Definition: StdExpansion.h:218
Array< OneD, LibUtilities::BasisSharedPtr > m_base
MatrixType GetMatrixType() const
Definition: StdMatrixKey.h:83
DNekMatSharedPtr v_CreateStdMatrix(const StdMatrixKey &mkey) override
void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Calculates the inner product of a given function f with the different modes of the expansion.
void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
DNekMatSharedPtr v_GenMatrix(const StdMatrixKey &mkey) override
void v_BwdTrans_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
LibUtilities::ShapeType v_DetShapeType() const override
Definition: StdPointExp.cpp:48
void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Definition: StdPointExp.cpp:61
void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Transform a given function from physical quadrature space to coefficient space.
Definition: StdPointExp.cpp:78
void v_GetCoords(Array< OneD, NekDouble > &coords_0, Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2) override
Definition: StdPointExp.cpp:53
void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
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 void Dcopy(const int &n, const double *x, const int &incx, double *y, const int &incy)
BLAS level 1: Copy x to y.
Definition: Blas.hpp:128
std::vector< double > w(NPUPPER)
std::vector< double > z(NPUPPER)
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
STL namespace.