Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Types | Private Member Functions | Private Attributes | Friends | List of all members
Nektar::AdjointAdvection Class Reference

#include <AdjointAdvection.h>

Inheritance diagram for Nektar::AdjointAdvection:
Inheritance graph
[legend]
Collaboration diagram for Nektar::AdjointAdvection:
Collaboration graph
[legend]

Static Public Member Functions

static AdvectionTermSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Creates an instance of this class.

Static Public Attributes

static std::string className = GetAdvectionTermFactory().RegisterCreatorFunction("Adjoint", AdjointAdvection::create)
 Name of class.

Protected Member Functions

DNekBlkMatSharedPtr GetFloquetBlockMatrix (FloquetMatType mattype, bool UseContCoeffs=false) const
DNekBlkMatSharedPtr GenFloquetBlockMatrix (FloquetMatType mattype, bool UseContCoeffs=false) const
 AdjointAdvection (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
virtual ~AdjointAdvection ()
virtual void v_InitObject ()
void SetUpBaseFields (SpatialDomains::MeshGraphSharedPtr &mesh)
void UpdateBase (const NekDouble m_slices, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble m_time, const NekDouble m_period)
void DFT (const string file, const NekDouble m_slices)
void ImportFldBase (std::string pInfile, SpatialDomains::MeshGraphSharedPtr pGraph, int cnt)
 Import Base flow.
void ImportFldBase (std::string pInfile, SpatialDomains::MeshGraphSharedPtr pGraph)
void WriteFldBase (std::string &outname)
 Write field data to the given filename.
void WriteFldBase (std::string &outname, MultiRegions::ExpListSharedPtr &field, Array< OneD, Array< OneD, NekDouble > > &fieldcoeffs, Array< OneD, std::string > &variables)
 Write input fields to the given filename.
- Protected Member Functions inherited from Nektar::AdvectionTerm
 AdvectionTerm (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Constructor.
virtual void v_ComputeAdvectionTerm (Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pV, const Array< OneD, const NekDouble > &pU, Array< OneD, NekDouble > &pOutarray, int pVelocityComponent, NekDouble m_time, Array< OneD, NekDouble > &pWk)

Protected Attributes

Array< OneD,
MultiRegions::ExpListSharedPtr
m_base
int m_slices
NekDouble m_period
Array< OneD, Array< OneD,
NekDouble > > 
m_interp
LibUtilities::NektarFFTSharedPtr m_FFT
Array< OneD, NekDoublem_tmpIN
Array< OneD, NekDoublem_tmpOUT
bool m_useFFTW
bool m_SingleMode
 flag to determine if use single mode or not
bool m_HalfMode
 flag to determine if use half mode or not
bool m_MultipleModes
 flag to determine if use multiple mode or not
FloquetBlockMatrixMapShPtr m_FloquetBlockMat
- Protected Attributes inherited from Nektar::AdvectionTerm
LibUtilities::SessionReaderSharedPtr m_session
std::string m_sessionName
 Name of the session.
SpatialDomains::MeshGraphSharedPtr m_graph
 Pointer to mesh graph.
bool m_homogen_dealiasing
 flag to determine if use Fourier dealising or not
bool m_specHP_dealiasing
 flag to determine if use Spectral/hp element dealising or not
bool m_SingleMode
 Flag to determine if use single mode or not.
bool m_HalfMode
 Flag to determine if use half mode or not.
MultiRegions::CoeffState m_CoeffState
enum MultiRegions::ProjectionType m_projectionType
 Type of projection, i.e. Galerkin or DG.
int m_spacedim
 Spatial dimension (> expansion dim)
int m_expdim
 Dimension of the expansion.
int nvariables
 Number of variables.
int m_nConvectiveFields
int m_slices
 Number of fields to be convected;.
NekDouble m_period
Array< OneD, Array< OneD,
NekDouble > > 
m_interp
LibUtilities::NektarFFTSharedPtr m_FFT
Array< OneD, NekDoublem_tmpIN
Array< OneD, NekDoublem_tmpOUT
bool m_useFFTW

Private Types

enum  FloquetMatType { eForwardsCoeff, eForwardsPhys }
enum  HomogeneousType { eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, eNotHomogeneous }
 Parameter for homogeneous expansions. More...
typedef map< FloquetMatType,
DNekBlkMatSharedPtr
FloquetBlockMatrixMap
 A map between matrix keys and their associated block matrices.
typedef boost::shared_ptr
< FloquetBlockMatrixMap
FloquetBlockMatrixMapShPtr
 A shared pointer to a BlockMatrixMap.

Private Member Functions

virtual void v_ComputeAdvectionTerm (Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pV, const Array< OneD, const NekDouble > &pU, Array< OneD, NekDouble > &pOutarray, int pVelocityComponent, NekDouble m_time, Array< OneD, NekDouble > &pWk)

Private Attributes

bool m_useFFT
 flag to determine if use or not the FFT for transformations
enum HomogeneousType m_HomogeneousType
NekDouble m_LhomX
 physical length in X direction (if homogeneous)
NekDouble m_LhomY
 physical length in Y direction (if homogeneous)
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous)
int m_npointsX
 number of points in X direction (if homogeneous)
int m_npointsY
 number of points in Y direction (if homogeneous)
int m_npointsZ
 number of points in Z direction (if homogeneous)
int m_HomoDirec
 number of homogenous directions
int m_NumMode
 Mode to use in case of single mode analysis.
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions

Friends

class MemoryManager< AdjointAdvection >

Additional Inherited Members

- Public Member Functions inherited from Nektar::AdvectionTerm
virtual ~AdvectionTerm ()
 Destructor.
void InitObject ()
void DoAdvection (Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int nConvectiveFields, const Array< OneD, int > &vel_loc, const Array< OneD, const Array< OneD, NekDouble > > &pInarray, Array< OneD, Array< OneD, NekDouble > > &pOutarray, NekDouble m_time, Array< OneD, NekDouble > &pWk=NullNekDouble1DArray)
 Compute advection term.
void DoAdvection (Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, const Array< OneD, NekDouble > > &Velocity, const Array< OneD, const Array< OneD, NekDouble > > &pInarray, Array< OneD, Array< OneD, NekDouble > > &pOutarray, NekDouble m_time, Array< OneD, NekDouble > &pWk=NullNekDouble1DArray)
bool GetSpecHPDealiasing (void)
void SetSpecHPDealiasing (bool value)

Detailed Description

Definition at line 53 of file AdjointAdvection.h.

Member Typedef Documentation

A map between matrix keys and their associated block matrices.

Definition at line 63 of file AdjointAdvection.h.

A shared pointer to a BlockMatrixMap.

Definition at line 65 of file AdjointAdvection.h.

Member Enumeration Documentation

Enumerator:
eForwardsCoeff 
eForwardsPhys 

Definition at line 55 of file AdjointAdvection.h.

Parameter for homogeneous expansions.

Enumerator:
eHomogeneous1D 
eHomogeneous2D 
eHomogeneous3D 
eNotHomogeneous 

Definition at line 151 of file AdjointAdvection.h.

Constructor & Destructor Documentation

Nektar::AdjointAdvection::AdjointAdvection ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)
protected

