Nektar++
Dummy.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: Dummy.cpp
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2017 Kilian Lackhove
10// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
11// Department of Aeronautics, Imperial College London (UK), and Scientific
12// Computing and Imaging Institute, University of Utah (USA).
13//
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: Dummy Equation System that only sends/receives fields
33//
34///////////////////////////////////////////////////////////////////////////////
35
36#include <iostream>
37
40
41namespace Nektar
42{
43
44std::string Dummy::className =
46 "Dummy", Dummy::create,
47 "Dummy Equation System that only sends/receives fields");
48
51 : UnsteadySystem(pSession, pGraph)
52{
53}
54
55/**
56 * @brief Initialization object for the Dummy class.
57 */
58void Dummy::v_InitObject(bool DeclareFields)
59{
60 UnsteadySystem::v_InitObject(DeclareFields);
61
64
66 m_fields, m_fields.size());
67
68 if (m_session->DefinesElement("Nektar/Coupling"))
69 {
70 TiXmlElement *vCoupling = m_session->GetElement("Nektar/Coupling");
71
72 ASSERTL0(vCoupling->Attribute("TYPE"),
73 "Missing TYPE attribute in Coupling");
74 std::string vType = vCoupling->Attribute("TYPE");
75 ASSERTL0(!vType.empty(),
76 "TYPE attribute must be non-empty in Coupling");
77
78 m_coupling = GetCouplingFactory().CreateInstance(vType, m_fields[0]);
79
80 auto sV = m_session->GetVariables();
81 for (auto const &sendVar : m_coupling->GetSendFieldNames())
82 {
83 auto match = find(sV.begin(), sV.end(), sendVar);
84 if (match != sV.end())
85 {
86 int id = distance(sV.begin(), match);
87 m_intVariables.push_back(id);
88 }
89 }
90 }
91}
92
93/**
94 * @brief Compute the right-hand side.
95 */
98 [[maybe_unused]] const NekDouble time)
99{
100 int nVariables = inarray.size();
101 int nq = GetTotPoints();
102
103 for (int i = 0; i < nVariables; ++i)
104 {
105 Vmath::Zero(nq, outarray[i], 1);
106 }
107}
108
109/**
110 * @brief Compute the projection and call the method for imposing the
111 * boundary conditions in case of discontinuous projection.
112 */
114 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
116 [[maybe_unused]] const NekDouble time)
117{
118 int nvariables = inarray.size();
119 int nq = m_fields[0]->GetNpoints();
120
121 // deep copy
122 if (inarray != outarray)
123 {
124 for (int i = 0; i < nvariables; ++i)
125 {
126 Vmath::Vcopy(nq, inarray[i], 1, outarray[i], 1);
127 }
128 }
129}
130
131/**
132 * @brief v_PreIntegrate
133 */
135{
136 if (m_coupling)
137 {
138 int numForceFields = 0;
139 for (auto &x : m_forcing)
140 {
141 numForceFields += x->GetForces().size();
142 }
143 std::vector<std::string> varNames;
145 numForceFields);
146 for (int i = 0; i < m_fields.size(); ++i)
147 {
148 varNames.push_back(m_session->GetVariable(i));
149 phys[i] = m_fields[i]->UpdatePhys();
150 }
151
152 int f = 0;
153 for (auto &x : m_forcing)
154 {
155 for (int i = 0; i < x->GetForces().size(); ++i)
156 {
157 phys[m_fields.size() + f + i] = x->GetForces()[i];
158 varNames.push_back("F_" + std::to_string(f) + "_" +
159 m_session->GetVariable(i));
160 }
161 f++;
162 }
163
164 m_coupling->Send(step, m_time, phys, varNames);
165 m_coupling->Receive(step, m_time, phys, varNames);
166 }
167
169}
170
171/**
172 * @brief v_PostIntegrate
173 */
175{
176 if (m_coupling && m_coupling->GetSendFieldNames().size() > 0)
177 {
178 LibUtilities::Timer timer1;
179 timer1.Start();
180
181 auto sV = m_session->GetVariables();
182 for (auto const &sendVar : m_coupling->GetSendFieldNames())
183 {
184 auto match = find(sV.begin(), sV.end(), sendVar);
185 if (match != sV.end())
186 {
187 int id = distance(sV.begin(), match);
188 GetFunction("SendFields", m_fields[id])
189 ->Evaluate(sendVar, m_fields[id]->UpdatePhys(), m_time);
190 }
191 }
192
193 timer1.Stop();
194 if (m_session->DefinesCmdLineArgument("verbose"))
195 {
196 std::cout << "Field evaluation time: " << timer1.TimePerTest(1)
197 << std::endl;
198 }
199 }
200
201 for (int i = 0; i < m_session->GetVariables().size(); ++i)
202 {
203
204 m_fields[i]->FwdTransLocalElmt(m_fields[i]->UpdatePhys(),
205 m_fields[i]->UpdateCoeffs());
206 m_fields[i]->SetPhysState(false);
207 }
208
210}
211
213{
214 if (m_coupling)
215 {
216 m_coupling->Finalize();
217 }
218
220
221 int f = 0;
222 for (auto &x : m_forcing)
223 {
224 for (int i = 0; i < x->GetForces().size(); ++i)
225 {
226 int npts = GetTotPoints();
227
228 NekDouble l2err = 0.0;
229 NekDouble linferr = 0.0;
230 for (int j = 0; j < npts; ++j)
231 {
232 l2err += x->GetForces()[i][j] * x->GetForces()[i][j];
233 linferr = std::max(linferr, fabs(x->GetForces()[i][j]));
234 }
235
236 m_comm->AllReduce(l2err, LibUtilities::ReduceSum);
237 m_comm->AllReduce(npts, LibUtilities::ReduceSum);
238 m_comm->AllReduce(linferr, LibUtilities::ReduceMax);
239
240 l2err /= npts;
241 l2err = sqrt(l2err);
242
243 if (m_comm->TreatAsRankZero())
244 {
245 std::cout << "L 2 error (variable "
246 << "F_" + std::to_string(f) + "_" +
247 m_session->GetVariable(i)
248 << ") : " << l2err << std::endl;
249
250 std::cout << "L inf error (variable "
251 << "F_" + std::to_string(f) + "_" +
252 m_session->GetVariable(i)
253 << ") : " << linferr << std::endl;
254 }
255 }
256 f++;
257 }
258}
259
260} // namespace Nektar
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
bool v_PreIntegrate(int step) override
v_PreIntegrate
Definition: Dummy.cpp:134
SolverUtils::CouplingSharedPtr m_coupling
Definition: Dummy.h:68
static EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
Definition: Dummy.h:54
void v_Output() override
Definition: Dummy.cpp:212
bool v_PostIntegrate(int step) override
v_PostIntegrate
Definition: Dummy.cpp:174
std::vector< SolverUtils::ForcingSharedPtr > m_forcing
Definition: Dummy.h:69
Dummy(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Definition: Dummy.cpp:49
void DoOdeProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the projection and call the method for imposing the boundary conditions in case of discontinu...
Definition: Dummy.cpp:113
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the right-hand side.
Definition: Dummy.cpp:96
static std::string className
Name of class.
Definition: Dummy.h:65
void v_InitObject(bool DeclareFields=true) override
Initialization object for the Dummy class.
Definition: Dummy.cpp:58
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
NekDouble TimePerTest(unsigned int n)
Returns amount of seconds per iteration in a test with n iterations.
Definition: Timer.cpp:65
LibUtilities::CommSharedPtr m_comm
Communicator.
NekDouble m_time
Current time of simulation.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
SOLVER_UTILS_EXPORT int GetTotPoints()
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
Get a SessionFunction by name.
virtual SOLVER_UTILS_EXPORT void v_Output(void)
static SOLVER_UTILS_EXPORT std::vector< ForcingSharedPtr > Load(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields=0)
Definition: Forcing.cpp:76
Base class for unsteady solvers.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate(int step)
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate(int step)
SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true) override
Init object for UnsteadySystem class.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
CouplingFactory & GetCouplingFactory()
Declaration of the Coupling factory singleton.
Definition: Coupling.cpp:40
EquationSystemFactory & GetEquationSystemFactory()
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:174
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:475
double NekDouble
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285