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

#include <FilterAverageFields.h>

Inheritance diagram for Nektar::SolverUtils::FilterAverageFields:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::FilterAverageFields:
Collaboration graph
[legend]

Public Member Functions

SOLVER_UTILS_EXPORT FilterAverageFields (const LibUtilities::SessionReaderSharedPtr &pSession, const std::map< std::string, std::string > &pParams)
SOLVER_UTILS_EXPORT ~FilterAverageFields ()
- Public Member Functions inherited from Nektar::SolverUtils::Filter
SOLVER_UTILS_EXPORT Filter (const LibUtilities::SessionReaderSharedPtr &pSession)
virtual SOLVER_UTILS_EXPORT ~Filter ()
SOLVER_UTILS_EXPORT void Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
SOLVER_UTILS_EXPORT void Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
SOLVER_UTILS_EXPORT void Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
SOLVER_UTILS_EXPORT bool IsTimeDependent ()

Static Public Member Functions

static FilterSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const std::map< std::string, std::string > &pParams)
 Creates an instance of this class.

Static Public Attributes

static std::string className = GetFilterFactory().RegisterCreatorFunction("AverageFields", FilterAverageFields::create)
 Name of the class.

Protected Member Functions

virtual void v_Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
virtual void v_Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
virtual void v_Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
virtual bool v_IsTimeDependent ()
void OutputAvgField (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, int dump=-1)

Private Attributes

unsigned int m_numAverages
unsigned int m_outputFrequency
unsigned int m_sampleFrequency
unsigned int m_index
unsigned int m_outputIndex
std::string m_outputFile
LibUtilities::FieldIOSharedPtr m_fld
LibUtilities::FieldMetaDataMap m_avgFieldMetaData
Array< OneD, Array< OneD,
NekDouble > > 
m_avgFields

Friends

class MemoryManager< FilterAverageFields >

Additional Inherited Members

- Protected Attributes inherited from Nektar::SolverUtils::Filter
LibUtilities::SessionReaderSharedPtr m_session

Detailed Description

Definition at line 45 of file FilterAverageFields.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::FilterAverageFields::FilterAverageFields ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::map< std::string, std::string > &  pParams 
)

Definition at line 44 of file FilterAverageFields.cpp.

References ASSERTL0, m_fld, m_index, m_numAverages, m_outputFile, m_outputFrequency, m_outputIndex, m_sampleFrequency, and Nektar::SolverUtils::Filter::m_session.

:
Filter(pSession)
{
if (pParams.find("OutputFile") == pParams.end())
{
m_outputFile = m_session->GetSessionName();
}
else
{
ASSERTL0(!(pParams.find("OutputFile")->second.empty()),
"Missing parameter 'OutputFile'.");
m_outputFile = pParams.find("OutputFile")->second;
}
if(pParams.find("SampleFrequency") == pParams.end())
{
}
else
{
m_sampleFrequency = atoi(pParams.find("SampleFrequency")->second.c_str());
}
if(pParams.find("OutputFrequency") == pParams.end())
{
m_outputFrequency = m_session->GetParameter("NumSteps");
}
else
{
m_outputFrequency = atoi(pParams.find("OutputFrequency")->second.c_str());
}
m_index = 0;
}
Nektar::SolverUtils::FilterAverageFields::~FilterAverageFields ( )

Definition at line 85 of file FilterAverageFields.cpp.

{
}

Member Function Documentation

static FilterSharedPtr Nektar::SolverUtils::FilterAverageFields::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::map< std::string, std::string > &  pParams 
)
inlinestatic

Creates an instance of this class.

Definition at line 51 of file FilterAverageFields.h.

void Nektar::SolverUtils::FilterAverageFields::OutputAvgField ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
int  dump = -1 
)
protected

Definition at line 129 of file FilterAverageFields.cpp.

References m_avgFieldMetaData, m_avgFields, m_fld, m_numAverages, m_outputFile, Nektar::SolverUtils::Filter::m_session, and Vmath::Smul().

Referenced by v_Finalise(), and v_Update().