Constructor. Creates ...

Parameters
\param

Definition at line 61 of file AdjointAdvection.cpp.

:
AdvectionTerm(pSession, pGraph)
{
}
Nektar::AdjointAdvection::~AdjointAdvection ( )
protectedvirtual

Definition at line 252 of file AdjointAdvection.cpp.

{
}

Member Function Documentation

static AdvectionTermSharedPtr Nektar::AdjointAdvection::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)
inlinestatic

Creates an instance of this class.

Definition at line 71 of file AdjointAdvection.h.

{
p->InitObject();
return p;
}
void Nektar::AdjointAdvection::DFT ( const string  file,
const NekDouble  m_slices 
)
protected

Definition at line 838 of file AdjointAdvection.cpp.

References Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), eForwardsPhys, Nektar::eWrapper, GetFloquetBlockMatrix(), Nektar::LibUtilities::GetNektarFFTFactory(), ImportFldBase(), m_base, m_FFT, Nektar::AdvectionTerm::m_graph, m_interp, m_period, Nektar::AdvectionTerm::m_session, m_slices, m_tmpIN, m_tmpOUT, Vmath::Smul(), Vmath::Vcopy(), and Vmath::Zero().

Referenced by v_InitObject().

{
int npoints=m_base[0]->GetTotPoints();
//Convected fields
int ConvectedFields=m_base.num_elements()-1;
m_interp= Array<OneD, Array<OneD, NekDouble> > (ConvectedFields);
for(int i=0; i<ConvectedFields;++i)
{
m_interp[i]=Array<OneD,NekDouble>(npoints*m_slices);
}
//Import the slides into auxiliary vector
//The base flow should be stored in the form filename_i.bse
for (int i=0; i< m_slices; ++i)
{
char chkout[16] = "";
sprintf(chkout, "%d", i);
ImportFldBase(file+"_"+chkout+".bse",m_graph,i);
}
// Discrete Fourier Transform of the fields
for(int k=0; k< ConvectedFields;++k)
{
#ifdef NEKTAR_USING_FFTW
//Discrete Fourier Transform using FFTW
Array<OneD, NekDouble> fft_in(npoints*m_slices);
Array<OneD, NekDouble> fft_out(npoints*m_slices);
Array<OneD, NekDouble> m_tmpIN(m_slices);
Array<OneD, NekDouble> m_tmpOUT(m_slices);
//Shuffle the data
for(int j= 0; j < m_slices; ++j)
{
Vmath::Vcopy(npoints,&m_interp[k][j*npoints],1,&(fft_in[j]),m_slices);
}
//FFT Transform
for(int i=0; i<npoints; i++)
{
m_FFT->FFTFwdTrans(m_tmpIN =fft_in + i*m_slices, m_tmpOUT =fft_out + i*m_slices);
}
//Reshuffle data
for(int s = 0; s < m_slices; ++s)
{
Vmath::Vcopy(npoints,&fft_out[s],m_slices,&m_interp[k][s*npoints],1);
}
Vmath::Zero(fft_in.num_elements(),&fft_in[0],1);
Vmath::Zero(fft_out.num_elements(),&fft_out[0],1);
#else
//Discrete Fourier Transform using MVM
int nrows = blkmat->GetRows();
int ncols = blkmat->GetColumns();
Array<OneD, NekDouble> sortedinarray(ncols);
Array<OneD, NekDouble> sortedoutarray(nrows);
//Shuffle the data
for(int j= 0; j < m_slices; ++j)
{
Vmath::Vcopy(npoints,&m_interp[k][j*npoints],1,&(sortedinarray[j]),m_slices);
}
// Create NekVectors from the given data arrays
NekVector<NekDouble> in (ncols,sortedinarray,eWrapper);
NekVector<NekDouble> out(nrows,sortedoutarray,eWrapper);
// Perform matrix-vector multiply.
out = (*blkmat)*in;
//Reshuffle data
for(int s = 0; s < m_slices; ++s)
{
Vmath::Vcopy(npoints,&sortedoutarray[s],m_slices,&m_interp[k][s*npoints],1);
}
Vmath::Zero(sortedinarray.num_elements(),&sortedinarray[0],1);
Vmath::Zero(sortedoutarray.num_elements(),&sortedoutarray[0],1);
#endif
//scaling of the Fourier coefficients
NekDouble j=-1;
for (int i = 2; i < m_slices; i += 2)
{
Vmath::Smul(2*npoints,j,&m_interp[k][i*npoints],1,&m_interp[k][i*npoints],1);
j=-j;
}
}
if(m_session->DefinesParameter("period"))
{
m_period=m_session->GetParameter("period");
}
else
{
m_period=(m_session->GetParameter("TimeStep")*m_slices)/(m_slices-1.);
}
}
DNekBlkMatSharedPtr Nektar::AdjointAdvection::GenFloquetBlockMatrix ( FloquetMatType  mattype,
bool  UseContCoeffs = false 
) const
protected
DNekBlkMatSharedPtr Nektar::AdjointAdvection::GetFloquetBlockMatrix ( FloquetMatType  mattype,
bool  UseContCoeffs = false 
) const
protected

