Nektar++
FilterEnergy1D.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: FilterEnergy1D.cpp
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description: Outputs orthogonal expansion of 1D elements.
32//
33///////////////////////////////////////////////////////////////////////////////
34
37
38using namespace std;
39
40namespace Nektar::SolverUtils
41{
42std::string FilterEnergy1D::className =
45
46/**
47 * @brief Set up filter with output file and frequency parameters.
48 *
49 * @param pSession Current session.
50 * @param pParams Map of parameters defined in XML file.
51 */
54 const std::shared_ptr<EquationSystem> &pEquation, const ParamMap &pParams)
55 : Filter(pSession, pEquation), m_index(0)
56{
57 ASSERTL0(pSession->GetComm()->GetSize() == 1,
58 "The 1D energy filter currently only works in serial.");
59
60 // OutputFile
61 std::string outName;
62 std::string ext = ".eny";
63 outName = Filter::SetupOutput(ext, pParams);
64
65 m_out.open(outName.c_str());
66
67 // OutputFrequency
68 auto it = pParams.find("OutputFrequency");
69 if (it == pParams.end())
70 {
72 }
73 else
74 {
75 LibUtilities::Equation equ(m_session->GetInterpreter(), it->second);
76 m_outputFrequency = round(equ.Evaluate());
77 }
78}
79
80/**
81 * @brief Destructor.
82 */
84{
85}
86
87/**
88 * @brief Initialize filter.
89 */
92 [[maybe_unused]] const NekDouble &time)
93{
94 ASSERTL0(pFields[0]->GetExp(0)->GetNumBases() == 1,
95 "The Energy 1D filter is only valid in 1D.");
96}
97
98/**
99 * @brief Update filter output with the current timestep's orthogonal
100 * coefficients.
101 */
104 const NekDouble &time)
105{
106 // Only output every m_outputFrequency
107 if ((m_index++) % m_outputFrequency)
108 {
109 return;
110 }
111
112 int nElmt = pFields[0]->GetExpSize();
113
114 // Loop over all elements
115 m_out << "##" << endl;
116 m_out << "## Time = " << time << endl;
117 m_out << "##" << endl;
118
119 for (int i = 0; i < nElmt; ++i)
120 {
121 // Figure out number of modes in this expansion.
122 LocalRegions::ExpansionSharedPtr exp = pFields[0]->GetExp(i);
123 int nModes = exp->GetBasis(0)->GetNumModes();
124
125 // Set uo basis key for orthogonal basis
127 LibUtilities::BasisKey bkeyOrth(btype, nModes,
128 exp->GetBasis(0)->GetPointsKey());
129
130 // Get basis key for existing expansion
131 LibUtilities::BasisKey bkey(exp->GetBasis(0)->GetBasisType(),
132 exp->GetBasis(0)->GetNumModes(),
133 exp->GetBasis(0)->GetPointsKey());
134
135 // Find coeffs for this element in the list of all coefficients
137 pFields[0]->GetCoeffs() + pFields[0]->GetCoeff_Offset(i);
138
139 // Storage for orthogonal coefficients
140 Array<OneD, NekDouble> coeffsOrth(nModes);
141
142 // Project from coeffs -> orthogonal coeffs
143 LibUtilities::InterpCoeff1D(bkey, coeffs, bkeyOrth, coeffsOrth);
144
145 // Write coeffs to file
146 m_out << "# Element " << i << " (ID " << exp->GetGeom()->GetGlobalID()
147 << ")" << endl;
148 for (int j = 0; j < nModes; ++j)
149 {
150 m_out << coeffsOrth[j] << endl;
151 }
152 }
153 m_out << endl;
154}
155
158 &pFields,
159 [[maybe_unused]] const NekDouble &time)
160{
161 m_out.close();
162}
163
165{
166 return true;
167}
168} // namespace Nektar::SolverUtils
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
Describes the specification for a Basis.
Definition: Basis.h:45
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
Initialize filter.
static std::string className
Name of the class.
SOLVER_UTILS_EXPORT ~FilterEnergy1D() override
Destructor.
unsigned int m_outputFrequency
Output frequency.
SOLVER_UTILS_EXPORT FilterEnergy1D(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const ParamMap &pParams)
Set up filter with output file and frequency parameters.
void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
Update filter output with the current timestep's orthogonal coefficients.
static FilterSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const std::map< std::string, std::string > &pParams)
Creates an instance of this class.
unsigned int m_index
Current index counter.
void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
std::ofstream m_out
Output file.
SOLVER_UTILS_EXPORT std::string SetupOutput(const std::string ext, const ParamMap &pParams)
Definition: Filter.h:139
LibUtilities::SessionReaderSharedPtr m_session
Definition: Filter.h:93
std::map< std::string, std::string > ParamMap
Definition: Filter.h:66
std::shared_ptr< SessionReader > SessionReaderSharedPtr
void InterpCoeff1D(const BasisKey &fbasis0, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, Array< OneD, NekDouble > &to)
Definition: InterpCoeff.cpp:42
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:42
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:66
FilterFactory & GetFilterFactory()
double NekDouble
STL namespace.