Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
SplitFld.cpp File Reference
#include <cstdio>
#include <cstdlib>
#include <LibUtilities/BasicUtils/SessionReader.h>
#include <MultiRegions/ExpList.h>
#include <MultiRegions/ContField1D.h>
#include <MultiRegions/ContField2D.h>
#include <MultiRegions/ContField3D.h>
#include <MultiRegions/ContField3DHomogeneous1D.h>
#include <MultiRegions/ContField3DHomogeneous2D.h>
Include dependency graph for SplitFld.cpp:

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
void SetFields (SpatialDomains::MeshGraphSharedPtr &mesh, SpatialDomains::BoundaryConditionsSharedPtr &boundaryConditions, LibUtilities::SessionReaderSharedPtr &session, Array< OneD, MultiRegions::ExpListSharedPtr > &Exp, int nvariables)

Function Documentation

int main ( int  argc,
char *  argv[] 
)

Definition at line 14 of file SplitFld.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::SessionReader::CreateInstance(), Nektar::LibUtilities::ePolyEvenlySpaced, GetReflectionIndex(), Nektar::LibUtilities::Import(), Nektar::SpatialDomains::MeshGraph::Read(), SetFields(), and Nektar::LibUtilities::Write().

{
Array<OneD,MultiRegions::ExpListSharedPtr> &Exp,int nvariables);
Array<OneD, int> GetReflectionIndex(Array<OneD,MultiRegions::ExpListSharedPtr> &Exp,
int Ireg);
if(argc != 3)
{
fprintf(stderr,"Usage: SplitFld meshfile fieldfile\n");
exit(1);
}
= LibUtilities::SessionReader::CreateInstance(argc, argv);
//----------------------------------------------
// Read in mesh from input file
string meshfile(argv[argc-2]);
SpatialDomains::MeshGraphSharedPtr graphShPt = SpatialDomains::MeshGraph::Read(vSession);//meshfile);
//----------------------------------------------
// Also read and store the boundary conditions
::AllocateSharedPtr(vSession, graphShPt);
SpatialDomains::BoundaryConditions bcs(vSession, graphShPt);
//----------------------------------------------
//----------------------------------------------
// Import field file.
string fieldfile(argv[argc-1]);
vector<LibUtilities::FieldDefinitionsSharedPtr> fielddef;
vector<vector<NekDouble> > fielddata;
LibUtilities::Import(fieldfile,fielddef,fielddata);
//----------------------------------------------
// Define Expansion
int nfields;
nfields = fielddef[0]->m_fields.size();
Array<OneD, MultiRegions::ExpListSharedPtr> Exp;
Exp = Array<OneD, MultiRegions::ExpListSharedPtr>(nfields);
std::string solvtype = vSession->GetSolverInfo("SOLVERTYPE");
if(solvtype == "CoupledLinearisedNS" && vSession->DefinesSolverInfo("HOMOGENEOUS") )
{
SetFields(graphShPt,boundaryConditions,vSession,Exp,nfields-1);
//decomment
//nfields = nfields-1;
//start
int lastfield = nfields-1;
cout<<"Set pressure: "<<lastfield<<endl;
int nplanes = fielddef[0]->m_numModes[2];
const LibUtilities::BasisKey Bkey(fielddef[0]->m_basis[2],nplanes,Pkey);
NekDouble lz = fielddef[0]->m_homogeneousLengths[0];
Exp3DH1 = MemoryManager<MultiRegions::ExpList3DHomogeneous1D>::AllocateSharedPtr(vSession,Bkey,lz,false,false,graphShPt,fielddef[0]->m_fields[0]);
Exp[lastfield] = Exp3DH1;
//end
}
else
{
SetFields(graphShPt,boundaryConditions,vSession,Exp,nfields);
}
//----------------------------------------------
//----------------------------------------------
// Copy data from field file
for(int j = 0; j < nfields; ++j)
{
for(int i = 0; i < fielddef.size(); ++i)
{
Exp[j]->ExtractDataToCoeffs(fielddef [i],
fielddata[i],
fielddef [i]->m_fields[j],
Exp[j]->UpdateCoeffs());
}
Exp[j]->BwdTrans_IterPerExp(Exp[j]->GetCoeffs(),Exp[j]->UpdatePhys());
}
//----------------------------------------------
// Write solution to file with additional computed fields
string fldfilename(argv[2]);
string out = fldfilename.substr(0, fldfilename.find_last_of("."));
string endfile("split.fld");
//Array<OneD, Array<OneD, NekDouble> > fieldcoeffs(1);
string outfile;
string var;
Array<OneD, Array<OneD, NekDouble> > fieldcoeffs(Exp.num_elements());
//NB in case of homo fields you CANNOT use the BwdTrans
//because the Im comp is set to 0
for(int i = 0; i < Exp.num_elements(); ++i)
{
fieldcoeffs[i] = Exp[i]->UpdateCoeffs();
}
// copy Data into FieldData and set variable
/*
for(int g=0; g<Exp[0]->GetPlane(1)->GetNcoeffs(); g++)
{
cout<<"g="<<g<<" coeff f0="<<Exp[lastfield]->GetPlane(0)->GetCoeff(g)<<" f1="<<Exp[lastfield]->GetPlane(1)->GetCoeff(g)<<endl;
}
*/
for(int j =0; j<nfields; j++)
{
outfile = out;
std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
= Exp[j]->GetFieldDefinitions();
std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
//fieldcoeffs[j] = fields[j]->UpdateCoeffs();
for(int i = 0; i < FieldDef.size(); ++i)
{
var = fielddef[i]->m_fields[j];
// Could do a search here to find correct variable
FieldDef[i]->m_fields.push_back(var);
Exp[j]->AppendFieldData(FieldDef[i], FieldData[i],fieldcoeffs[j]);
}
outfile += "_"+var+"_"+endfile;
LibUtilities::Write(outfile,FieldDef,FieldData);
}
//-----------------------------------------------
return 0;
}
void SetFields ( SpatialDomains::MeshGraphSharedPtr mesh,
SpatialDomains::BoundaryConditionsSharedPtr boundaryConditions,
LibUtilities::SessionReaderSharedPtr session,
Array< OneD, MultiRegions::ExpListSharedPtr > &  Exp,
int  nvariables 
)