Definition at line 799 of file AdjointAdvection.cpp.

References Nektar::StdRegions::StdExpansion::DetShapeType(), Nektar::eDIAGONAL, Nektar::LibUtilities::eFourier, Nektar::LibUtilities::eFourierEvenlySpaced, Nektar::StdRegions::eFwdTrans, Nektar::StdRegions::StdExpansion::GetStdMatrix(), m_base, and m_slices.

Referenced by DFT().

{
int n_exp = 0;
n_exp = m_base[0]->GetTotPoints(); // will operatore on m_phys
Array<OneD,unsigned int> nrows(n_exp);
Array<OneD,unsigned int> ncols(n_exp);
nrows = Array<OneD, unsigned int>(n_exp,m_slices);
ncols = Array<OneD, unsigned int>(n_exp,m_slices);
MatrixStorage blkmatStorage = eDIAGONAL;
::AllocateSharedPtr(nrows,ncols,blkmatStorage);
StdSeg.DetShapeType(),
StdSeg);
loc_mat = StdSeg.GetStdMatrix(matkey);
// set up array of block matrices.
for(int i = 0; i < n_exp; ++i)
{
BlkMatrix->SetBlock(i,i,loc_mat);
}
return BlkMatrix;
}
void Nektar::AdjointAdvection::ImportFldBase ( std::string  pInfile,
SpatialDomains::MeshGraphSharedPtr  pGraph,
int  cnt 
)
protected

Import Base flow.

Import field from infile and load into m_fields. This routine will also perform a BwdTrans to ensure data is in both the physical and coefficient storage.

Parameters
infileFilename to read.

Definition at line 459 of file AdjointAdvection.cpp.

References ASSERTL1, m_base, m_HalfMode, m_interp, Nektar::AdvectionTerm::m_nConvectiveFields, m_npointsZ, Nektar::AdvectionTerm::m_session, m_SingleMode, Vmath::Vcopy(), and Vmath::Zero().

Referenced by DFT(), and v_InitObject().

{
std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
int nqtot = m_base[0]->GetTotPoints();
//Get Homogeneous
m_session->GetComm());
fld->Import(pInfile, FieldDef, FieldData);
int nvar = m_session->GetVariables().size();
int s;
if(m_session->DefinesSolverInfo("HOMOGENEOUS"))
{
std::string HomoStr = m_session->GetSolverInfo("HOMOGENEOUS");
}
// copy FieldData into m_fields
for(int j = 0; j < nvar; ++j)
{
for(int i = 0; i < FieldDef.size(); ++i)
{
if ((m_session->DefinesSolverInfo("HOMOGENEOUS") &&
(m_session->GetSolverInfo("HOMOGENEOUS")=="HOMOGENEOUS1D" ||
m_session->GetSolverInfo("HOMOGENEOUS")=="1D" ||
m_session->GetSolverInfo("HOMOGENEOUS")=="Homo1D")) &&
nvar==4)
{
// w-component must be ignored and set to zero.
if (j != nvar - 2)
{
// p component (it is 4th variable of the 3D and corresponds 3nd variable of 2D)
s = (j == nvar - 1) ? 2 : j;
//extraction of the 2D
m_base[j]->ExtractDataToCoeffs(
FieldDef[i],
FieldData[i],
FieldDef[i]->m_fields[s],
m_base[j]->UpdateCoeffs());
}
// Put zero on higher modes
int ncplane = (m_base[0]->GetNcoeffs()) / m_npointsZ;
if (m_npointsZ > 2)
{
Vmath::Zero(ncplane*(m_npointsZ-2),
&m_base[j]->UpdateCoeffs()[2*ncplane],1);
}
}
//2D cases and Homogeneous1D Base Flows
else
{
bool flag = FieldDef[i]->m_fields[j] ==
m_session->GetVariable(j);
ASSERTL1(flag, (std::string("Order of ") + pInfile
+ std::string(" data and that defined in "
"m_boundaryconditions differs")).c_str());
m_base[j]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
FieldDef[i]->m_fields[j],
m_base[j]->UpdateCoeffs());
}
}
{
m_base[j]->SetWaveSpace(true);
m_base[j]->BwdTrans(m_base[j]->GetCoeffs(),
m_base[j]->UpdatePhys());
{
//copy the bwd into the second plane for single Mode Analysis
int ncplane=(m_base[0]->GetNpoints())/m_npointsZ;
Vmath::Vcopy(ncplane,&m_base[j]->GetPhys()[0],1,&m_base[j]->UpdatePhys()[ncplane],1);
}
}
else
{
m_base[j]->BwdTrans(m_base[j]->GetCoeffs(),
m_base[j]->UpdatePhys());
}
}
//std::string outname ="BaseFlow.bse";
//WriteFldBase(outname);
if(m_session->DefinesParameter("N_slices"))
{
m_nConvectiveFields = m_base.num_elements()-1;
for(int i=0; i<m_nConvectiveFields;++i)
{
Vmath::Vcopy(nqtot, &m_base[i]->GetPhys()[0], 1, &m_interp[i][cnt*nqtot], 1);
}
}
}
void Nektar::AdjointAdvection::ImportFldBase ( std::string  pInfile,
SpatialDomains::MeshGraphSharedPtr  pGraph 
)
protected
void Nektar::AdjointAdvection::SetUpBaseFields ( SpatialDomains::MeshGraphSharedPtr mesh)
protected

