47 namespace MultiRegions
50 ExpListHomogeneous1D::ExpListHomogeneous1D():
63 const bool dealiasing):
68 m_dealiasing(dealiasing)
71 "Homogeneous Basis is a null basis");
76 m_comm->GetColumnComm()->GetColumnComm() :
81 "HomModesZ should be a multiple of npz.");
91 "HomModesZ/npz should be an even integer.");
120 ASSERTL0(
false,
"Dealiasing available just in combination "
132 m_transposition(In.m_transposition),
133 m_StripZcomm(In.m_StripZcomm),
134 m_useFFT(In.m_useFFT),
137 m_tmpOUT(In.m_tmpOUT),
138 m_homogeneousBasis(In.m_homogeneousBasis),
140 m_homogeneous1DBlockMat(In.m_homogeneous1DBlockMat),
141 m_dealiasing(In.m_dealiasing),
142 m_padsize(In.m_padsize)
148 const std::vector<unsigned int> &eIDs):
150 m_transposition(In.m_transposition),
151 m_useFFT(In.m_useFFT),
154 m_tmpOUT(In.m_tmpOUT),
155 m_homogeneousBasis(In.m_homogeneousBasis),
158 m_dealiasing(In.m_dealiasing),
159 m_padsize(In.m_padsize)
203 int num_dofs = inarray1.num_elements();
221 int num_points_per_plane = num_dofs/
m_planes.num_elements();
223 if(!
m_session->DefinesSolverInfo(
"HomoStrip"))
225 num_proc =
m_comm->GetColumnComm()->GetSize();
231 int num_dfts_per_proc = num_points_per_plane / num_proc
232 + (num_points_per_plane % num_proc > 0);
250 for(
int i = 0 ; i < num_dfts_per_proc ; i++)
254 Vmath::Vcopy(N, &(ShufV1[i*N]), 1, &(ShufV1_PAD_coef[0]), 1);
255 Vmath::Vcopy(N, &(ShufV2[i*N]), 1, &(ShufV2_PAD_coef[0]), 1);
258 m_FFT_deal->FFTBwdTrans(ShufV1_PAD_coef, ShufV1_PAD_phys);
259 m_FFT_deal->FFTBwdTrans(ShufV2_PAD_coef, ShufV2_PAD_phys);
265 ShufV1V2_PAD_phys, 1);
269 m_FFT_deal->FFTFwdTrans(ShufV1V2_PAD_phys, ShufV1V2_PAD_coef);
274 &(ShufV1V2[i*N]), 1);
306 int ndim = inarray1.num_elements();
307 ASSERTL1( inarray2.num_elements() % ndim == 0,
308 "Wrong dimensions for DealiasedDotProd.");
309 int nvec = inarray2.num_elements() / ndim;
311 int num_dofs = inarray1[0].num_elements();
314 int num_points_per_plane = num_dofs/
m_planes.num_elements();
316 if(!
m_session->DefinesSolverInfo(
"HomoStrip"))
318 num_proc =
m_comm->GetColumnComm()->GetSize();
324 int num_dfts_per_proc = num_points_per_plane / num_proc
325 + (num_points_per_plane % num_proc > 0);
332 for (
int i = 0; i < ndim; i++)
336 for (
int i = 0; i < ndim*nvec; i++)
343 for (
int i = 0; i < ndim; i++)
348 for (
int i = 0; i < ndim*nvec; i++)
359 for (
int i = 0; i < ndim; i++)
362 (num_dfts_per_proc*N,0.0);
370 for (
int i = 0; i < ndim*nvec; i++)
373 (num_dfts_per_proc*N,0.0);
381 for (
int i = 0; i < nvec; i++)
384 (num_dfts_per_proc*N,0.0);
388 for (
int i = 0; i < ndim; i++)
393 for (
int i = 0; i < ndim*nvec; i++)
400 for(
int i = 0 ; i < num_dfts_per_proc ; i++)
402 for (
int j = 0; j < ndim; j++)
407 &(ShufV1_PAD_coef[0]), 1);
409 m_FFT_deal->FFTBwdTrans(ShufV1_PAD_coef, ShufV1_PAD_phys[j]);
411 for (
int j = 0; j < ndim*nvec; j++)
414 &(ShufV2_PAD_coef[0]), 1);
415 m_FFT_deal->FFTBwdTrans(ShufV2_PAD_coef, ShufV2_PAD_phys[j]);
420 for (
int j = 0; j < nvec; j++)
423 for (
int k = 0; k < ndim; k++)
426 ShufV2_PAD_phys[j*ndim+k], 1,
427 ShufV1V2_PAD_phys, 1,
428 ShufV1V2_PAD_phys, 1);
432 m_FFT_deal->FFTFwdTrans(ShufV1V2_PAD_phys, ShufV1V2_PAD_coef);
436 &(ShufV1V2[j][i*N]), 1);
443 for (
int j = 0; j < nvec; j++)
453 for (
int j = 0; j < nvec; j++)
467 int cnt = 0, cnt1 = 0;
470 for(
int n = 0; n <
m_planes.num_elements(); ++n)
472 m_planes[n]->FwdTrans(inarray+cnt, tmparray = outarray + cnt1,
489 int cnt = 0, cnt1 = 0;
493 for(
int n = 0; n <
m_planes.num_elements(); ++n)
495 m_planes[n]->FwdTrans_IterPerExp(inarray+cnt, tmparray = outarray + cnt1);
510 int cnt = 0, cnt1 = 0;
513 for(
int n = 0; n <
m_planes.num_elements(); ++n)
515 m_planes[n]->BwdTrans(inarray+cnt, tmparray = outarray + cnt1,
518 cnt1 +=
m_planes[n]->GetTotPoints();
531 int cnt = 0, cnt1 = 0;
534 for(
int n = 0; n <
m_planes.num_elements(); ++n)
536 m_planes[n]->BwdTrans_IterPerExp(inarray+cnt, tmparray = outarray + cnt1);
539 cnt1 +=
m_planes[n]->GetTotPoints();
552 int cnt = 0, cnt1 = 0;
565 for(
int n = 0; n <
m_planes.num_elements(); ++n)
567 m_planes[n]->IProductWRTBase(tmpIn+cnt, tmparray = outarray + cnt1,coeffstate);
579 int cnt = 0, cnt1 = 0;
592 for(
int n = 0; n <
m_planes.num_elements(); ++n)
594 m_planes[n]->IProductWRTBase_IterPerExp(tmpIn+cnt, tmparray = outarray + cnt1);
614 num_dofs = inarray.num_elements();
618 num_dofs = outarray.num_elements();
623 int num_points_per_plane = num_dofs/
m_planes.num_elements();
624 int num_dfts_per_proc;
625 if(!
m_session->DefinesSolverInfo(
"HomoStrip"))
627 int nP =
m_comm->GetColumnComm()->GetSize();
628 num_dfts_per_proc = num_points_per_plane / nP
629 + (num_points_per_plane % nP > 0);
634 num_dfts_per_proc = num_points_per_plane / nP
635 + (num_points_per_plane % nP > 0);
653 for(
int i = 0 ; i < num_dfts_per_proc ; i++)
660 for(
int i = 0 ; i < num_dfts_per_proc ; i++)
673 fft_out,1,outarray,1);
703 int nrows = blkmat->GetRows();
704 int ncols = blkmat->GetColumns();
748 return matrixIter->second;
758 int num_trans_per_proc = 0;
767 n_exp =
m_planes[0]->GetTotPoints();
770 num_trans_per_proc = n_exp/
m_comm->GetColumnComm()->GetSize() + (n_exp%
m_comm->GetColumnComm()->GetSize() > 0);
798 StdPoint.DetShapeType(),
801 loc_mat = StdPoint.GetStdMatrix(matkey);
806 StdPoint.DetShapeType(),
809 loc_mat = StdPoint.GetStdMatrix(matkey);
820 StdSeg.DetShapeType(),
823 loc_mat = StdSeg.GetStdMatrix(matkey);
828 StdSeg.DetShapeType(),
831 loc_mat = StdSeg.GetStdMatrix(matkey);
837 for(
int i = 0; i < num_trans_per_proc; ++i)
839 BlkMatrix->SetBlock(i,i,loc_mat);
847 std::vector<LibUtilities::FieldDefinitionsSharedPtr> returnval;
852 std::vector<NekDouble> HomoLen;
853 HomoLen.push_back(
m_lhom);
855 std::vector<unsigned int> StripsIDs;
858 m_session->MatchSolverInfo(
"HomoStrip",
"True",strips,
false);
864 std::vector<unsigned int> PlanesIDs;
874 for(
int i = 0; i <
m_planes.num_elements(); i++)
879 m_planes[0]->GeneralGetFieldDefinitions(returnval, 1, HomoBasis,
880 HomoLen, strips, StripsIDs, PlanesIDs);
889 std::vector<NekDouble> HomoLen;
890 HomoLen.push_back(
m_lhom);
892 std::vector<unsigned int> StripsIDs;
895 m_session->MatchSolverInfo(
"HomoStrip",
"True",strips,
false);
901 std::vector<unsigned int> PlanesIDs;
909 for(
int i = 0; i <
m_planes.num_elements(); i++)
915 m_planes[0]->GeneralGetFieldDefinitions(fielddef, 1, HomoBasis,
916 HomoLen, strips, StripsIDs, PlanesIDs);
929 int ncoeffs_per_plane =
m_planes[0]->GetNcoeffs();
935 for(i = 0; i <
m_planes[0]->GetExpSize(); ++i)
941 for(i = 0; i < fielddef->m_elementIDs.size(); ++i)
944 int datalen = (*m_exp)[eid]->GetNcoeffs();
946 for(n = 0; n <
m_planes.num_elements(); ++n)
948 fielddata.insert(fielddata.end(),&coeffs[
m_coeff_offset[eid]+n*ncoeffs_per_plane],&coeffs[m_coeff_offset[eid]+n*ncoeffs_per_plane]+datalen);
961 std::vector<NekDouble> &fielddata,
968 int datalen = fielddata.size()/fielddef->m_fields.size();
969 std::vector<unsigned int> fieldDefHomoZids;
973 for(i = 0; i < fielddef->m_fields.size(); ++i)
975 if(fielddef->m_fields[i] == field)
982 if(i == fielddef->m_fields.size())
984 cout <<
"Field "<< field<<
"not found in data file. " << endl;
989 int modes_offset = 0;
990 int planes_offset = 0;
998 for (i = 0; i <
m_planes.num_elements(); ++i)
1003 for (i = 0; i <
m_planes[0]->GetExpSize(); ++i)
1009 if(fielddef->m_numHomogeneousDir)
1011 nzmodes = fielddef->m_homogeneousZIDs.size();
1012 fieldDefHomoZids = fielddef->m_homogeneousZIDs;
1017 fieldDefHomoZids.push_back(0);
1021 int ncoeffs_per_plane =
m_planes[0]->GetNcoeffs();
1023 for(i = 0; i < fielddef->m_elementIDs.size(); ++i)
1025 if(fielddef->m_uniOrder ==
true)
1031 fielddef->m_numModes,
1040 offset += datalen*nzmodes;
1041 modes_offset += (*m_exp)[0]->GetNumBases() +
1042 fielddef->m_numHomogeneousDir;
1046 int eid = it->second;
1049 for(n = 0; n < nzmodes; ++n, offset += datalen)
1060 planes_offset = it->second;
1061 if(datalen == (*
m_exp)[eid]->GetNcoeffs())
1067 (*m_exp)[eid]->ExtractDataToCoeffs(&fielddata[offset], fielddef->m_numModes,modes_offset,&coeffs[
m_coeff_offset[eid] + planes_offset*ncoeffs_per_plane], fielddef->m_basis);
1070 modes_offset += (*m_exp)[0]->GetNumBases() + fielddef->m_numHomogeneousDir;
1080 int fromNcoeffs_per_plane = fromExpList->GetPlane(0)->GetNcoeffs();
1081 int toNcoeffs_per_plane =
m_planes[0]->GetNcoeffs();
1084 for(i = 0; i <
m_planes.num_elements(); ++i)
1086 m_planes[i]->ExtractCoeffsToCoeffs(fromExpList->GetPlane(i),fromcoeffs_tmp = fromCoeffs + fromNcoeffs_per_plane*i, tocoeffs_tmp = toCoeffs + toNcoeffs_per_plane*i);
1096 m_planes[0]->WriteVtkPieceData(outfile, expansion, var);
1101 int nq = (*m_exp)[expansion]->GetTotPoints();
1102 int npoints_per_plane =
m_planes[0]->GetTotPoints();
1105 int outputExtraPlane = 0;
1111 outputExtraPlane = 1;
1122 int fromRank = (rank+1) % size;
1123 int toRank = (rank == 0) ? size-1 : rank-1;
1129 fromRank, extraPlane);
1134 outfile <<
" <DataArray type=\"Float64\" Name=\""
1135 << var <<
"\">" << endl;
1137 for (
int n = 0; n <
m_planes.num_elements(); ++n)
1140 for(i = 0; i < nq; ++i)
1145 if (outputExtraPlane)
1147 for(i = 0; i < nq; ++i)
1150 0 : extraPlane[i]) <<
" ";
1154 outfile <<
" </DataArray>" << endl;
1162 cnt1 =
m_planes[0]->Get1DScaledTotPoints(scale);
1164 ASSERTL1(
m_planes.num_elements()*cnt1 <= outarray.num_elements(),
"size of outarray does not match internal estimage");
1167 for(
int i = 0; i <
m_planes.num_elements(); i++)
1170 m_planes[i]->PhysInterp1DScaled(scale,inarray+i*cnt,
1171 tmparray = outarray+i*cnt1);
1180 cnt =
m_planes[0]->Get1DScaledTotPoints(scale);
1181 cnt1 =
m_planes[0]->GetTotPoints();
1183 ASSERTL1(
m_planes.num_elements()*cnt <= inarray.num_elements(),
"size of outarray does not match internal estimage");
1186 for(
int i = 0; i <
m_planes.num_elements(); i++)
1188 m_planes[i]->PhysGalerkinProjection1DScaled(scale,inarray+i*cnt,
1189 tmparray = outarray+i*cnt1);
1198 int nT_pts = inarray.num_elements();
1199 int nP_pts = nT_pts/
m_planes.num_elements();
1207 for(
int i = 0; i <
m_planes.num_elements(); i++)
1209 m_planes[i]->PhysDeriv(inarray + i*nP_pts ,tmp2 = out_d0 + i*nP_pts , tmp3 = out_d1 + i*nP_pts );
1219 temparray = inarray;
1246 for(
int i = 0; i <
m_planes.num_elements(); i++)
1250 Vmath::Smul(nP_pts,beta,tmp1 = temparray + i*nP_pts,1,tmp2 = outarray + (i-
int(sign))*nP_pts,1);
1267 if(!
m_session->DefinesSolverInfo(
"HomoStrip"))
1270 "Parallelisation in the homogeneous direction "
1271 "implemented just for Fourier basis");
1276 "Parallelisation in the homogeneous direction "
1277 "implemented just for Fourier basis");
1282 ASSERTL0(
false,
"Semi-phyisical time-stepping not "
1283 "implemented yet for non-Fourier "
1292 for(
int i = 0; i < nP_pts; i++)
1294 StdSeg.PhysDeriv(temparray + i*
m_planes.num_elements(), tmp2 = outarray + i*
m_planes.num_elements());
1309 int nT_pts = inarray.num_elements();
1310 int nP_pts = nT_pts/
m_planes.num_elements();
1321 for(
int i=0; i <
m_planes.num_elements(); i++)
1323 m_planes[i]->PhysDeriv(edir, inarray + i*nP_pts ,tmp2 = out_d + i*nP_pts);
1333 temparray = inarray;
1359 for(
int i = 0; i <
m_planes.num_elements(); i++)
1363 Vmath::Smul(nP_pts,beta,tmp1 = temparray + i*nP_pts,1,tmp2 = outarray + (i-
int(sign))*nP_pts,1);
1379 if(!
m_session->DefinesSolverInfo(
"HomoStrip"))
1382 "Parallelisation in the homogeneous direction "
1383 "implemented just for Fourier basis");
1388 "Parallelisation in the homogeneous direction "
1389 "implemented just for Fourier basis");
1394 ASSERTL0(
false,
"Semi-phyisical time-stepping not implemented yet for non-Fourier basis");
1402 for(
int i = 0; i < nP_pts; i++)
1404 StdSeg.PhysDeriv(temparray + i*
m_planes.num_elements(), tmp2 = outarray + i*
m_planes.num_elements());
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
void HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Homo1DBlockMatrixMapShPtr m_homogeneous1DBlockMat
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
#define ASSERTL0(condition, msg)
DNekBlkMatSharedPtr GetHomogeneous1DBlockMatrix(Homogeneous1DMatType mattype, CoeffState coeffstate=eLocal) const
virtual void v_PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
boost::shared_ptr< Transposition > TranspositionSharedPtr
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
static Array< OneD, NekDouble > NullNekDouble1DArray
boost::unordered_map< int, int > m_elmtToExpId
Mapping from geometry ID of element to index inside m_exp.
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
#define sign(a, b)
return the sign(b)*a
LibUtilities::TranspositionSharedPtr m_transposition
virtual void v_WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var)
virtual void v_PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
boost::shared_ptr< FieldDefinitions > FieldDefinitionsSharedPtr
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
LibUtilities::NektarFFTSharedPtr m_FFT_deal
static BasisSharedPtr NullBasisSharedPtr
Array< OneD, NekDouble > m_tmpOUT
NekDouble m_lhom
Width of homogeneous direction.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
virtual void v_HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
boost::shared_ptr< DNekMat > DNekMatSharedPtr
NektarFFTFactory & GetNektarFFTFactory()
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
virtual LibUtilities::TranspositionSharedPtr v_GetTransposition(void)
int GetNumberOfCoefficients(ShapeType shape, std::vector< unsigned int > &modes, int offset)
virtual void v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
BasisManagerT & BasisManager(void)
Class representing a segment element in reference space.
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
LibUtilities::BasisSharedPtr m_homogeneousBasis
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
Base class for all multi-elemental spectral/hp expansions.
boost::unordered_map< int, int > m_zIdToPlane
1D Evenly-spaced points using Fourier Fit
static const NekDouble kNekZeroTol
Fourier Modified expansions with just the real part of the first mode .
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
virtual void v_BwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
LibUtilities::SessionReaderSharedPtr m_session
Session.
Array< OneD, ExpListSharedPtr > m_planes
void HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
virtual void v_FwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
ExpListHomogeneous1D()
Default constructor.
Array< OneD, NekDouble > m_tmpIN
Fourier Modified expansions with just the imaginary part of the first mode .
virtual Array< OneD, const unsigned int > v_GetZIDs(void)
virtual void v_DealiasedDotProd(const Array< OneD, Array< OneD, NekDouble > > &inarray1, const Array< OneD, Array< OneD, NekDouble > > &inarray2, Array< OneD, Array< OneD, NekDouble > > &outarray, CoeffState coeffstate=eLocal)
virtual void v_ExtractDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract data from raw field data into expansion list.
boost::shared_ptr< DNekBlkMat > DNekBlkMatSharedPtr
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
virtual void v_DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
LibUtilities::CommSharedPtr m_comm
Communicator.
void Homogeneous1DTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool IsForwards, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
virtual std::vector< LibUtilities::FieldDefinitionsSharedPtr > v_GetFieldDefinitions(void)
DNekBlkMatSharedPtr GenHomogeneous1DBlockMatrix(Homogeneous1DMatType mattype, CoeffState coeffstate=eLocal) const
virtual void v_IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Fourier ModifiedExpansion with just the first mode .
virtual ~ExpListHomogeneous1D()
Destructor.
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
LibUtilities::CommSharedPtr m_StripZcomm
bool m_useFFT
FFT variables.
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
std::map< Homogeneous1DMatType, DNekBlkMatSharedPtr > Homo1DBlockMatrixMap
A map between homo matrix keys and their associated block matrices.
virtual NekDouble v_GetHomoLen(void)
void Zero(int n, T *x, const int incx)
Zero vector.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
Describes the specification for a Basis.
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.
LibUtilities::NektarFFTSharedPtr m_FFT