< physical length in Y direction (if homogeneous)

< physical length in Z direction (if homogeneous)

< number of points in Y direction (if homogeneous)

< number of points in Z direction (if homogeneous)

Parameter for homogeneous expansions

Definition at line 163 of file SplitFld.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::LibUtilities::eFourier, and Nektar::LibUtilities::eFourierEvenlySpaced.

{
// Setting parameteres for homogenous problems
//NekDouble LhomX; ///< physical length in X direction (if homogeneous)
NekDouble LhomY; ///< physical length in Y direction (if homogeneous)
NekDouble LhomZ; ///< physical length in Z direction (if homogeneous)
bool DeclareCoeffPhysArrays = true;
//int npointsX; ///< number of points in X direction (if homogeneous)
int npointsY; ///< number of points in Y direction (if homogeneous)
int npointsZ; ///< number of points in Z direction (if homogeneous)
bool useFFT = false;
bool deal = false;
///Parameter for homogeneous expansions
enum HomogeneousType
{
eHomogeneous1D,
eHomogeneous2D,
eHomogeneous3D,
eNotHomogeneous
};
enum HomogeneousType HomogeneousType = eNotHomogeneous;
if(session->DefinesSolverInfo("HOMOGENEOUS"))
{
std::string HomoStr = session->GetSolverInfo("HOMOGENEOUS");
//m_spacedim = 3;
if((HomoStr == "HOMOGENEOUS1D")||(HomoStr == "Homogeneous1D")||
(HomoStr == "1D")||(HomoStr == "Homo1D"))
{
HomogeneousType = eHomogeneous1D;
npointsZ = session->GetParameter("HomModesZ");
LhomZ = session->GetParameter("LZ");
}
if((HomoStr == "HOMOGENEOUS2D")||(HomoStr == "Homogeneous2D")||
(HomoStr == "2D")||(HomoStr == "Homo2D"))
{
HomogeneousType = eHomogeneous2D;
npointsY = session->GetParameter("HomModesY");
LhomY = session->GetParameter("LY");
npointsZ = session->GetParameter("HomModesZ");
LhomZ = session->GetParameter("LZ");
}
if((HomoStr == "HOMOGENEOUS3D")||(HomoStr == "Homogeneous3D")||
(HomoStr == "3D")||(HomoStr == "Homo3D"))
{
HomogeneousType = eHomogeneous3D;
//npointsX = session->GetParameter("HomModesX");
//LhomX = session->GetParameter("LX");
npointsY = session->GetParameter("HomModesY");
LhomY = session->GetParameter("LY");
npointsZ = session->GetParameter("HomModesZ");
LhomZ = session->GetParameter("LZ");
}
if(session->DefinesSolverInfo("USEFFT"))
{
useFFT = true;
}
}
int i;
int expdim = mesh->GetMeshDimension();
//Exp= Array<OneD, MultiRegions::ExpListSharedPtr>(nvariables);
// I can always have 3 variables in a 2D mesh (oech vel component i a function which can depend on 1-3 var)
// Continuous Galerkin projection
switch(expdim)
{
case 1:
{
if(HomogeneousType == eHomogeneous2D)
{
const LibUtilities::BasisKey BkeyY(LibUtilities::eFourier,npointsY,PkeyY);
const LibUtilities::BasisKey BkeyZ(LibUtilities::eFourier,npointsZ,PkeyZ);
for(i = 0 ; i < nvariables; i++)
{
::AllocateSharedPtr(session,BkeyY,BkeyZ,LhomY,LhomZ,useFFT,deal,mesh,session->GetVariable(i));
}
}
else
{
for(i = 0 ; i < nvariables; i++)
{
::AllocateSharedPtr(session,mesh,session->GetVariable(i));
}
}
break;
}
case 2:
{
if(HomogeneousType == eHomogeneous1D)
{
const LibUtilities::BasisKey BkeyZ(LibUtilities::eFourier,npointsZ,PkeyZ);
for(i = 0 ; i < nvariables; i++)
{
::AllocateSharedPtr(session,BkeyZ,LhomZ,useFFT,deal,mesh,session->GetVariable(i));
}
}
else
{
i = 0;
::AllocateSharedPtr(session,mesh,session->GetVariable(i),DeclareCoeffPhysArrays);
Exp[0] = firstfield;
for(i = 1 ; i < nvariables; i++)
{
::AllocateSharedPtr(*firstfield,mesh,session->GetVariable(i),DeclareCoeffPhysArrays);
}
}
break;
}
case 3:
{
if(HomogeneousType == eHomogeneous3D)
{
ASSERTL0(false,"3D fully periodic problems not implemented yet");
}
else
{
i = 0;
::AllocateSharedPtr(session,mesh,session->GetVariable(i));
Exp[0] = firstfield;
for(i = 1 ; i < nvariables; i++)
{
::AllocateSharedPtr(*firstfield,mesh,session->GetVariable(i));
}
}
break;
}
default:
ASSERTL0(false,"Expansion dimension not recognised");
break;
}
}