Definition at line 257 of file AdjointAdvection.cpp.

References ASSERTL0, Nektar::LibUtilities::eFourier, Nektar::LibUtilities::eFourierEvenlySpaced, Nektar::LibUtilities::eFourierHalfModeRe, Nektar::LibUtilities::eFourierSingleModeSpaced, Nektar::MultiRegions::eGalerkin, eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, m_base, Nektar::AdvectionTerm::m_expdim, Nektar::AdvectionTerm::m_graph, m_HalfMode, Nektar::AdvectionTerm::m_homogen_dealiasing, m_HomogeneousType, m_LhomY, m_LhomZ, m_npointsY, m_npointsZ, Nektar::AdvectionTerm::m_projectionType, Nektar::AdvectionTerm::m_session, m_SingleMode, m_useFFT, and Nektar::AdvectionTerm::nvariables.

Referenced by v_InitObject().

{
int nvariables = m_session->GetVariables().size();
int i;
m_base = Array<OneD, MultiRegions::ExpListSharedPtr>(nvariables);
{
switch (m_expdim)
{
case 1:
{
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
}
}
else {
for(i = 0 ; i < m_base.num_elements(); i++)
{
m_session->GetVariable(i));
}
}
}
break;
case 2:
{
{
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
}
}
else if(m_HalfMode)
{
//1 plane field (half mode expansion)
for(i = 0 ; i < m_base.num_elements(); i++)
{
}
}
else
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
m_base[i]->SetWaveSpace(false);
}
}
}
else
{
i = 0;
m_session->GetVariable(i));
m_base[0]=firstbase;
for(i = 1 ; i < m_base.num_elements(); i++)
{
::AllocateSharedPtr(*firstbase,mesh,
m_session->GetVariable(i));
}
}
}
break;
case 3:
{
{
ASSERTL0(false,"3D fully periodic problems not implemented yet");
}
else
{
m_session->GetVariable(0));
m_base[0] = firstbase;
for(i = 1 ; i < m_base.num_elements(); i++)
{
::AllocateSharedPtr(*firstbase,mesh,
m_session->GetVariable(i));
}
}
}
break;
default:
ASSERTL0(false,"Expansion dimension not recognised");
break;
}
}
else
{
switch(m_expdim)
{
case 1:
{
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
}
}
else
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
::DisContField1D>::AllocateSharedPtr(m_session,mesh,
m_session->GetVariable(i));
}
}
break;
}
case 2:
{
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
}
}
else
{
for(i = 0 ; i < m_base.num_elements(); i++)
{
::DisContField2D>::AllocateSharedPtr(m_session, mesh,
m_session->GetVariable(i));
}
}
break;
}
case 3:
ASSERTL0(false,"3 D not set up");
default:
ASSERTL0(false,"Expansion dimension not recognised");
break;
}
}
}
void Nektar::AdjointAdvection::UpdateBase ( const NekDouble  m_slices,
Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const NekDouble  m_time,
const NekDouble  m_period 
)
protected

Definition at line 773 of file AdjointAdvection.cpp.

References m_base, m_period, m_slices, Vmath::Svtvp(), and Vmath::Vcopy().

Referenced by v_ComputeAdvectionTerm().

