Nektar++
FilterOffsetPhase.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File FilterOffsetPhase.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// License for the specific language governing rights and limitations under
14// Permission is hereby granted, free of charge, to any person obtaining a
15// copy of this software and associated documentation files (the "Software"),
16// to deal in the Software without restriction, including without limitation
17// the rights to use, copy, modify, merge, publish, distribute, sublicense,
18// and/or sell copies of the Software, and to permit persons to whom the
19// Software is furnished to do so, subject to the following conditions:
20//
21// The above copyright notice and this permission notice shall be included
22// in all copies or substantial portions of the Software.
23//
24// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30// DEALINGS IN THE SOFTWARE.
31//
32// Description: Outputs phase fields during time-stepping.
33//
34///////////////////////////////////////////////////////////////////////////////
35
37
38namespace Nektar
39{
42 "OffsetPhase", FilterOffsetPhase::create);
43
46 const std::shared_ptr<SolverUtils::EquationSystem> &pEquation,
47 const ParamMap &pParams)
48 : Filter(pSession, pEquation)
49{
50 std::string ext = "";
51 m_outputFile = Filter::SetupOutput(ext, pParams);
52
53 ASSERTL0(pParams.find("OutputFrequency") != pParams.end(),
54 "Missing parameter 'OutputFrequency'.");
55 LibUtilities::Equation equ(m_session->GetInterpreter(),
56 pParams.find("OutputFrequency")->second);
57 m_outputFrequency = floor(equ.Evaluate());
58
59 if (pParams.find("MeanV") != pParams.end())
60 {
61 LibUtilities::Equation equ(m_session->GetInterpreter(),
62 pParams.find("MeanV")->second);
63 meanV = equ.Evaluate();
64 }
65
66 m_outputIndex = 0;
67 m_index = 0;
69}
70
72{
73}
74
77 [[maybe_unused]] const NekDouble &time)
78{
79 m_index = 0;
80 m_outputIndex = 0;
81
82 oldV = pFields[0]->GetPhys(); // Allocate storage for past data
83
84 v_Update(pFields, 0.0);
85}
86
89 const NekDouble &time)
90{
91 if (m_index++ % m_outputFrequency > 0)
92 {
93 return;
94 }
95
96 std::stringstream vTmpFilename;
97 std::string vOutputFilename;
98 std::string ext = ".chk";
99
100 vTmpFilename << fs::path(m_outputFile).replace_extension("").string() << "_"
101 << m_outputIndex << ext;
102 vOutputFilename = Filter::SetupOutput(ext, vTmpFilename.str());
103
104 // Get FieldDefinitions
105 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef =
106 pFields[0]->GetFieldDefinitions();
107 std::vector<std::vector<NekDouble>> FieldData(FieldDef.size());
108
109 Array<OneD, NekDouble> currentV = pFields[0]->GetPhys();
110 Array<OneD, NekDouble> phase(currentV.size());
111 Array<OneD, NekDouble> phase_coeff(pFields[0]->GetNcoeffs());
112
113 // Calculate phase angle in physical space:
114 for (int i = 0; i < currentV.size(); ++i)
115 {
116 phase[i] = atan2(currentV[i] - meanV,
117 oldV[i] - meanV); // subtract mean -57.2986
118 }
119
120 // Convert to coefficient space:
121 pFields[0]->FwdTransLocalElmt(phase, phase_coeff);
122
123 // for each item in FieldDef
124 // push "phase" to FieldDef[i]->m_fields
125 // call AppendFieldData(FieldDef[i], FieldData[i], data)
126 for (int i = 0; i < FieldDef.size(); ++i)
127 {
128 FieldDef[i]->m_fields.push_back("phase");
129 pFields[0]->AppendFieldData(FieldDef[i], FieldData[i], phase_coeff);
130 }
131
132 // Save currentV
133 oldV = currentV;
134
135 // Update time in field info if required
136 LibUtilities::FieldMetaDataMap fieldMetaDataMap;
137 fieldMetaDataMap["Time"] = boost::lexical_cast<std::string>(time);
138
139 m_fld->Write(vOutputFilename, FieldDef, FieldData, fieldMetaDataMap);
141}
142
145 &pFields,
146 [[maybe_unused]] const NekDouble &time)
147{
148}
149
151{
152 return true;
153}
154} // namespace Nektar
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
FilterOffsetPhase(const LibUtilities::SessionReaderSharedPtr &pSession, const std::shared_ptr< SolverUtils::EquationSystem > &pEquation, const ParamMap &pParams)
static std::string className
Name of the class.
Array< OneD, NekDouble > oldV
void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) 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.
LibUtilities::FieldIOSharedPtr m_fld
void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
static std::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
Definition: FieldIO.cpp:194
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
std::map< std::string, std::string > ParamMap
Definition: Filter.h:66
std::map< std::string, std::string > FieldMetaDataMap
Definition: FieldIO.h:50
std::shared_ptr< SessionReader > SessionReaderSharedPtr
FilterFactory & GetFilterFactory()
double NekDouble