44 namespace MultiRegions
47 ExpList3DHomogeneous1D::ExpList3DHomogeneous1D():
58 const bool dealiasing):
70 const bool dealiasing,
72 const std::string &var,
87 const bool dealiasing,
111 for(j = 0; j < nel; ++j)
116 for(n = 1; n <
m_planes.num_elements(); ++n)
119 for(j = 0; j < nel; ++j)
121 (*m_exp).push_back((*
m_exp)[j]);
142 if(DeclarePlanesSetCoeffPhys)
147 for(
int n = 0; n <
m_planes.num_elements(); ++n)
162 const std::vector<unsigned int> &eIDs,
163 bool DeclarePlanesSetCoeffPhys,
169 if(DeclarePlanesSetCoeffPhys)
172 std::vector<unsigned int> eIDsPlane;
173 int nel = eIDs.size()/
m_planes.num_elements();
175 for (
int i = 0; i < nel; ++i)
177 eIDsPlane.push_back(eIDs[i]);
186 for(
int n = 0; n <
m_planes.num_elements(); ++n)
205 int ncoeffs_per_plane =
m_planes[0]->GetNcoeffs();
206 int npoints_per_plane =
m_planes[0]->GetTotPoints();
208 int nzplanes =
m_planes.num_elements();
217 int nel =
m_planes[0]->GetExpSize();
222 for(cnt = n = 0; n < nzplanes; ++n)
225 m_planes[n]->SetPhysArray(tmparray =
m_phys + npoints_per_plane*n);
227 for(i = 0; i < nel; ++i)
230 m_phys_offset[cnt++] =
m_planes[n]->GetPhys_Offset(i) + n*npoints_per_plane;
242 int nzplanes =
m_planes.num_elements();
245 (*m_exp)[eid]->GetCoords(xc0,xc1);
251 for(n = 0; n <
m_planes.num_elements(); n++)
261 for(n = 0; n < nzplanes; ++n)
263 Vmath::Fill(npoints,z[n],tmp_xc = xc2 + npoints*n,1);
294 int nzplanes =
m_planes.num_elements();
295 int npoints =
m_planes[0]->GetTotPoints();
304 for(n = 0; n <
m_planes.num_elements(); n++)
314 for(n = 0; n < nzplanes; ++n)
316 Vmath::Fill(npoints,z[n],tmp_xc = xc2 + npoints*n,1);
327 ASSERTL0(expansion == -1,
"Multi-zone output not supported for homogeneous expansions.");
329 const int nPtsPlane =
m_planes[0]->GetNpoints();
330 const int nElmt =
m_planes[0]->GetExpSize();
331 const int nPlanes =
m_planes.num_elements();
335 for (
int i = 0; i < nElmt; ++i)
337 const int np0 = (*m_exp)[i]->GetNumPoints(0);
338 const int np1 = (*m_exp)[i]->GetNumPoints(1);
340 for (
int n = 1; n < nPlanes; ++n)
342 const int o1 = (n-1) * nPtsPlane;
343 const int o2 = n * nPtsPlane;
344 for (
int j = 1; j < np1; ++j)
346 for(
int k = 1; k < np0; ++k)
348 outfile << cnt + (j-1)*np0 + (k-1) + o1 + 1 <<
" ";
349 outfile << cnt + (j-1)*np0 + (k-1) + o2 + 1 <<
" ";
350 outfile << cnt + (j-1)*np0 + k + o2 + 1 <<
" ";
351 outfile << cnt + (j-1)*np0 + k + o1 + 1 <<
" ";
352 outfile << cnt + j *np0 + (k-1) + o1 + 1 <<
" ";
353 outfile << cnt + j *np0 + (k-1) + o2 + 1 <<
" ";
354 outfile << cnt + j *np0 + k + o2 + 1 <<
" ";
355 outfile << cnt + j *np0 + k + o1 + 1 << endl;
367 std::ostream &outfile,
374 m_planes[0]->WriteVtkPieceHeader(outfile, expansion);
379 int outputExtraPlane = 0;
384 outputExtraPlane = 1;
388 int nq0 = (*m_exp)[expansion]->GetNumPoints(0);
389 int nq1 = (*m_exp)[expansion]->GetNumPoints(1);
390 int nq2 =
m_planes.num_elements() + outputExtraPlane;
391 int ntot = nq0*nq1*nq2;
392 int ntotminus = (nq0-1)*(nq1-1)*(nq2-1);
398 GetCoords(expansion,coords[0],coords[1],coords[2]);
400 if (outputExtraPlane)
405 tmp = coords[0] + (nq2-1)*nq0*nq1, 1);
407 tmp = coords[1] + (nq2-1)*nq0*nq1, 1);
410 (coords[2][nq0*nq1] - coords[2][0]);
411 Vmath::Fill(nq0*nq1, z, tmp = coords[2] + (nq2-1)*nq0*nq1, 1);
415 m_session->LoadParameter(
"DistStrip", DistStrip, 0);
417 for(
int i = 0; i < ntot; i++)
419 coords[2][i] += istrip*DistStrip;
422 outfile <<
" <Piece NumberOfPoints=\""
423 << ntot <<
"\" NumberOfCells=\""
424 << ntotminus <<
"\">" << endl;
425 outfile <<
" <Points>" << endl;
426 outfile <<
" <DataArray type=\"Float64\" "
427 <<
"NumberOfComponents=\"3\" format=\"ascii\">" << endl;
429 for (i = 0; i < ntot; ++i)
431 for (j = 0; j < 3; ++j)
433 outfile << coords[j][i] <<
" ";
438 outfile <<
" </DataArray>" << endl;
439 outfile <<
" </Points>" << endl;
440 outfile <<
" <Cells>" << endl;
441 outfile <<
" <DataArray type=\"Int32\" "
442 <<
"Name=\"connectivity\" format=\"ascii\">" << endl;
443 for (i = 0; i < nq0-1; ++i)
445 for (j = 0; j < nq1-1; ++j)
447 for (k = 0; k < nq2-1; ++k)
449 outfile << k*nq0*nq1 + j*nq0 + i <<
" "
450 << k*nq0*nq1 + j*nq0 + i + 1 <<
" "
451 << k*nq0*nq1 + (j+1)*nq0 + i + 1 <<
" "
452 << k*nq0*nq1 + (j+1)*nq0 + i <<
" "
453 << (k+1)*nq0*nq1 + j*nq0 + i <<
" "
454 << (k+1)*nq0*nq1 + j*nq0 + i + 1 <<
" "
455 << (k+1)*nq0*nq1 + (j+1)*nq0 + i + 1 <<
" "
456 << (k+1)*nq0*nq1 + (j+1)*nq0 + i << endl;
461 outfile <<
" </DataArray>" << endl;
462 outfile <<
" <DataArray type=\"Int32\" "
463 <<
"Name=\"offsets\" format=\"ascii\">" << endl;
464 for (i = 0; i < ntotminus; ++i)
466 outfile << i*8+8 <<
" ";
469 outfile <<
" </DataArray>" << endl;
470 outfile <<
" <DataArray type=\"UInt8\" "
471 <<
"Name=\"types\" format=\"ascii\">" << endl;
472 for (i = 0; i < ntotminus; ++i)
477 outfile <<
" </DataArray>" << endl;
478 outfile <<
" </Cells>" << endl;
479 outfile <<
" <PointData>" << endl;
491 for(
int n = 0; n <
m_planes.num_elements(); n++)
498 for(
int n = 0; n <
m_planes.num_elements(); ++n)
500 errL2 =
m_planes[n]->L2(inarray + cnt);
502 err += errL2*errL2*local_w[n]*
m_lhom*0.5;
507 for(
int n = 0; n <
m_planes.num_elements(); ++n)
509 errL2 =
m_planes[n]->L2(inarray + cnt, soln + cnt);
511 err += errL2*errL2*local_w[n]*
m_lhom*0.5;
525 for (
int n = 0; n <
m_planes[0]->GetExpSize(); ++n)
528 area +=
m_planes[0]->GetExp(n)->Integral(inarray);
534 for (
int n = 0; n <
m_planes.num_elements(); n += 2)
540 for(
int i = 0; i <
m_planes[n]->GetExpSize(); ++i)
547 energy[n/2] += err*err;
553 energy[n/2] += err*err;
557 energy[n/2] /= 2.0*area;
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
const Array< OneD, const NekDouble > & GetCoeffs() const
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
virtual void v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
#define ASSERTL0(condition, msg)
int GetCoeff_Offset(int n) const
Get the start offset position for a global list of m_coeffs correspoinding to element n...
virtual Array< OneD, const NekDouble > v_HomogeneousEnergy(void)
static Array< OneD, NekDouble > NullNekDouble1DArray
virtual ~ExpList3DHomogeneous1D()
Destructor.
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
LibUtilities::TranspositionSharedPtr m_transposition
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
NekDouble m_lhom
Width of homogeneous direction.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
int GetExpSize(void)
This function returns the number of elements in the expansion.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
virtual void v_GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2)
void SetCoeffPhys(void)
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
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...
1D Evenly-spaced points using Fourier Fit
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
ExpList3DHomogeneous1D()
Default constructor.
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
int GetTotPoints(void) const
Returns the total number of quadrature points m_npoints .
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
LibUtilities::SessionReaderSharedPtr m_session
Session.
Array< OneD, ExpListSharedPtr > m_planes
boost::shared_ptr< ExpList2D > ExpList2DSharedPtr
Shared pointer to an ExpList2D object.
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
LibUtilities::CommSharedPtr m_comm
Communicator.
void GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1=NullNekDouble1DArray, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
This function calculates the coordinates of all the elemental quadrature points . ...
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
void GenExpList3DHomogeneous1D(const SpatialDomains::ExpansionMap &expansions, const Collections::ImplementationType ImpType)
virtual NekDouble v_L2(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
Describes the specification for a Basis.
std::map< int, ExpansionShPtr > ExpansionMap