{
int npoints=m_base[0]->GetTotPoints();
NekDouble BetaT=2*M_PI*fmod (m_time, m_period) / m_period;
NekDouble phase;
Array<OneD, NekDouble> auxiliary(npoints);
Vmath::Vcopy(npoints,&inarray[0],1,&outarray[0],1);
Vmath::Svtvp(npoints, cos(0.5*m_slices*BetaT),&inarray[npoints],1,&outarray[0],1,&outarray[0],1);
for (int i = 2; i < m_slices; i += 2)
{
phase = (i>>1) * BetaT;
Vmath::Svtvp(npoints, cos(phase),&inarray[i*npoints],1,&outarray[0],1,&outarray[0],1);
Vmath::Svtvp(npoints, sin(phase), &inarray[(i+1)*npoints], 1, &outarray[0], 1,&outarray[0],1);
}
}
void Nektar::AdjointAdvection::v_ComputeAdvectionTerm ( Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const Array< OneD, Array< OneD, NekDouble > > &  pV,
const Array< OneD, const NekDouble > &  pU,
Array< OneD, NekDouble > &  pOutarray,
int  pVelocityComponent,
NekDouble  m_time,
Array< OneD, NekDouble > &  pWk 
)
privatevirtual

Definition at line 571 of file AdjointAdvection.cpp.

References ASSERTL0, Nektar::LibUtilities::eFunctionTypeFile, m_base, m_HalfMode, m_interp, Nektar::AdvectionTerm::m_nConvectiveFields, m_period, Nektar::AdvectionTerm::m_session, m_slices, Vmath::Neg(), UpdateBase(), Vmath::Vmul(), and Vmath::Vvtvp().

