Nektar++
Static Public Member Functions | Static Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | Friends | List of all members
Nektar::SolverUtils::ForcingAbsorption Class Reference

#include <ForcingAbsorption.h>

Inheritance diagram for Nektar::SolverUtils::ForcingAbsorption:
[legend]

Static Public Member Functions

static SOLVER_UTILS_EXPORT ForcingSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields, const TiXmlElement *pForce)
 Creates an instance of this class. More...
 
- Static Public Member Functions inherited from Nektar::SolverUtils::Forcing
static SOLVER_UTILS_EXPORT std::vector< ForcingSharedPtrLoad (const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields=0)
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Types

typedef bg::model::point< NekDouble, 3, bg::cs::cartesian > BPoint
 
typedef std::pair< BPoint, unsigned int > BPointPair
 
typedef bgi::rtree< BPointPair, bgi::rstar< 16 > > BRTree
 

Protected Member Functions

void CalculateForcing (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, const NekDouble &time)
 
SOLVER_UTILS_EXPORT void v_InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields, const TiXmlElement *pForce) override
 
SOLVER_UTILS_EXPORT void v_Apply (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time) override
 
SOLVER_UTILS_EXPORT void v_ApplyCoeff (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time) override
 
- Protected Member Functions inherited from Nektar::SolverUtils::Forcing
SOLVER_UTILS_EXPORT Forcing (const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation)
 Constructor. More...
 
virtual SOLVER_UTILS_EXPORT ~Forcing ()=default
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields, const TiXmlElement *pForce)=0
 
virtual SOLVER_UTILS_EXPORT void v_Apply (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)=0
 
