35 #include <boost/core/ignore_unused.hpp>
47 namespace MultiRegions
53 m_lhom_z(1), m_homogeneous2DBlockMat(
63 const NekDouble lhom_z,
const bool useFFT,
const bool dealiasing)
64 :
ExpList(type), m_useFFT(useFFT), m_lhom_y(lhom_y), m_lhom_z(lhom_z),
65 m_homogeneous2DBlockMat(
67 m_dealiasing(dealiasing)
70 m_comm = pSession->GetComm();
73 "Homogeneous Basis in y direction is a null basis");
75 "Homogeneous Basis in z direction is a null basis");
82 HomoBasis_y, HomoBasis_z,
m_comm->GetColumnComm());
103 "Remove dealiasing if you want to run in parallel");
112 :
ExpList(In, false), m_useFFT(In.m_useFFT), m_FFT_y(In.m_FFT_y),
113 m_FFT_z(In.m_FFT_z), m_transposition(In.m_transposition),
114 m_Ycomm(In.m_Ycomm), m_Zcomm(In.m_Ycomm),
115 m_homogeneousBasis_y(In.m_homogeneousBasis_y),
116 m_homogeneousBasis_z(In.m_homogeneousBasis_z), m_lhom_y(In.m_lhom_y),
117 m_lhom_z(In.m_lhom_z),
118 m_homogeneous2DBlockMat(In.m_homogeneous2DBlockMat), m_ny(In.m_ny),
119 m_nz(In.m_nz), m_dealiasing(In.m_dealiasing), m_padsize_y(In.m_padsize_y),
120 m_padsize_z(In.m_padsize_z), MatFwdPAD(In.MatFwdPAD),
121 MatBwdPAD(In.MatBwdPAD)
128 :
ExpList(In, eIDs, false), m_useFFT(In.m_useFFT), m_FFT_y(In.m_FFT_y),
129 m_FFT_z(In.m_FFT_z), m_transposition(In.m_transposition),
130 m_Ycomm(In.m_Ycomm), m_Zcomm(In.m_Ycomm),
131 m_homogeneousBasis_y(In.m_homogeneousBasis_y),
132 m_homogeneousBasis_z(In.m_homogeneousBasis_z), m_lhom_y(In.m_lhom_y),
133 m_lhom_z(In.m_lhom_z),
134 m_homogeneous2DBlockMat(
136 m_ny(In.m_ny), m_nz(In.m_nz), m_dealiasing(In.m_dealiasing),
137 m_padsize_y(In.m_padsize_y), m_padsize_z(In.m_padsize_z),
138 MatFwdPAD(In.MatFwdPAD), MatBwdPAD(In.MatBwdPAD)
170 int npoints = outarray.size();
174 int nslabs = npoints /
222 for (
int j = 0; j < nslabs; j++)
226 for (
int i = 0; i <
m_nz; i++)
229 &(PadV1_slab_coeff[i * 2 *
m_ny]), 1);
231 &(PadV2_slab_coeff[i * 2 *
m_ny]), 1);
235 PadOUT_V1 = (*MatBwdPAD) * PadIN_V1;
236 PadOUT_V2 = (*MatBwdPAD) * PadIN_V2;
241 PadV2_slab_phys, 1, PadRe_slab_phys, 1);
245 PadOUT_Re = (*MatFwdPAD) * PadIN_Re;
249 for (
int i = 0; i <
m_nz; i++)
252 &(ShufV1V2[i *
m_ny + j * nlines]), 1);
276 int ndim = inarray1.size();
277 ASSERTL1(inarray2.size() % ndim == 0,
278 "Wrong dimensions for DealiasedDotProd.");
279 int nvec = inarray2.size() % ndim;
280 int npts = inarray1[0].size();
283 for (
int i = 0; i < nvec; i++)
286 for (
int j = 0; j < ndim; j++)
289 Vmath::Vadd(npts, outarray[i], 1, out, 1, outarray[i], 1);
298 int cnt = 0, cnt1 = 0;
302 for (
int n = 0; n < nlines; ++n)
304 m_lines[n]->FwdTrans(inarray + cnt, tmparray = outarray + cnt1);
305 cnt +=
m_lines[n]->GetTotPoints();
306 cnt1 +=
m_lines[n]->GetNcoeffs();
318 int cnt = 0, cnt1 = 0;
322 for (
int n = 0; n < nlines; ++n)
324 m_lines[n]->FwdTransLocalElmt(inarray + cnt,
325 tmparray = outarray + cnt1);
327 cnt +=
m_lines[n]->GetTotPoints();
328 cnt1 +=
m_lines[n]->GetNcoeffs();
340 int cnt = 0, cnt1 = 0;
344 for (
int n = 0; n < nlines; ++n)
346 m_lines[n]->BwdTrans(inarray + cnt, tmparray = outarray + cnt1);
347 cnt +=
m_lines[n]->GetNcoeffs();
348 cnt1 +=
m_lines[n]->GetTotPoints();
360 int cnt = 0, cnt1 = 0;
364 for (
int n = 0; n < nlines; ++n)
366 m_lines[n]->IProductWRTBase(inarray + cnt, tmparray = outarray + cnt1);
368 cnt +=
m_lines[n]->GetNcoeffs();
369 cnt1 +=
m_lines[n]->GetTotPoints();
377 boost::ignore_unused(Shuff, UnShuff);
384 int s = inarray.size();
398 for (
int i = 0; i < (
p *
m_nz); i++)
406 for (
int i = 0; i < (
p *
m_nz); i++)
418 for (
int i = 0; i < (
p *
m_ny); i++)
426 for (
int i = 0; i < (
p *
m_ny); i++)
472 int nrowsY = blkmatY->GetRows();
473 int ncolsY = blkmatY->GetColumns();
478 int nrowsZ = blkmatZ->GetRows();
479 int ncolsZ = blkmatZ->GetColumns();
493 outY = (*blkmatY) * inY;
498 outZ = (*blkmatZ) * inZ;
520 return matrixIter->second;
562 n_exp =
m_lines[0]->GetNcoeffs();
566 n_exp =
m_lines[0]->GetTotPoints();
611 for (i = 0; i < (n_exp * NumPencils); ++i)
613 BlkMatrix->SetBlock(i, i, loc_mat);
622 std::vector<LibUtilities::FieldDefinitionsSharedPtr> returnval;
628 std::vector<NekDouble> HomoLen(2);
637 std::vector<unsigned int> yIDs;
638 std::vector<unsigned int> zIDs;
640 for (
int n = 0; n < nhom_modes_z; ++n)
642 for (
int m = 0; m < nhom_modes_y; ++m)
649 m_lines[0]->GeneralGetFieldDefinitions(returnval, 2, HomoBasis, HomoLen,
650 false, sIDs, zIDs, yIDs);
655 std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef)
661 std::vector<NekDouble> HomoLen(2);
670 std::vector<unsigned int> yIDs;
671 std::vector<unsigned int> zIDs;
673 for (
int n = 0; n < nhom_modes_z; ++n)
675 for (
int m = 0; m < nhom_modes_y; ++m)
683 m_lines[0]->GeneralGetFieldDefinitions(fielddef, 2, HomoBasis, HomoLen,
684 false, sIDs, zIDs, yIDs);
696 int ncoeffs_per_line =
m_lines[0]->GetNcoeffs();
700 map<int, int> ElmtID_to_ExpID;
701 for (i = 0; i <
m_lines[0]->GetExpSize(); ++i)
703 ElmtID_to_ExpID[(*m_exp)[i]->GetGeom()->GetGlobalID()] = i;
706 for (i = 0; i < fielddef->m_elementIDs.size(); ++i)
708 int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
709 int datalen = (*m_exp)[eid]->GetNcoeffs();
711 for (k = 0; k < (NumMod_y * NumMod_z); ++k)
723 std::vector<NekDouble> &fielddata)
731 std::vector<NekDouble> &fielddata, std::string &field,
736 int datalen = fielddata.size() / fielddef->m_fields.size();
737 int ncoeffs_per_line =
m_lines[0]->GetNcoeffs();
742 for (i = 0; i < fielddef->m_fields.size(); ++i)
744 if (fielddef->m_fields[i] == field)
751 ASSERTL0(i != fielddef->m_fields.size(),
"Field not found in data file");
755 map<int, int> ElmtID_to_ExpID;
756 for (i = 0; i <
m_lines[0]->GetExpSize(); ++i)
758 ElmtID_to_ExpID[(*m_exp)[i]->GetGeom()->GetGlobalID()] = i;
761 for (i = 0; i < fielddef->m_elementIDs.size(); ++i)
763 int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
764 int datalen = (*m_exp)[eid]->GetNcoeffs();
766 for (k = 0; k < (NumMod_y * NumMod_z); ++k)
777 int expansion, std::string var)
780 int nq = (*m_exp)[expansion]->GetTotPoints();
781 int npoints_per_line =
m_lines[0]->GetTotPoints();
784 outfile <<
" <DataArray type=\"Float64\" Name=\"" << var <<
"\">"
787 for (
int n = 0; n <
m_lines.size(); ++n)
791 for (i = 0; i < nq; ++i)
798 outfile <<
" </DataArray>" << endl;
808 int npoints = inarray.size();
810 int n_points_line = npoints / nyzlines;
819 for (
int i = 0; i < nyzlines; i++)
821 m_lines[i]->PhysDeriv(tmp1 = inarray + i * n_points_line,
822 tmp2 = out_d0 + i * n_points_line);
840 for (
int i = 0; i <
m_ny; i++)
844 for (
int j = 0; j <
m_nz; j++)
847 tmp1 = temparray + n_points_line * (i + j *
m_ny),
850 n_points_line * ((i -
int(
sign)) + j *
m_ny),
859 for (
int i = 0; i <
m_nz; i++)
864 tmp1 = temparray + i *
m_ny * n_points_line, 1,
865 tmp2 = temparray2 + (i -
int(
sign)) *
m_ny * n_points_line, 1);
883 ASSERTL0(
false,
"Semi-phyisical time-stepping not implemented yet "
884 "for non-Fourier basis")
894 for (
int i = 0; i < n_points_line; i++)
896 StdQuad.
PhysDeriv(tmp1 = temparray + i * nyzlines,
897 tmp2 = temparray1 + i * nyzlines,
898 tmp3 = temparray2 + i * nyzlines);
918 int npoints = inarray.size();
920 int n_points_line = npoints / nyzlines;
933 for (
int i = 0; i < nyzlines; i++)
935 m_lines[i]->PhysDeriv(tmp1 = inarray + i * n_points_line,
936 tmp2 = out_d + i * n_points_line);
958 for (
int i = 0; i <
m_ny; i++)
962 for (
int j = 0; j <
m_nz; j++)
966 tmp1 = temparray + n_points_line * (i + j *
m_ny),
969 n_points_line * ((i -
int(
sign)) + j *
m_ny),
986 for (
int i = 0; i <
m_nz; i++)
990 tmp1 = temparray + i *
m_ny * n_points_line, 1,
992 (i -
int(
sign)) *
m_ny * n_points_line,
1010 ASSERTL0(
false,
"Semi-phyisical time-stepping not implemented "
1011 "yet for non-Fourier basis")
1022 for (
int i = 0; i < n_points_line; i++)
1024 StdQuad.
PhysDeriv(tmp1 = temparray + i * nyzlines,
1025 tmp2 = temparray1 + i * nyzlines,
1026 tmp3 = temparray2 + i * nyzlines);
#define ASSERTL0(condition, msg)
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
#define sign(a, b)
return the sign(b)*a
Represents a basis of a given type.
const BasisKey GetBasisKey() const
Returns the specification for the Basis.
Describes the specification for a Basis.
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Defines a specification for a set of points.
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
virtual void v_FwdTransLocalElmt(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool Shuff=true, bool UnShuff=true)
DNekBlkMatSharedPtr GenHomogeneous2DBlockMatrix(Homogeneous2DMatType mattype) const
virtual std::vector< LibUtilities::FieldDefinitionsSharedPtr > v_GetFieldDefinitions(void)
DNekMatSharedPtr MatFwdPAD
LibUtilities::NektarFFTSharedPtr m_FFT_y
virtual void v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
ExpListHomogeneous2D(const ExpansionType type)
Default constructor.
int m_nz
Number of modes = number of poitns in z direction.
Homo2DBlockMatrixMapShPtr m_homogeneous2DBlockMat
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.
Array< OneD, NekDouble > m_tmpOUT
LibUtilities::NektarFFTSharedPtr m_FFT_z
LibUtilities::BasisSharedPtr m_homogeneousBasis_y
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
NekDouble m_lhom_z
Width of homogeneous direction z.
DNekMatSharedPtr MatBwdPAD
Array< OneD, NekDouble > m_tmpIN
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)
void HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool Shuff=true, bool UnShuff=true)
int m_ny
Number of modes = number of poitns in y direction.
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
LibUtilities::CommSharedPtr m_Ycomm
LibUtilities::CommSharedPtr m_Zcomm
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void Homogeneous2DTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool IsForwards, bool Shuff=true, bool UnShuff=true)
virtual void v_WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var)
virtual void v_DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray)
bool m_useFFT
FFT variables.
virtual void v_HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool Shuff=true, bool UnShuff=true)
LibUtilities::TranspositionSharedPtr m_transposition
Array< OneD, ExpListSharedPtr > m_lines
Vector of ExpList, will be filled with ExpList1D.
LibUtilities::BasisSharedPtr m_homogeneousBasis_z
Base expansion in z direction.
LibUtilities::BasisSharedPtr m_paddingBasis_z
Base expansion in z direction.
void SetPaddingBase(void)
DNekBlkMatSharedPtr GetHomogeneous2DBlockMatrix(Homogeneous2DMatType mattype) const
virtual void v_DealiasedDotProd(const Array< OneD, Array< OneD, NekDouble >> &inarray1, const Array< OneD, Array< OneD, NekDouble >> &inarray2, Array< OneD, Array< OneD, NekDouble >> &outarray)
LibUtilities::BasisSharedPtr m_paddingBasis_y
Base expansion in y direction.
NekDouble m_lhom_y
Width of homogeneous direction y.
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool Shuff=true, bool UnShuff=true)
virtual ~ExpListHomogeneous2D()
Destructor.
Base class for all multi-elemental spectral/hp expansions.
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
LibUtilities::CommSharedPtr m_comm
Communicator.
void DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray)
LibUtilities::SessionReaderSharedPtr m_session
Session.
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Class representing a segment element in reference space.
static BasisSharedPtr NullBasisSharedPtr
static std::vector< unsigned int > NullUnsignedIntVector
BasisManagerT & BasisManager(void)
std::shared_ptr< Basis > BasisSharedPtr
NektarFFTFactory & GetNektarFFTFactory()
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< FieldDefinitions > FieldDefinitionsSharedPtr
@ beta
Gauss Radau pinned at x=-1,.
@ eFourierEvenlySpaced
1D Evenly-spaced points using Fourier Fit
@ eFourier
Fourier Expansion .
std::map< Homogeneous2DMatType, DNekBlkMatSharedPtr > Homo2DBlockMatrixMap
A map between homo matrix keys and their associated block matrices.
@ eBackwardsCoeffSpaceY1D
@ eBackwardsCoeffSpaceZ1D
static const NekDouble kNekZeroTol
The above copyright notice and this permission notice shall be included.
std::shared_ptr< DNekBlkMat > DNekBlkMatSharedPtr
std::shared_ptr< DNekMat > DNekMatSharedPtr
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.
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.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Zero(int n, T *x, const int incx)
Zero vector.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)