49 const std::shared_ptr<EquationSystem> &pEquation,
51 :
Filter(pSession, pEquation), m_index(-1), m_homogeneous(false), m_planes()
54 std::string ext =
".eny";
57 m_comm = pSession->GetComm();
58 if (
m_comm->GetRank() == 0)
62 m_outFile.setf(ios::scientific, ios::floatfield);
64 <<
"Enstrophy" << endl
65 <<
"# ---------------------------------------------"
66 <<
"--------------" << endl;
71 auto it = pParams.find(
"OutputFrequency");
72 ASSERTL0(it != pParams.end(),
"Missing parameter 'OutputFrequency'.");
89 "1D expansion not supported for energy filter");
92 "2D expansion not supported for energy filter");
95 "Homogeneous 2D expansion not supported for energy filter");
106 areaField = pFields[0]->GetPlane(0);
110 areaField = pFields[0];
114 m_area = areaField->Integral(inarray);
130 int i, nPoints = pFields[0]->GetNpoints();
140 auto equ =
m_equ.lock();
141 ASSERTL0(equ,
"Weak pointer expired");
143 auto fluidEqu = std::dynamic_pointer_cast<FluidInterface>(equ);
144 ASSERTL0(fluidEqu,
"Energy filter is incompatible with this solver.");
148 for (i = 0; i < pFields.size(); ++i)
150 physfields[i] = pFields[i]->GetPhys();
158 for (i = 0; i < 3; ++i)
162 fluidEqu->GetVelocity(physfields, u);
164 for (i = 0; i < 3; ++i)
168 pFields[i]->HomogeneousBwdTrans(nPoints, u[i], u[i]);
171 Vmath::Vvtvp(nPoints, u[i], 1, u[i], 1, tmp, 1, tmp, 1);
174 if (!fluidEqu->HasConstantDensity())
177 fluidEqu->GetDensity(physfields, density);
184 pFields[0]->HomogeneousFwdTrans(nPoints, tmp, tmp2);
189 Ek = pFields[0]->Integral(tmp);
194 if (
m_comm->GetRank() == 0)
196 m_outFile << setw(17) << setprecision(8) << time << setw(22)
197 << setprecision(11) << Ek;
200 bool waveSpace[3] = {pFields[0]->GetWaveSpace(), pFields[1]->GetWaveSpace(),
201 pFields[2]->GetWaveSpace()};
205 for (i = 0; i < 3; ++i)
207 pFields[i]->SetWaveSpace(
false);
214 for (i = 0; i < 3; ++i)
216 int f1 = (i + 2) % 3, c2 = f1;
217 int c1 = (i + 1) % 3, f2 = c1;
218 pFields[f1]->PhysDeriv(c1, u[f1], tmp2);
219 pFields[f2]->PhysDeriv(c2, u[f2], tmp3);
221 Vmath::Vvtvp(nPoints, tmp2, 1, tmp2, 1, tmp, 1, tmp, 1);
224 if (!fluidEqu->HasConstantDensity())
231 for (i = 0; i < 3; ++i)
233 pFields[i]->SetWaveSpace(waveSpace[i]);
235 pFields[0]->HomogeneousFwdTrans(nPoints, tmp, tmp);
240 Ek = pFields[0]->Integral(tmp);
245 if (
m_comm->GetRank() == 0)
247 m_outFile << setw(22) << setprecision(11) << Ek << endl;
#define ASSERTL0(condition, msg)
NekDouble Evaluate() const
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
SOLVER_UTILS_EXPORT void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
SOLVER_UTILS_EXPORT void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
static std::string className
Name of the class.
SOLVER_UTILS_EXPORT void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
LibUtilities::CommSharedPtr m_comm
SOLVER_UTILS_EXPORT FilterEnergy(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const ParamMap &pParams)
Array< OneD, unsigned int > m_planes
NekDouble m_homogeneousLength
SOLVER_UTILS_EXPORT ~FilterEnergy() override
static SolverUtils::FilterSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< SolverUtils::EquationSystem > &pEquation, const ParamMap &pParams)
Creates an instance of this class.
unsigned int m_outputFrequency
SOLVER_UTILS_EXPORT bool v_IsTimeDependent() override
SOLVER_UTILS_EXPORT std::string SetupOutput(const std::string ext, const ParamMap &pParams)
LibUtilities::SessionReaderSharedPtr m_session
const std::weak_ptr< EquationSystem > m_equ
std::map< std::string, std::string > ParamMap
bool m_updateOnInitialise
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
FilterFactory & GetFilterFactory()
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 Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
void Zero(int n, T *x, const int incx)
Zero vector.
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.