virtual SOLVER_UTILS_EXPORT void v_PreApply (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 
virtual SOLVER_UTILS_EXPORT void v_ApplyCoeff (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, std::string pName, bool pCache=false)
 Get a SessionFunction by name. More...
 
SOLVER_UTILS_EXPORT void EvaluateTimeFunction (LibUtilities::SessionReaderSharedPtr pSession, std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, NekDouble pTime=NekDouble(0))
 
SOLVER_UTILS_EXPORT void EvaluateTimeFunction (const NekDouble pTime, const LibUtilities::EquationSharedPtr &pEqn, Array< OneD, NekDouble > &pArray)
 

Protected Attributes

bool m_hasRefFlow
 
bool m_hasRefFlowTime
 
Array< OneD, Array< OneD, NekDouble > > m_Absorption
 
Array< OneD, Array< OneD, NekDouble > > m_Refflow
 
std::string m_funcNameTime
 
std::vector< unsigned int > m_bRegions
 
std::shared_ptr< BRTreem_rtree
 
- Protected Attributes inherited from Nektar::SolverUtils::Forcing
LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
const std::weak_ptr< EquationSystemm_equ
 Weak pointer to equation system using this forcing. More...
 
Array< OneD, Array< OneD, NekDouble > > m_Forcing
 Evaluated forcing function. More...
 
int m_NumVariable
 Number of variables. More...
 
std::map< std::string, SolverUtils::SessionFunctionSharedPtrm_sessionFunctions
 Map of known SessionFunctions. More...
 

Private Member Functions

 ForcingAbsorption (const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation)
 
 ~ForcingAbsorption (void) override=default
 
void CalcAbsorption (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pForce)
 

Private Attributes

bool m_homogeneous
 

Friends

class MemoryManager< ForcingAbsorption >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::Forcing
SOLVER_UTILS_EXPORT void InitObject (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields, const TiXmlElement *pForce)
 Initialise the forcing object. More...
 
SOLVER_UTILS_EXPORT void Apply (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 Apply the forcing. More...
 
SOLVER_UTILS_EXPORT void PreApply (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 Change the advection velocity before applying the forcing. For example, subtracting the frame velocity from the advection velocity in the MovingRefercenceFrame. More...
 
SOLVER_UTILS_EXPORT void ApplyCoeff (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 Apply the forcing. More...
 
SOLVER_UTILS_EXPORT const Array< OneD, const Array< OneD, NekDouble > > & GetForces ()
 
SOLVER_UTILS_EXPORT Array< OneD, Array< OneD, NekDouble > > & UpdateForces ()
 

Detailed Description

Definition at line 55 of file ForcingAbsorption.h.

Member Typedef Documentation

◆ BPoint

typedef bg::model::point<NekDouble, 3, bg::cs::cartesian> Nektar::SolverUtils::ForcingAbsorption::BPoint
protected

Definition at line 78 of file ForcingAbsorption.h.

◆ BPointPair

typedef std::pair<BPoint, unsigned int> Nektar::SolverUtils::ForcingAbsorption::BPointPair
protected

Definition at line 79 of file ForcingAbsorption.h.

◆ BRTree

typedef bgi::rtree<BPointPair, bgi::rstar<16> > Nektar::SolverUtils::ForcingAbsorption::BRTree
protected

Definition at line 80 of file ForcingAbsorption.h.

Constructor & Destructor Documentation

◆ ForcingAbsorption()

Nektar::SolverUtils::ForcingAbsorption::ForcingAbsorption ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::weak_ptr< EquationSystem > &  pEquation 
)
private

Definition at line 47 of file ForcingAbsorption.cpp.

50 : Forcing(pSession, pEquation), m_hasRefFlow(false), m_hasRefFlowTime(false)
51{
52}
SOLVER_UTILS_EXPORT Forcing(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation)
Constructor.
Definition: Forcing.cpp:47

◆ ~ForcingAbsorption()

Nektar::SolverUtils::ForcingAbsorption::~ForcingAbsorption ( void  )
overrideprivatedefault

Member Function Documentation

◆ CalcAbsorption()

void Nektar::SolverUtils::ForcingAbsorption::CalcAbsorption ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const TiXmlElement *  pForce 
)
private

Definition at line 110 of file ForcingAbsorption.cpp.

114{
115 const TiXmlElement *funcNameElmt = pForce->FirstChildElement("COEFF");
116 ASSERTL0(funcNameElmt,
117 "Requires COEFF tag, specifying function "
118 "name which prescribes absorption layer coefficient.");
119 std::string funcName = funcNameElmt->GetText();
120 ASSERTL0(m_session->DefinesFunction(funcName),
121 "Function '" + funcName + "' not defined.");
122
123 int npts = pFields[0]->GetTotPoints();
124
125 m_Absorption = Array<OneD, Array<OneD, NekDouble>>(m_NumVariable);
126 for (int i = 0; i < m_NumVariable; ++i)
127 {
128 m_Absorption[i] = Array<OneD, NekDouble>(npts, 0.0);
129 }
130
131 funcNameElmt = pForce->FirstChildElement("BOUNDARYREGIONS");
132 if (funcNameElmt)
133 {
134 ASSERTL0(
135 ParseUtils::GenerateVector(funcNameElmt->GetText(), m_bRegions),
136 "Unable to process list of BOUNDARYREGIONS in Absorption "
137 "Forcing: " +
138 std::string(funcNameElmt->GetText()));
139
140 // alter m_bRegions so that it contains the boundaryRegions of this rank
141 std::vector<unsigned int> localBRegions;
142 SpatialDomains::BoundaryConditions bcs(m_session,
143 pFields[0]->GetGraph());
145 bcs.GetBoundaryRegions();
146 SpatialDomains::BoundaryRegionCollection::iterator it1;
147 int n = 0;
148 for (it1 = regions.begin(); it1 != regions.end(); ++it1)
149 {
150 if (std::find(m_bRegions.begin(), m_bRegions.end(), it1->first) !=
151 m_bRegions.end())
152 {
153 localBRegions.push_back(n);
154 }
155 n++;
156 }
157 m_bRegions = localBRegions;
158
159 if (m_bRegions.size() == 0)
160 {
161 return;
162 }
163
164 std::vector<Array<OneD, const NekDouble>> points;
165
166 Array<OneD, Array<OneD, NekDouble>> x(3);
167 for (int i = 0; i < 3; i++)
168 {
169 x[i] = Array<OneD, NekDouble>(npts, 0.0);
170 }
171 pFields[0]->GetCoords(x[0], x[1], x[2]);
172 for (int i = 0; i < 3; i++)
173 {
174 points.push_back(x[i]);
175 }
176
177 Array<OneD, NekDouble> t(npts, 0.0);
178 points.push_back(t);
179
180 Array<OneD, NekDouble> r(npts, 0.0);
181 std::vector<unsigned int>::iterator it;
182 std::vector<BPointPair> inPoints;
183 Array<OneD, Array<OneD, NekDouble>> b(3);
184 for (it = m_bRegions.begin(); it != m_bRegions.end(); ++it)
185 {
186 int bpts = pFields[0]->GetBndCondExpansions()[*it]->GetNpoints();
187 for (int i = 0; i < 3; i++)
188 {
189 b[i] = Array<OneD, NekDouble>(bpts, 0.0);
190 }
191 pFields[0]->GetBndCondExpansions()[*it]->GetCoords(b[0], b[1],
192 b[2]);
193 for (int i = 0;
194 i < pFields[0]->GetBndCondExpansions()[*it]->GetNpoints(); ++i)
195 {
196 inPoints.push_back(
197 BPointPair(BPoint(b[0][i], b[1][i], b[2][i]), i));
198 }
199 }
201 m_rtree->insert(inPoints.begin(), inPoints.end());
202
203 for (int i = 0; i < npts; ++i)
204 {
205 std::vector<BPointPair> result;
206 BPoint sPoint(x[0][i], x[1][i], x[2][i]);
207 m_rtree->query(bgi::nearest(sPoint, 1), std::back_inserter(result));
208 r[i] = bg::distance(sPoint, result[0].first);
209 }
210 points.push_back(r);
211
212 std::string s_FieldStr;
213 for (int i = 0; i < m_NumVariable; ++i)
214 {
215 s_FieldStr = m_session->GetVariable(i);
216 ASSERTL0(m_session->DefinesFunction(funcName, s_FieldStr),
217 "Variable '" + s_FieldStr + "' not defined.");
218
220 m_session->GetFunction(funcName, s_FieldStr);
221 ASSERTL0(ffunc->GetVlist() == "x y z t r",
222 "EVARS of " + funcName + " must be 'r'");
223
224 ffunc->Evaluate(points, m_Absorption[i]);
225 }
226 }
227 else
228 {
229 for (int i = 0; i < m_NumVariable; ++i)
230 {
231 std::string s_FieldStr = m_session->GetVariable(i);
232 GetFunction(pFields, m_session, funcName)
233 ->Evaluate(s_FieldStr, m_Absorption[i]);
234 }
235 }
236
237 // If homogeneous expansion is used, transform the forcing term to
238 // be in the Fourier space
239 if (m_homogeneous)
240 {
241 for (int i = 0; i < m_NumVariable; ++i)
242 {
243 pFields[i]->HomogeneousFwdTrans(pFields[i]->GetTotPoints(),
245 }
246 }
247}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
static bool GenerateVector(const std::string &str, std::vector< T > &out)
Takes a comma-separated string and converts it to entries in a vector.
Definition: ParseUtils.cpp:130
bg::model::point< NekDouble, 3, bg::cs::cartesian > BPoint
std::pair< BPoint, unsigned int > BPointPair
std::vector< unsigned int > m_bRegions
Array< OneD, Array< OneD, NekDouble > > m_Absorption
int m_NumVariable
Number of variables.
Definition: Forcing.h:129
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, std::string pName, bool pCache=false)
Get a SessionFunction by name.
Definition: Forcing.cpp:150
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: Forcing.h:123
std::shared_ptr< Equation > EquationSharedPtr
Definition: Equation.h:125
std::map< int, BoundaryRegionShPtr > BoundaryRegionCollection
Definition: Conditions.h:211
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:475

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::StdRegions::find(), Nektar::ParseUtils::GenerateVector(), Nektar::SpatialDomains::BoundaryConditions::GetBoundaryRegions(), Nektar::SolverUtils::Forcing::GetFunction(), m_Absorption, m_bRegions, m_homogeneous, Nektar::SolverUtils::Forcing::m_NumVariable, m_rtree, and Nektar::SolverUtils::Forcing::m_session.

Referenced by v_InitObject().

◆ CalculateForcing()

void Nektar::SolverUtils::ForcingAbsorption::CalculateForcing ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
const NekDouble time 
)
protected

Definition at line 278 of file ForcingAbsorption.cpp.

281{
282 int nq = m_Forcing[0].size();
283
284 std::string s_FieldStr;
285 Array<OneD, NekDouble> TimeScale(1);
286 Array<OneD, Array<OneD, NekDouble>> RefflowScaled(m_NumVariable);
287
288 if (m_hasRefFlow)
289 {
290 for (int i = 0; i < m_NumVariable; i++)
291 {
292 RefflowScaled[i] = Array<OneD, NekDouble>(nq);
294 {
295 s_FieldStr = m_session->GetVariable(i);
296
297 std::string s_FieldStr = m_session->GetVariable(i);
299 ->Evaluate(s_FieldStr, m_Refflow[i], time);
300 Vmath::Vcopy(nq, m_Refflow[i], 1, RefflowScaled[i], 1);
301 }
302 else
303 {
304 Vmath::Vcopy(nq, m_Refflow[i], 1, RefflowScaled[i], 1);
305 }
306
307 Vmath::Vsub(nq, inarray[i], 1, RefflowScaled[i], 1, m_Forcing[i],
308 1);
309 Vmath::Vmul(nq, m_Absorption[i], 1, m_Forcing[i], 1, m_Forcing[i],
310 1);
311 }
312 }
313 else
314 {
315 for (int i = 0; i < m_NumVariable; i++)
316 {
317 Vmath::Vmul(nq, m_Absorption[i], 1, inarray[i], 1, m_Forcing[i], 1);
318 }
319 }
320}
Array< OneD, Array< OneD, NekDouble > > m_Refflow
Array< OneD, Array< OneD, NekDouble > > m_Forcing
Evaluated forcing function.
Definition: Forcing.h:127
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.
Definition: Vmath.hpp:72
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
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.
Definition: Vmath.hpp:220

References Nektar::SolverUtils::Forcing::GetFunction(), m_Absorption, Nektar::SolverUtils::Forcing::m_Forcing, m_funcNameTime, m_hasRefFlow, m_hasRefFlowTime, Nektar::SolverUtils::Forcing::m_NumVariable, m_Refflow, Nektar::SolverUtils::Forcing::m_session, Vmath::Vcopy(), Vmath::Vmul(), and Vmath::Vsub().

Referenced by v_Apply(), and v_ApplyCoeff().

◆ create()

static SOLVER_UTILS_EXPORT ForcingSharedPtr Nektar::SolverUtils::ForcingAbsorption::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::weak_ptr< EquationSystem > &  pEquation,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const unsigned int &  pNumForcingFields,
const TiXmlElement *  pForce 
)
inlinestatic