{
int ndim = m_nConvectiveFields;
int nPointsTot = pFields[0]->GetNpoints();
Array<OneD, NekDouble> grad0,grad1,grad2;
//Evaluation of the gradiend of each component of the base flow
//\nabla U
Array<OneD, NekDouble> grad_base_u0,grad_base_u1,grad_base_u2;
// \nabla V
Array<OneD, NekDouble> grad_base_v0,grad_base_v1,grad_base_v2;
// \nabla W
Array<OneD, NekDouble> grad_base_w0,grad_base_w1,grad_base_w2;
grad0 = Array<OneD, NekDouble> (nPointsTot);
grad_base_u0 = Array<OneD, NekDouble> (nPointsTot);
grad_base_v0 = Array<OneD, NekDouble> (nPointsTot);
grad_base_w0 = Array<OneD, NekDouble> (nPointsTot);
//Evaluation of the base flow for periodic cases
//(it requires fld files)
if(m_slices>1)
{
if (m_session->GetFunctionType("BaseFlow", 0)
{
for(int i=0; i<m_nConvectiveFields;++i)
{
UpdateBase(m_slices,m_interp[i],m_base[i]->UpdatePhys(),m_time,m_period);
}
}
else
{
ASSERTL0(false, "Periodic Base flow requires .fld files");
}
}
//Evaluate the Adjoint advection term
switch(ndim)
{
// 1D
case 1:
pFields[0]->PhysDeriv(pU,grad0);
pFields[0]->PhysDeriv(m_base[0]->GetPhys(),grad_base_u0);
//Evaluate U du'/dx
Vmath::Vmul(nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U du'/dx+ u' dU/dx
Vmath::Vvtvp(nPointsTot,grad_base_u0,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
break;
//2D
case 2:
grad1 = Array<OneD, NekDouble> (nPointsTot);
grad_base_u1 = Array<OneD, NekDouble> (nPointsTot);
grad_base_v1 = Array<OneD, NekDouble> (nPointsTot);
pFields[0]->PhysDeriv(pU,grad0,grad1);
//Derivates of the base flow
pFields[0]-> PhysDeriv(m_base[0]->GetPhys(), grad_base_u0, grad_base_u1);
pFields[0]-> PhysDeriv(m_base[1]->GetPhys(), grad_base_v0, grad_base_v1);
//Since the components of the velocity are passed one by one, it is necessary to distinguish which
//term is consider
switch (pVelocityComponent)
{
//x-equation
case 0:
// Evaluate U du'/dx
Vmath::Vmul (nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U du'/dx+ V du'/dy
Vmath::Vvtvp(nPointsTot,grad1,1,m_base[1]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate - (U du'/dx+ V du'/dy)
Vmath::Neg(nPointsTot,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy)+u' dU/dx
Vmath::Vvtvp(nPointsTot,grad_base_u0,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy) +u' dU/dx +v' dV/dx
Vmath::Vvtvp(nPointsTot,grad_base_v0,1,pVelocity[1],1,pOutarray,1,pOutarray,1);
break;
//y-equation
case 1:
// Evaluate U dv'/dx
Vmath::Vmul (nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U dv'/dx+ V dv'/dy
Vmath::Vvtvp(nPointsTot,grad1,1,m_base[1]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate -(U dv'/dx+ V dv'/dy)
Vmath::Neg(nPointsTot,pOutarray,1);
//Evaluate (U dv'/dx+ V dv'/dy)+u' dU/dy
Vmath::Vvtvp(nPointsTot,grad_base_u1,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
//Evaluate (U dv'/dx+ V dv'/dy +u' dv/dx)+v' dV/dy
Vmath::Vvtvp(nPointsTot,grad_base_v1,1,pVelocity[1],1,pOutarray,1,pOutarray,1);
break;
}
break;
//3D
case 3:
grad1 = Array<OneD, NekDouble> (nPointsTot);
grad2 = Array<OneD, NekDouble> (nPointsTot);
grad_base_u1 = Array<OneD, NekDouble> (nPointsTot);
grad_base_v1 = Array<OneD, NekDouble> (nPointsTot);
grad_base_w1 = Array<OneD, NekDouble> (nPointsTot);
grad_base_u2 = Array<OneD, NekDouble> (nPointsTot);
grad_base_v2 = Array<OneD, NekDouble> (nPointsTot);
grad_base_w2 = Array<OneD, NekDouble> (nPointsTot);
m_base[0]->PhysDeriv(m_base[0]->GetPhys(), grad_base_u0, grad_base_u1,grad_base_u2);
m_base[0]->PhysDeriv(m_base[1]->GetPhys(), grad_base_v0, grad_base_v1,grad_base_v2);
m_base[0]->PhysDeriv(m_base[2]->GetPhys(), grad_base_w0, grad_base_w1, grad_base_w2);
//HalfMode has W(x,y,t)=0
{
for(int i=0; i<grad_base_u2.num_elements();++i)
{
grad_base_u2[i]=0;
grad_base_v2[i]=0;
grad_base_w2[i]=0;
}
}
pFields[0]->PhysDeriv(pU, grad0, grad1, grad2);
switch (pVelocityComponent)
{
//x-equation
case 0:
//Evaluate U du'/dx
Vmath::Vmul (nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U du'/dx+ V du'/dy
Vmath::Vvtvp(nPointsTot,grad1,1,m_base[1]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate U du'/dx+ V du'/dy+W du'/dz
Vmath::Vvtvp(nPointsTot,grad2,1,m_base[2]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy+W du'/dz)
Vmath::Neg(nPointsTot,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy+W du'/dz)+u' dU/dx
Vmath::Vvtvp(nPointsTot,grad_base_u0,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy+W du'/dz)+u'dU/dx+ v' dV/dx
Vmath::Vvtvp(nPointsTot,grad_base_v0,1,pVelocity[1],1,pOutarray,1,pOutarray,1);
//Evaluate -(U du'/dx+ V du'/dy+W du'/dz)+u'dU/dx+ v' dV/dx+ w' dW/dz
Vmath::Vvtvp(nPointsTot,grad_base_w0,1,pVelocity[2],1,pOutarray,1,pOutarray,1);
break;
//y-equation
case 1:
//Evaluate U dv'/dx
Vmath::Vmul (nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U dv'/dx+ V dv'/dy
Vmath::Vvtvp(nPointsTot,grad1,1,m_base[1]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate U dv'/dx+ V dv'/dy+W dv'/dz
Vmath::Vvtvp(nPointsTot,grad2,1,m_base[2]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate -(U dv'/dx+ V dv'/dy+W dv'/dz)
Vmath::Neg(nPointsTot,pOutarray,1);
//Evaluate -(U dv'/dx+ V dv'/dy+W dv'/dz)+u' dU/dy
Vmath::Vvtvp(nPointsTot,grad_base_u1,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
//Evaluate -(U dv'/dx+ V dv'/dy+W dv'/dz)+u' dU/dy +v' dV/dy
Vmath::Vvtvp(nPointsTot,grad_base_v1,1,pVelocity[1],1,pOutarray,1,pOutarray,1);
//Evaluate -(U dv'/dx+ V dv'/dy+W dv'/dz)+u' dU/dy +v' dV/dy+ w' dW/dy
Vmath::Vvtvp(nPointsTot,grad_base_w1,1,pVelocity[2],1,pOutarray,1,pOutarray,1);
break;
//z-equation
case 2:
//Evaluate U dw'/dx
Vmath::Vmul (nPointsTot,grad0,1,m_base[0]->GetPhys(),1,pOutarray,1);
//Evaluate U dw'/dx+ V dw'/dx
Vmath::Vvtvp(nPointsTot,grad1,1,m_base[1]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate U dw'/dx+ V dw'/dx+ W dw'/dz
Vmath::Vvtvp(nPointsTot,grad2,1,m_base[2]->GetPhys(),1,pOutarray,1,pOutarray,1);
//Evaluate -(U dw'/dx+ V dw'/dx+ W dw'/dz)
Vmath::Neg(nPointsTot,pOutarray,1);
//Evaluate -(U dw'/dx+ V dw'/dx+ W dw'/dz)+u' dU/dz
Vmath::Vvtvp(nPointsTot,grad_base_u2,1,pVelocity[0],1,pOutarray,1,pOutarray,1);
//Evaluate -(U dw'/dx+ V dw'/dx+ W dw'/dz)+u' dU/dz+v'dV/dz
Vmath::Vvtvp(nPointsTot,grad_base_v2,1,pVelocity[1],1,pOutarray,1,pOutarray,1);
//Evaluate -(U dw'/dx+ V dw'/dx+ W dw'/dz)+u' dU/dz+v'dV/dz + w' dW/dz
Vmath::Vvtvp(nPointsTot,grad_base_w2,1,pVelocity[2],1,pOutarray,1,pOutarray,1);
break;
}
break;
default:
ASSERTL0(false,"dimension unknown");
}
}
void Nektar::AdjointAdvection::v_InitObject ( )
protectedvirtual

Reimplemented from Nektar::AdvectionTerm.

Definition at line 69 of file AdjointAdvection.cpp.

References ASSERTL0, DFT(), Nektar::MultiRegions::eDiscontinuous, Nektar::LibUtilities::eFunctionTypeFile, Nektar::MultiRegions::eGalerkin, eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, eNotHomogeneous, ImportFldBase(), m_base, m_boundaryConditions, Nektar::AdvectionTerm::m_graph, m_HalfMode, m_HomoDirec, m_HomogeneousType, m_LhomX, m_LhomY, m_LhomZ, m_MultipleModes, m_npointsX, m_npointsY, m_npointsZ, Nektar::AdvectionTerm::m_projectionType, Nektar::AdvectionTerm::m_session, m_SingleMode, m_slices, Nektar::AdvectionTerm::m_spacedim, m_useFFT, and SetUpBaseFields().

{
//Setting parameters for homogeneous problems
m_useFFT = false;
m_SingleMode =false;
m_HalfMode =false;
if(m_session->DefinesSolverInfo("HOMOGENEOUS"))
{
std::string HomoStr = m_session->GetSolverInfo("HOMOGENEOUS");
if((HomoStr == "HOMOGENEOUS1D")||(HomoStr == "Homogeneous1D")||
(HomoStr == "1D")||(HomoStr == "Homo1D"))
{
m_LhomZ = m_session->GetParameter("LZ");
if(m_session->DefinesSolverInfo("ModeType"))
{
m_session->MatchSolverInfo("ModeType","SingleMode",m_SingleMode,false);
m_session->MatchSolverInfo("ModeType","HalfMode",m_HalfMode,false);
m_session->MatchSolverInfo("ModeType","MultipleModes",m_MultipleModes,false);
}
if(m_session->DefinesSolverInfo("ModeType"))
{
{
}
else if(m_HalfMode)
{
}
else if(m_MultipleModes)
{
m_npointsZ = m_session->GetParameter("HomModesZ");
}
else
{
ASSERTL0(false, "SolverInfo ModeType not valid");
}
}
else
{
m_session->LoadParameter("HomModesZ",m_npointsZ);
}
}
if((HomoStr == "HOMOGENEOUS2D")||(HomoStr == "Homogeneous2D")||
(HomoStr == "2D")||(HomoStr == "Homo2D"))
{
m_session->LoadParameter("HomModesY", m_npointsY);
m_session->LoadParameter("LY", m_LhomY);
m_session->LoadParameter("HomModesZ", m_npointsZ);
m_session->LoadParameter("LZ", m_LhomZ);
}
if((HomoStr == "HOMOGENEOUS3D")||(HomoStr == "Homogeneous3D")||
(HomoStr == "3D")||(HomoStr == "Homo3D"))
{
m_session->LoadParameter("HomModesX",m_npointsX);
m_session->LoadParameter("LX", m_LhomX );
m_session->LoadParameter("HomModesY",m_npointsY);
m_session->LoadParameter("LY", m_LhomY );
m_session->LoadParameter("HomModesZ",m_npointsZ);
m_session->LoadParameter("LZ", m_LhomZ );
}
if(m_session->DefinesSolverInfo("USEFFT"))
{
m_useFFT = true;
}
}
else
{
m_npointsZ = 1; // set to default value so can use to identify 2d or 3D (homogeneous) expansions
}
if(m_session->DefinesSolverInfo("PROJECTION"))
{
std::string ProjectStr
= m_session->GetSolverInfo("PROJECTION");
if((ProjectStr == "Continuous")||(ProjectStr == "Galerkin")||
(ProjectStr == "CONTINUOUS")||(ProjectStr == "GALERKIN"))
{
}
else if(ProjectStr == "DisContinuous")
{
}
else
{
ASSERTL0(false,"PROJECTION value not recognised");
}
}
else
{
cerr << "Projection type not specified in SOLVERINFO,"
"defaulting to continuous Galerkin" << endl;
}
ASSERTL0(m_session->DefinesFunction("BaseFlow"),
"Base flow must be defined for linearised forms.");
string file = m_session->GetFunctionFilename("BaseFlow", 0);
//Periodic base flows
if(m_session->DefinesParameter("N_slices"))
{
m_session->LoadParameter("N_slices",m_slices);
if(m_slices>1)
{
DFT(file,m_slices);
}
else
{
ASSERTL0(false,"Number of slices must be a positive number greater than 1");
}
}
//Steady base-flow
else
{
//BaseFlow from file
if (m_session->GetFunctionType("BaseFlow", m_session->GetVariable(0))
{
}
//analytic base flow
else
{
int nq = m_base[0]->GetNpoints();
Array<OneD,NekDouble> x0(nq);
Array<OneD,NekDouble> x1(nq);
Array<OneD,NekDouble> x2(nq);
// get the coordinates (assuming all fields have the same
// discretisation)
m_base[0]->GetCoords(x0,x1,x2);
for(unsigned int i = 0 ; i < m_base.num_elements(); i++)
{
= m_session->GetFunction("BaseFlow", i);
ifunc->Evaluate(x0,x1,x2,m_base[i]->UpdatePhys());
m_base[i]->SetPhysState(true);
m_base[i]->FwdTrans_IterPerExp(m_base[i]->GetPhys(),
m_base[i]->UpdateCoeffs());
}
}
}
}
void Nektar::AdjointAdvection::WriteFldBase ( std::string &  outname)
protected

Write field data to the given filename.

void Nektar::AdjointAdvection::WriteFldBase ( std::string &  outname,
MultiRegions::ExpListSharedPtr field,
Array< OneD, Array< OneD, NekDouble > > &  fieldcoeffs,
Array< OneD, std::string > &  variables 
)
protected

Write input fields to the given filename.

Friends And Related Function Documentation

friend class MemoryManager< AdjointAdvection >
friend

Definition at line 68 of file AdjointAdvection.h.

Member Data Documentation

string Nektar::AdjointAdvection::className = GetAdvectionTermFactory().RegisterCreatorFunction("Adjoint", AdjointAdvection::create)
static

Name of class.

Definition at line 79 of file AdjointAdvection.h.

Array<OneD, MultiRegions::ExpListSharedPtr> Nektar::AdjointAdvection::m_base
protected

Definition at line 83 of file AdjointAdvection.h.

Referenced by DFT(), GetFloquetBlockMatrix(), ImportFldBase(), SetUpBaseFields(), UpdateBase(), v_ComputeAdvectionTerm(), and v_InitObject().

SpatialDomains::BoundaryConditionsSharedPtr Nektar::AdjointAdvection::m_boundaryConditions
private

Definition at line 178 of file AdjointAdvection.h.

Referenced by v_InitObject().

LibUtilities::NektarFFTSharedPtr Nektar::AdjointAdvection::m_FFT
protected

Definition at line 92 of file AdjointAdvection.h.

Referenced by DFT().

FloquetBlockMatrixMapShPtr Nektar::AdjointAdvection::m_FloquetBlockMat
protected

Definition at line 102 of file AdjointAdvection.h.

bool Nektar::AdjointAdvection::m_HalfMode
protected

flag to determine if use half mode or not

Definition at line 97 of file AdjointAdvection.h.

Referenced by ImportFldBase(), SetUpBaseFields(), v_ComputeAdvectionTerm(), and v_InitObject().

int Nektar::AdjointAdvection::m_HomoDirec
private

number of homogenous directions

Definition at line 172 of file AdjointAdvection.h.

Referenced by v_InitObject().

enum HomogeneousType Nektar::AdjointAdvection::m_HomogeneousType
private

Definition at line 162 of file AdjointAdvection.h.

Referenced by SetUpBaseFields(), and v_InitObject().

Array<OneD, Array<OneD, NekDouble> > Nektar::AdjointAdvection::m_interp
protected

Definition at line 90 of file AdjointAdvection.h.

Referenced by DFT(), ImportFldBase(), and v_ComputeAdvectionTerm().

NekDouble Nektar::AdjointAdvection::m_LhomX
private

physical length in X direction (if homogeneous)

Definition at line 164 of file AdjointAdvection.h.

Referenced by v_InitObject().

NekDouble Nektar::AdjointAdvection::m_LhomY
private

physical length in Y direction (if homogeneous)

Definition at line 165 of file AdjointAdvection.h.

Referenced by SetUpBaseFields(), and v_InitObject().

NekDouble Nektar::AdjointAdvection::m_LhomZ
private

physical length in Z direction (if homogeneous)

Definition at line 166 of file AdjointAdvection.h.

Referenced by SetUpBaseFields(), and v_InitObject().

bool Nektar::AdjointAdvection::m_MultipleModes
protected

flag to determine if use multiple mode or not

Definition at line 98 of file AdjointAdvection.h.

Referenced by v_InitObject().

int Nektar::AdjointAdvection::m_npointsX
private

number of points in X direction (if homogeneous)

Definition at line 168 of file AdjointAdvection.h.

Referenced by v_InitObject().

int Nektar::AdjointAdvection::m_npointsY
private

number of points in Y direction (if homogeneous)

Definition at line 169 of file AdjointAdvection.h.

Referenced by SetUpBaseFields(), and v_InitObject().

int Nektar::AdjointAdvection::m_npointsZ
private

number of points in Z direction (if homogeneous)

Definition at line 170 of file AdjointAdvection.h.

Referenced by ImportFldBase(), SetUpBaseFields(), and v_InitObject().

int Nektar::AdjointAdvection::m_NumMode
private

Mode to use in case of single mode analysis.

Definition at line 174 of file AdjointAdvection.h.

NekDouble Nektar::AdjointAdvection::m_period
protected

Definition at line 88 of file AdjointAdvection.h.

Referenced by DFT(), UpdateBase(), and v_ComputeAdvectionTerm().

bool Nektar::AdjointAdvection::m_SingleMode
protected

flag to determine if use single mode or not

Definition at line 96 of file AdjointAdvection.h.

Referenced by ImportFldBase(), SetUpBaseFields(), and v_InitObject().

int Nektar::AdjointAdvection::m_slices
protected

Definition at line 86 of file AdjointAdvection.h.

Referenced by DFT(), GetFloquetBlockMatrix(), UpdateBase(), v_ComputeAdvectionTerm(), and v_InitObject().

Array<OneD,NekDouble> Nektar::AdjointAdvection::m_tmpIN
protected

Definition at line 93 of file AdjointAdvection.h.

Referenced by DFT().

Array<OneD,NekDouble> Nektar::AdjointAdvection::m_tmpOUT
protected

Definition at line 94 of file AdjointAdvection.h.

Referenced by DFT().

bool Nektar::AdjointAdvection::m_useFFT
private

flag to determine if use or not the FFT for transformations

Definition at line 159 of file AdjointAdvection.h.

Referenced by SetUpBaseFields(), and v_InitObject().

bool Nektar::AdjointAdvection::m_useFFTW
protected

Definition at line 95 of file AdjointAdvection.h.