{
for(int n = 0; n < m_avgFields.num_elements(); ++n)
{
Vmath::Smul(m_avgFields[n].num_elements(), 1.0/m_numAverages,m_avgFields[n],1,m_avgFields[n],1);
}
std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
= pFields[0]->GetFieldDefinitions();
std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
Array<OneD, NekDouble> fieldcoeffs;
int ncoeffs = pFields[0]->GetNcoeffs();
// copy Data into FieldData and set variable
for(int j = 0; j < pFields.num_elements(); ++j)
{
// check to see if field is same order a zeroth field
if(m_avgFields[j].num_elements() == ncoeffs)
{
fieldcoeffs = m_avgFields[j];
}
else
{
fieldcoeffs = Array<OneD,NekDouble>(ncoeffs);
pFields[0]->ExtractCoeffsToCoeffs(pFields[j],m_avgFields[j],fieldcoeffs);
}
for(int i = 0; i < FieldDef.size(); ++i)
{
// Could do a search here to find correct variable
FieldDef[i]->m_fields.push_back(m_session->GetVariable(j));
pFields[0]->AppendFieldData(FieldDef[i], FieldData[i], fieldcoeffs);
}
}
m_avgFieldMetaData["NumberOfFieldDumps"] = boost::lexical_cast<std::string>(m_numAverages);
std::stringstream outname;
if(dump == -1) // final dump
{
outname << m_outputFile << "_avg.fld";
}
else
{
outname << m_outputFile<< "_" << dump << "_avg.fld";
}
m_fld->Write(outname.str(),FieldDef,FieldData,m_avgFieldMetaData);
if(dump != -1) // not final dump so rescale cummulative average
{
for(int n = 0; n < m_avgFields.num_elements(); ++n)
{
}
}
}
void Nektar::SolverUtils::FilterAverageFields::v_Finalise ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Implements Nektar::SolverUtils::Filter.

Definition at line 124 of file FilterAverageFields.cpp.

References OutputAvgField().

{
OutputAvgField(pFields);
}
void Nektar::SolverUtils::FilterAverageFields::v_Initialise ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Implements Nektar::SolverUtils::Filter.

Definition at line 89 of file FilterAverageFields.cpp.

References m_avgFieldMetaData, and m_avgFields.

{
m_avgFields = Array<OneD, Array<OneD, NekDouble> >(pFields.num_elements());
for(int n =0; n < pFields.num_elements(); ++n)
{
m_avgFields[n] = Array<OneD, NekDouble>(pFields[n]->GetNcoeffs(),0.0);
}
m_avgFieldMetaData["InitialTime"] = boost::lexical_cast<std::string>(time);
}
bool Nektar::SolverUtils::FilterAverageFields::v_IsTimeDependent ( )
protectedvirtual

Implements Nektar::SolverUtils::Filter.

Definition at line 189 of file FilterAverageFields.cpp.

{
return true;
}
void Nektar::SolverUtils::FilterAverageFields::v_Update ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Implements Nektar::SolverUtils::Filter.

Definition at line 99 of file FilterAverageFields.cpp.

References m_avgFieldMetaData, m_avgFields, m_index, m_numAverages, m_outputFrequency, m_outputIndex, m_sampleFrequency, OutputAvgField(), and Vmath::Vadd().

{
{
return;
}
for(int n = 0; n < pFields.num_elements(); ++n)
{
Vmath::Vadd(m_avgFields[n].num_elements(),
pFields[n]->GetCoeffs(),1,m_avgFields[n],1,
m_avgFields[n],1);
}
// update FinalTime here since this is when last field will be added.
m_avgFieldMetaData["FinalTime"] = boost::lexical_cast<std::string>(time);
{
}
}

Friends And Related Function Documentation

friend class MemoryManager< FilterAverageFields >
friend

Definition at line 48 of file FilterAverageFields.h.

Member Data Documentation

std::string Nektar::SolverUtils::FilterAverageFields::className = GetFilterFactory().RegisterCreatorFunction("AverageFields", FilterAverageFields::create)
static

Name of the class.

Definition at line 59 of file FilterAverageFields.h.

LibUtilities::FieldMetaDataMap Nektar::SolverUtils::FilterAverageFields::m_avgFieldMetaData
private

Definition at line 82 of file FilterAverageFields.h.

Referenced by OutputAvgField(), v_Initialise(), and v_Update().

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::FilterAverageFields::m_avgFields
private

Definition at line 83 of file FilterAverageFields.h.

Referenced by OutputAvgField(), v_Initialise(), and v_Update().

LibUtilities::FieldIOSharedPtr Nektar::SolverUtils::FilterAverageFields::m_fld
private

Definition at line 81 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and OutputAvgField().

unsigned int Nektar::SolverUtils::FilterAverageFields::m_index
private

Definition at line 78 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and v_Update().

unsigned int Nektar::SolverUtils::FilterAverageFields::m_numAverages
private

Definition at line 75 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), OutputAvgField(), and v_Update().

std::string Nektar::SolverUtils::FilterAverageFields::m_outputFile
private

Definition at line 80 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and OutputAvgField().

unsigned int Nektar::SolverUtils::FilterAverageFields::m_outputFrequency
private

Definition at line 76 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and v_Update().

unsigned int Nektar::SolverUtils::FilterAverageFields::m_outputIndex
private

Definition at line 79 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and v_Update().

unsigned int Nektar::SolverUtils::FilterAverageFields::m_sampleFrequency
private

Definition at line 77 of file FilterAverageFields.h.

Referenced by FilterAverageFields(), and v_Update().