Creates an instance of this class.

Definition at line 61 of file ForcingAbsorption.h.

66 {
69 pEquation);
70 p->InitObject(pFields, pNumForcingFields, pForce);
71 return p;
72 }
SOLVER_UTILS_EXPORT typedef std::shared_ptr< Forcing > ForcingSharedPtr
A shared pointer to an EquationSystem object.
Definition: Forcing.h:53

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::SolverUtils::ForcingSharedPtr, and CellMLToNektar.cellml_metadata::p.

◆ v_Apply()

void Nektar::SolverUtils::ForcingAbsorption::v_Apply ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble time 
)
overrideprotectedvirtual

Implements Nektar::SolverUtils::Forcing.

Definition at line 249 of file ForcingAbsorption.cpp.

253{
254 int nq = m_Forcing[0].size();
255 CalculateForcing(fields, inarray, time);
256 for (int i = 0; i < m_NumVariable; ++i)
257 {
258 Vmath::Vadd(nq, m_Forcing[i], 1, outarray[i], 1, outarray[i], 1);
259 }
260}
void CalculateForcing(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, const NekDouble &time)
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.hpp:180

References CalculateForcing(), Nektar::SolverUtils::Forcing::m_Forcing, Nektar::SolverUtils::Forcing::m_NumVariable, and Vmath::Vadd().

