Nektar++
FilterMean.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: FilterMean.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: Output mean.
32//
33///////////////////////////////////////////////////////////////////////////////
34
35#include <iomanip>
36
39
40using namespace std;
41
42namespace Nektar::SolverUtils
43{
44std::string FilterMean::className =
47
49 const std::shared_ptr<EquationSystem> &pEquation,
50 const ParamMap &pParams)
51 : Filter(pSession, pEquation), m_index(-1), m_homogeneous(false), m_planes()
52{
53 // OutputFile
54 std::string ext = ".avg";
55 m_outputFile = Filter::SetupOutput(ext, pParams);
56
57 // OutputFrequency
58 auto it = pParams.find("OutputFrequency");
59 ASSERTL0(it != pParams.end(), "Missing parameter 'OutputFrequency'.");
60 LibUtilities::Equation equ(m_session->GetInterpreter(), it->second);
61 m_outputFrequency = round(equ.Evaluate());
62
63 pSession->LoadParameter("LZ", m_homogeneousLength, 0.0);
64}
65
67{
68}
69
72 [[maybe_unused]] const NekDouble &time)
73{
75 int spacedim = 2;
76
77 ASSERTL0(pFields[0]->GetExpType() != MultiRegions::e1D,
78 "1D expansion not supported for mean filter");
79
80 ASSERTL0(pFields[0]->GetExpType() != MultiRegions::e3DH2D,
81 "Homogeneous 2D expansion not supported for mean filter");
82
83 // Lock equation system pointer
84 auto equ = m_equ.lock();
85 ASSERTL0(equ, "Weak pointer expired");
86
87 auto fluidEqu = std::dynamic_pointer_cast<FluidInterface>(equ);
88 ASSERTL0(fluidEqu, "Mean filter is incompatible with this solver.");
89
90 if (pFields[0]->GetExpType() == MultiRegions::e3DH1D)
91 {
92 m_homogeneous = true;
93 spacedim = 3;
94 }
95
96 // Calculate area/volume of domain.
97 if (m_homogeneous)
98 {
99 m_planes = pFields[0]->GetZIDs();
100 areaField = pFields[0]->GetPlane(0);
101 }
102 else
103 {
104 areaField = pFields[0];
105 }
106
107 Array<OneD, NekDouble> inarray(areaField->GetNpoints(), 1.0);
108 m_area = areaField->Integral(inarray);
109
110 if (m_homogeneous)
111 {
113 }
114
115 // Open OutputFile
116 std::string volname[3] = {"length", "area", "volume"};
117 LibUtilities::CommSharedPtr vComm = pFields[0]->GetComm();
118 if (vComm->GetRank() == 0)
119 {
120 m_outputStream.open(m_outputFile.c_str());
122 "Unable to open: '" + m_outputFile + "'");
123 m_outputStream.setf(ios::scientific, ios::floatfield);
124 m_outputStream << "# Time";
125 for (int i = 0; i < pFields.size(); ++i)
126 {
127 m_outputStream << setw(22) << equ->GetVariable(i);
128 }
129 m_outputStream << setw(22) << volname[spacedim - 1] << " " << m_area;
130 m_outputStream << endl;
131 }
132
133 // Output values at initial time.
134 m_index = 0;
136 {
137 v_Update(pFields, time);
138 }
139}
140
143 const NekDouble &time)
144{
145 // Only output every m_outputFrequency
146 if ((m_index++) % m_outputFrequency)
147 {
148 return;
149 }
150
151 int i;
152 LibUtilities::CommSharedPtr vComm = pFields[0]->GetComm();
153
154 // Calculate average values.
155 Array<OneD, NekDouble> avg(pFields.size());
156 for (i = 0; i < pFields.size(); ++i)
157 {
158 avg[i] = 0.0;
159 }
160
161 if (m_homogeneous)
162 {
163 for (i = 0; i < pFields.size(); ++i)
164 {
165 avg[i] = pFields[0]->GetPlane(0)->Integral(pFields[i]->GetPhys()) *
167 }
168 }
169 else
170 {
171 for (i = 0; i < pFields.size(); ++i)
172 {
173 avg[i] = pFields[0]->Integral(pFields[i]->GetPhys());
174 }
175 }
176
177 for (i = 0; i < pFields.size(); ++i)
178 {
179 avg[i] /= m_area;
180 }
181
182 if (vComm->GetRank() == 0)
183 {
184 m_outputStream << setw(17) << setprecision(8) << time;
185 for (int i = 0; i < pFields.size(); ++i)
186 {
187 m_outputStream << setw(22) << setprecision(11) << avg[i];
188 }
189 m_outputStream << endl;
190 }
191}
192
195 &pFields,
196 [[maybe_unused]] const NekDouble &time)
197{
198 if (pFields[0]->GetComm()->GetRank() == 0)
199 {
200 m_outputStream.close();
201 }
202}
203
205{
206 return true;
207}
208
209} // namespace Nektar::SolverUtils
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
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
const std::weak_ptr< EquationSystem > m_equ
Definition: Filter.h:94
std::map< std::string, std::string > ParamMap
Definition: Filter.h:66
SOLVER_UTILS_EXPORT void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
Definition: FilterMean.cpp:141
SOLVER_UTILS_EXPORT void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
Definition: FilterMean.cpp:193
SOLVER_UTILS_EXPORT bool v_IsTimeDependent() override
Definition: FilterMean.cpp:204
SOLVER_UTILS_EXPORT void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pField, const NekDouble &time) override
Definition: FilterMean.cpp:70
Array< OneD, unsigned int > m_planes
Definition: FilterMean.h:86
static SolverUtils::FilterSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< SolverUtils::EquationSystem > &pEquation, const ParamMap &pParams)
Creates an instance of this class.
Definition: FilterMean.h:46
SOLVER_UTILS_EXPORT ~FilterMean() override
Definition: FilterMean.cpp:66
SOLVER_UTILS_EXPORT FilterMean(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< EquationSystem > &pEquation, const ParamMap &pParams)
Definition: FilterMean.cpp:48
static std::string className
Name of the class.
Definition: FilterMean.h:58
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:55
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
FilterFactory & GetFilterFactory()
double NekDouble
STL namespace.