◆ v_ApplyCoeff()

void Nektar::SolverUtils::ForcingAbsorption::v_ApplyCoeff ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble time 
)
overrideprotectedvirtual

Reimplemented from Nektar::SolverUtils::Forcing.

Definition at line 262 of file ForcingAbsorption.cpp.

266{
267 int ncoeff = outarray[m_NumVariable - 1].size();
268 Array<OneD, NekDouble> tmp(ncoeff, 0.0);
269 CalculateForcing(fields, inarray, time);
270
271 for (int i = 0; i < m_NumVariable; ++i)
272 {
273 fields[i]->FwdTrans(m_Forcing[i], tmp);
274 Vmath::Vadd(ncoeff, tmp, 1, outarray[i], 1, outarray[i], 1);
275 }
276}

References CalculateForcing(), Nektar::SolverUtils::Forcing::m_Forcing, Nektar::SolverUtils::Forcing::m_NumVariable, and Vmath::Vadd().

◆ v_InitObject()

void Nektar::SolverUtils::ForcingAbsorption::v_InitObject ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const unsigned int &  pNumForcingFields,
const TiXmlElement *  pForce 
)
overrideprotectedvirtual

Implements Nektar::SolverUtils::Forcing.

Definition at line 54 of file ForcingAbsorption.cpp.

57{
58 m_NumVariable = pNumForcingFields;
59 int npts = pFields[0]->GetTotPoints();
60
61 // Check for homogeneous expansion
62 // for later transformation of evaluated Absorption (SessionFunction)
63 m_homogeneous = pFields[0]->GetExpType() == MultiRegions::e3DH1D ||
64 pFields[0]->GetExpType() == MultiRegions::e3DH2D;
65
66 CalcAbsorption(pFields, pForce);
67
68 m_Forcing = Array<OneD, Array<OneD, NekDouble>>(m_NumVariable);
69 for (int i = 0; i < m_NumVariable; ++i)
70 {
71 m_Forcing[i] = Array<OneD, NekDouble>(npts, 0.0);
72 }
73
74 const TiXmlElement *funcNameElmt = pForce->FirstChildElement("REFFLOW");
75 if (funcNameElmt)
76 {
77 std::string funcName = funcNameElmt->GetText();
78 ASSERTL0(m_session->DefinesFunction(funcName),
79 "Function '" + funcName + "' not defined.");
80 m_Refflow = Array<OneD, Array<OneD, NekDouble>>(m_NumVariable);
81 for (int i = 0; i < m_NumVariable; ++i)
82 {
83 std::string s_FieldStr = m_session->GetVariable(i);
84 ASSERTL0(m_session->DefinesFunction(funcName, s_FieldStr),
85 "Variable '" + s_FieldStr + "' not defined.");
86 m_Refflow[i] = Array<OneD, NekDouble>(npts, 0.0);
87 GetFunction(pFields, m_session, funcName)
88 ->Evaluate(s_FieldStr, m_Refflow[i]);
89 }
90 m_hasRefFlow = true;
91 }
92
93 funcNameElmt = pForce->FirstChildElement("REFFLOWTIME");
94 if (funcNameElmt)
95 {
96 m_funcNameTime = funcNameElmt->GetText();
97 m_hasRefFlowTime = true;
98 m_hasRefFlow = true;
99 m_Refflow = Array<OneD, Array<OneD, NekDouble>>(m_NumVariable);
100 for (int i = 0; i < m_NumVariable; ++i)
101 {
102 std::string s_FieldStr = m_session->GetVariable(i);
103 ASSERTL0(m_session->DefinesFunction(m_funcNameTime, s_FieldStr),
104 "Variable '" + s_FieldStr + "' not defined.");
105 m_Refflow[i] = Array<OneD, NekDouble>(npts, 0.0);
106 }
107 }
108}
void CalcAbsorption(const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const TiXmlElement *pForce)

References ASSERTL0, CalcAbsorption(), Nektar::MultiRegions::e3DH1D, Nektar::MultiRegions::e3DH2D, Nektar::SolverUtils::Forcing::GetFunction(), Nektar::SolverUtils::Forcing::m_Forcing, m_funcNameTime, m_hasRefFlow, m_hasRefFlowTime, m_homogeneous, Nektar::SolverUtils::Forcing::m_NumVariable, m_Refflow, and Nektar::SolverUtils::Forcing::m_session.

Friends And Related Function Documentation

◆ MemoryManager< ForcingAbsorption >

friend class MemoryManager< ForcingAbsorption >
friend

Definition at line 1 of file ForcingAbsorption.h.

Member Data Documentation

◆ className

std::string Nektar::SolverUtils::ForcingAbsorption::className
static
Initial value:
=
"Absorption", ForcingAbsorption::create, "Forcing Absorption")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static SOLVER_UTILS_EXPORT ForcingSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const unsigned int &pNumForcingFields, const TiXmlElement *pForce)
Creates an instance of this class.
ForcingFactory & GetForcingFactory()
Declaration of the forcing factory singleton.
Definition: Forcing.cpp:41

Name of the class.

Definition at line 75 of file ForcingAbsorption.h.

◆ m_Absorption

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::ForcingAbsorption::m_Absorption
protected

Definition at line 84 of file ForcingAbsorption.h.

Referenced by CalcAbsorption(), and CalculateForcing().

◆ m_bRegions

std::vector<unsigned int> Nektar::SolverUtils::ForcingAbsorption::m_bRegions
protected

Definition at line 87 of file ForcingAbsorption.h.

Referenced by CalcAbsorption().

◆ m_funcNameTime

std::string Nektar::SolverUtils::ForcingAbsorption::m_funcNameTime
protected

Definition at line 86 of file ForcingAbsorption.h.

Referenced by CalculateForcing(), and v_InitObject().

◆ m_hasRefFlow

bool Nektar::SolverUtils::ForcingAbsorption::m_hasRefFlow
protected

Definition at line 82 of file ForcingAbsorption.h.

Referenced by CalculateForcing(), and v_InitObject().

◆ m_hasRefFlowTime

bool Nektar::SolverUtils::ForcingAbsorption::m_hasRefFlowTime
protected

Definition at line 83 of file ForcingAbsorption.h.

Referenced by CalculateForcing(), and v_InitObject().

◆ m_homogeneous

bool Nektar::SolverUtils::ForcingAbsorption::m_homogeneous
private

Definition at line 114 of file ForcingAbsorption.h.

Referenced by CalcAbsorption(), and v_InitObject().

◆ m_Refflow

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::ForcingAbsorption::m_Refflow
protected

Definition at line 85 of file ForcingAbsorption.h.

Referenced by CalculateForcing(), and v_InitObject().

◆ m_rtree

std::shared_ptr<BRTree> Nektar::SolverUtils::ForcingAbsorption::m_rtree
protected

Definition at line 88 of file ForcingAbsorption.h.

Referenced by CalcAbsorption().