Nektar++
SmoothedProfileMethod.cpp
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: SmoothedProfileMethod.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: Smoothed Profile Method for the Incompressible
33// Navier Stokes equations
34///////////////////////////////////////////////////////////////////////////////
35
41
42namespace Nektar
43{
44
45using namespace MultiRegions;
46
49 "SmoothedProfileMethod", SmoothedProfileMethod::create);
50
53 "SolverType", "SmoothedProfileMethod", eSmoothedProfileMethod);
54
55/**
56 * @brief Construct a new Smoothed Profile Method object
57 *
58 * @param pSession
59 * @param pGraph
60 */
64 : UnsteadySystem(pSession, pGraph),
65 VelocityCorrectionScheme(pSession, pGraph)
66{
67}
68
70{
72
73 // Update implicit time-intregration class operators
76
77 // Number of dims as number of velocity vectors
78 size_t nvel = m_velocity.size();
79
80 // Initialization of correction pressure and shape function
81 switch (nvel)
82 {
83 case 1:
85 {
86 SetUpExpansions<ContField>(nvel);
87 }
89 {
90 SetUpExpansions<DisContField>(nvel);
91 }
92 break;
93
94 case 2:
96 {
97 SetUpExpansions<ContField>(nvel);
98 }
100 {
101 SetUpExpansions<DisContField>(nvel);
102 }
103 break;
104
105 case 3:
107 {
109 {
110 SetUpExpansions<ContField>(nvel);
111 }
113 {
114 SetUpExpansions<ContField3DHomogeneous1D>(nvel);
115 }
118 {
119 SetUpExpansions<ContField3DHomogeneous2D>(nvel);
120 }
121 }
123 {
125 {
126 SetUpExpansions<DisContField>(nvel);
127 }
129 {
130 SetUpExpansions<DisContField3DHomogeneous1D>(nvel);
131 }
134 {
135 SetUpExpansions<DisContField3DHomogeneous2D>(nvel);
136 }
137 }
138 break;
139 }
140
141 // Read 'm_phi' and its velocity
142 ASSERTL0(m_session->DefinesFunction("ShapeFunction"),
143 "ShapeFunction must be defined in the session file.")
144 ReadPhi();
145
146 // Allocate the vector 'm_up'
147 size_t physTot = m_phi->GetTotPoints();
148 m_velName.push_back("Up");
149 if (nvel > 1)
150 {
151 m_velName.push_back("Vp");
152 }
153 if (nvel == 3)
154 {
155 m_velName.push_back("Wp");
156 }
157
159 for (size_t i = 0; i < nvel; ++i)
160 {
161 m_up[i] = Array<OneD, NekDouble>(physTot, 0.0);
162 }
163
164 // Make sure that m_phi and m_up are defined
165 UpdatePhiUp(0.0);
166
167 // Get the time integration scheme.
169 if (m_session->DefinesTimeIntScheme())
170 {
171 timeInt = m_session->GetTimeIntScheme();
172 }
173 else
174 {
175 timeInt.method = m_session->GetSolverInfo("TimeIntegrationMethod");
176 timeInt.order = timeInt.method.back() - '0';
177
178 // Remove everything past the IMEX.
179 timeInt.method = timeInt.method.substr(0, 4);
180 }
181
182 // Select 'm_gamma0' depending on IMEX order
183 ASSERTL0(
184 boost::iequals(timeInt.method, "IMEX") && 1 <= timeInt.order &&
185 timeInt.order <= 4,
186 "The TimeIntegrationMethod scheme must be IMEX with order '1' to '4'.")
187
188 switch (timeInt.order)
189 {
190 case 1:
191 m_gamma0 = 1.0;
192 break;
193
194 case 2:
195 m_gamma0 = 3.0 / 2.0;
196 break;
197
198 case 3:
199 m_gamma0 = 11.0 / 6.0;
200 break;
201
202 case 4:
203 m_gamma0 = 25.0 / 12.0;
204 break;
205 }
206
207 // Check if the aeroforces filter is active, negative if inactive
208 m_forcesFilter = -1;
209 for (size_t i = 0; i < m_session->GetFilters().size(); ++i)
210 {
211 if (boost::iequals(m_session->GetFilters()[i].name, "AeroForcesSPM"))
212 {
213 m_forcesFilter = i;
214 break;
215 }
216 }
217}
218
219/**
220 * @brief Generates the summary of the current simulation
221 *
222 * @param s
223 */
225{
227 SolverUtils::AddSummaryItem(s, "IB formulation",
228 "Smoothed Profile Method (SPM)");
229}
230
231/**
232 * @brief Linear terms due to pressure and visosity are calculated here.
233 * After solving the velocity filed without taking into account the
234 * immersed boundaries, a new correction is applied through the force
235 * \f$f_s\f$:
236 *
237 * \f[ \mathbf{f_s} = \frac{\Phi^{n+1}(\mathbf{u_p}-\mathbf{u^*})}
238 * {\Delta t} \f]
239 *
240 * @param inarray
241 * @param outarray
242 * @param time
243 * @param a_iixDt
244 */
246 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
247 Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble time,
248 const NekDouble a_iixDt)
249{
251 time, a_iixDt);
252
253 size_t physTot = m_pressureP->GetNpoints();
254
255 /* SPM correction of velocity */
256 // Update 'm_phi' and 'm_up' if needed (evaluated at next time step)
257 UpdatePhiUp(time + a_iixDt);
258 // Update calculation of IB forcing 'm_fs'
259 UpdateForcing(outarray, a_iixDt);
260 // Estimate forces only if requested
261 if (m_forcesFilter >= 0)
262 {
263 std::static_pointer_cast<FilterAeroForcesSPM>(
265 ->CalculateForces(outarray, m_upPrev, m_phi, time, a_iixDt);
266 }
267 // Set BC conditions for pressure p_p
268 SetUpCorrectionPressure(outarray, m_F);
269 // Solve Poisson equation for pressure p_p
271 // Solve velocity in the next step with IB
272 SolveCorrectedVelocity(m_F, outarray, a_iixDt);
273
274 // Add pressures to get final value
275 Vmath::Vadd(physTot, m_pressure->GetPhys(), 1, m_pressureP->GetPhys(), 1,
276 m_pressure->UpdatePhys(), 1);
277 m_pressure->FwdTrans(m_pressure->GetPhys(), m_pressure->UpdateCoeffs());
278
279 // Add presure to outflow bc if using convective like BCs
280 m_extrapolation->AddPressureToOutflowBCs(m_kinvis);
281}
282
283/**
284 * @brief Sets the forcing term of the equation for the correction pressure
285 * \f$p_p\f$:
286 *
287 * \f[ \nabla\cdot\mathbf{f_s} \f]
288 *
289 * @param fields
290 * @param Forcing
291 */
293 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>> &fields,
295{
296 size_t physTot = m_fs[0]->GetNpoints();
297 size_t nvel = m_velocity.size();
298
299 // Set boundary conditions
301
302 // Divergence of 'fs'
303 m_fields[m_velocity[0]]->PhysDeriv(eX, m_fs[0]->GetPhys(), Forcing[0]);
304
305 // Using 'Forcing[1]' as storage
306 for (size_t i = 1; i < nvel; ++i)
307 {
308 size_t ind = m_velocity[i];
309 m_fields[ind]->PhysDeriv(DirCartesianMap[i], m_fs[i]->GetPhys(),
310 Forcing[1]);
311 Vmath::Vadd(physTot, Forcing[1], 1, Forcing[0], 1, Forcing[0], 1);
312 }
313}
314
315/**
316 * @brief Solves the Poisson equation for the correction pressure
317 * \f$p_p\f$:
318 *
319 * \f[ \nabla^2 p_p = \nabla\cdot\mathbf{f_s} \f]
320 *
321 * @param Forcing
322 */
325{
327 // Factor 'lambda=0' in Helmholtz equation to get the Poisson form
329
330 // Solve the Poisson equation
331 m_pressureP->HelmSolve(Forcing, m_pressureP->UpdateCoeffs(), factors);
332
333 // Update node values from coefficients
334 m_pressureP->BwdTrans(m_pressureP->GetCoeffs(), m_pressureP->UpdatePhys());
335}
336
337/**
338 * @brief Corrects the velocity field so that the IBs are taken into
339 * account. Solves the explicit equation:
340 *
341 * \f[ \frac{\gamma_0(\mathbf{u_p}^{n+1} - \mathbf{u}^*)}{\Delta t} =
342 * \mathbf{f_s} - \nabla p_p \f]
343 *
344 * @param Forcing
345 * @param fields
346 * @param dt
347 */
351{
352 size_t physTot = m_phi->GetNpoints();
353
354 // Gradient of p_p
355 size_t nvel = m_velocity.size();
356 if (nvel == 2)
357 {
358 m_pressureP->PhysDeriv(m_pressureP->GetPhys(), Forcing[0], Forcing[1]);
359 }
360 else
361 {
362 m_pressureP->PhysDeriv(m_pressureP->GetPhys(), Forcing[0], Forcing[1],
363 Forcing[2]);
364 }
365
366 // Velocity correction
367 for (size_t i = 0; i < nvel; ++i)
368 {
369 // Adding -(1-m_phi)*grad(p_p) instead of -grad(p_p) reduces the
370 // flux through the walls, but the flow is not incompressible
371 if (m_session->DefinesSolverInfo("ForceBoundary") &&
372 boost::iequals(m_session->GetSolverInfo("ForceBoundary"), "True"))
373 {
374 Vmath::Vvtvm(physTot, m_phi->GetPhys(), 1, Forcing[i], 1,
375 Forcing[i], 1, Forcing[i], 1);
376 Vmath::Vadd(physTot, m_fs[i]->GetPhys(), 1, Forcing[i], 1,
377 Forcing[i], 1);
378 }
379 else
380 {
381 Vmath::Vsub(physTot, m_fs[i]->GetPhys(), 1, Forcing[i], 1,
382 Forcing[i], 1);
383 }
384 Blas::Daxpy(physTot, dt / m_gamma0, Forcing[i], 1, fields[i], 1);
385 }
386}
387
388/**
389 * @brief Updates the BCs for boundaries with Neumann or periodic BCs in
390 * the pressure:
391 *
392 * \f[ \frac{\partial p_p}{\partial\mathbf{n}} =
393 * \mathbf{f_s}\cdot\mathbf{n} \f]
394 */
396{
397 size_t nvel = m_velocity.size();
400
401 // Get the BC expansions
402 BndExp = m_pressureP->GetBndCondExpansions();
403 BndCond = m_pressureP->GetBndConditions();
404
405 // For each boundary...
406 for (size_t b = 0; b < BndExp.size(); ++b)
407 {
408 // Only for BCs based on the derivative
409 if (BndCond[b]->GetBoundaryConditionType() ==
411 BndCond[b]->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
412 {
413 // Calculate f_s values
415 for (size_t i = 0; i < nvel; ++i)
416 {
417 f_s[i] = m_fs[0]->GetBndCondExpansions()[b]->GetPhys();
418 }
419
420 // BC is f_s * n
421 BndExp[b]->NormVectorIProductWRTBase(f_s, BndExp[b]->UpdatePhys());
422 }
423 }
424}
425
426/**
427 * @brief Calculates the values of the shape function
428 *
429 * @param t
430 */
432{
433 // Initialise 'm_up' and 'm_phi' during first step
434 if (t <= 0.0)
435 {
436 if (!m_filePhi)
437 {
438 // Update 'm_phi' only if it was provided as a function
439 m_phiEvaluator->Evaluate("Phi", m_phi->UpdatePhys(), t);
440 }
441
442 // Initialize both variables for the first step
443 m_phiEvaluator->Evaluate(m_velName, m_up, t);
444
445 // Initialise 'm_upPrev' in all cases
446 m_upPrev = m_up;
447 }
448 // If timedependent 'm_phi'
449 // Phi functions from files are not timedependent
450 else if (m_timeDependentPhi)
451 {
452 m_phiEvaluator->Evaluate("Phi", m_phi->UpdatePhys(), t);
453
454 // And if velocities are timedependent as well
456 {
457 // Store previous value of u_p during simulation
458 m_upPrev = m_up;
459 m_phiEvaluator->Evaluate(m_velName, m_up, t);
460 }
461 }
462}
463
464/**
465 * @brief For a body with a velocity \f$\mathbf{u_p}\f$, the force
466 * \f$\mathbf{f_s}\f$ applied to the fluid ensures that the IBC are met:
467 *
468 * \f[ \mathbf{f_s} = \frac{\Phi^{n+1}\left(\mathbf{u_p}^{n+1} -
469 * \mathbf{u^*}\right)}{\Delta t} \f]
470 *
471 * @param fields
472 * @param dt
473 * @param f_s
474 */
476 const Array<OneD, const Array<OneD, NekDouble>> &fields, NekDouble dt)
477{
478 size_t nvel = m_velocity.size();
479 size_t nq = m_phi->GetNpoints();
480
481 for (size_t i = 0; i < nvel; ++i)
482 {
483 // In homogeneous cases, switch out of wave space
484 Array<OneD, NekDouble> tmpField(nq);
485 size_t ind = m_velocity[i];
486
488 m_fields[ind]->GetWaveSpace())
489 {
490 m_fields[ind]->HomogeneousBwdTrans(nq, fields[i], tmpField);
491 m_fs[i]->HomogeneousBwdTrans(nq, m_fs[i]->GetPhys(),
492 m_fs[i]->UpdatePhys());
493 }
494 else
495 {
496 tmpField = fields[i];
497 }
498
499 Vmath::Vsub(nq, m_up[i], 1, tmpField, 1, m_fs[i]->UpdatePhys(), 1);
500 Vmath::Vmul(nq, m_phi->GetPhys(), 1, m_fs[i]->GetPhys(), 1,
501 m_fs[i]->UpdatePhys(), 1);
502 Vmath::Smul(nq, m_gamma0 / dt, m_fs[i]->GetPhys(), 1,
503 m_fs[i]->UpdatePhys(), 1);
504
505 // And go back to wave space if the 'if' was executed
507 m_fields[ind]->GetWaveSpace())
508 {
509 m_fs[i]->HomogeneousFwdTrans(nq, m_fs[i]->GetPhys(),
510 m_fs[i]->UpdatePhys());
511 }
512 }
513}
514
515/**
516 * @brief True if the function is timedependent, false otherwise
517 *
518 * @param name
519 * @param type
520 * @param attribute
521 * @return string
522 */
524 std::string elemName)
525{
526 // Get the handler of the function
527 TiXmlElement *function = GetFunctionHdl(funcName);
528
529 // Go to the first element
530 TiXmlElement *functionDef = function->FirstChildElement();
531 ASSERTL0(functionDef, "At least one element must be defined in " + funcName)
532
533 // And search the element with name 'elemName'
534 std::string varName = functionDef->Attribute("VAR");
535 while (functionDef && !boost::iequals(varName, elemName))
536 {
537 functionDef = functionDef->NextSiblingElement();
538 varName = functionDef->Attribute("VAR");
539 }
540
541 ASSERTL0(functionDef,
542 "Variable " + elemName + " must be defined in " + funcName + ".");
543
544 // And return the value of USERDEFINEDTYPE
545 std::string attr;
546 int err = functionDef->QueryStringAttribute("USERDEFINEDTYPE", &attr);
547 bool output = boost::iequals(attr, "TimeDependent");
548
549 ASSERTL0((err == TIXML_NO_ATTRIBUTE) || (err == TIXML_SUCCESS && output),
550 "USERDEFINEDTYPE in " + elemName +
551 " must be TimeDependent if defined");
552
553 return output;
554}
555
556/**
557 * @brief Returns a handle to the requested function. Returns NULL if it
558 * does not exist
559 *
560 * @param functionName
561 * @return TiXmlElement*
562 */
563TiXmlElement *SmoothedProfileMethod::GetFunctionHdl(std::string functionName)
564{
565 // Get the handler of first function block
566 TiXmlElement *conds = m_session->GetElement("Nektar/Conditions");
567 TiXmlElement *function = conds->FirstChildElement("FUNCTION");
568
569 // Loop over functions until the block 'name' is found
570 std::string functionType = function->Attribute("NAME");
571 while (function && !boost::iequals(functionType, functionName))
572 {
573 function = function->NextSiblingElement("FUNCTION");
574 functionType = function->Attribute("NAME");
575 }
576
577 return function;
578}
579
581{
582 // Function evaluator for Phi and Up
583 m_phiEvaluator = GetFunction("ShapeFunction");
584
585 TiXmlElement *function = GetFunctionHdl("ShapeFunction");
586 TiXmlElement *child = function->FirstChildElement();
587 m_filePhi = false;
588
589 // If defined by using a file
590 if (boost::iequals(child->ValueStr(), "F"))
591 {
592 // Get name of STL file
593 std::string fileName;
594 int status = child->QueryStringAttribute("FILE", &fileName);
595 ASSERTL0(status == TIXML_SUCCESS,
596 "An FLD file with the values "
597 "of the phi function has to be supplied.")
598 ASSERTL0(boost::iequals(fileName.substr(fileName.length() - 4), ".fld"),
599 "A valid FLD file must be supplied in the "
600 "'ShapeFunction' field.")
601
602 // Get phi values from XML file (after "FieldConvert" the STL file)
603 // First, load the data
604 std::vector<LibUtilities::FieldDefinitionsSharedPtr> fieldDef;
605 std::vector<std::vector<NekDouble>> fieldData;
606 LibUtilities::FieldMetaDataMap fieldMetaData;
609 phiFile->Import(fileName, fieldDef, fieldData, fieldMetaData);
610
611 // Only Phi field should be defined in the file
612 ASSERTL0(fieldData.size() == 1, "Only one field (phi) must be "
613 "defined in the FLD file.")
614
615 // Extract Phi field to output
616 std::string tmp("phi");
617 m_phi->ExtractDataToCoeffs(fieldDef[0], fieldData[0], tmp,
618 m_phi->UpdateCoeffs());
619 m_phi->BwdTrans(m_phi->GetCoeffs(), m_phi->UpdatePhys());
620 m_filePhi = true;
621 m_timeDependentPhi = false;
622 m_timeDependentUp = false;
623 }
624 else
625 {
626 // Check if Phi is timedependent
627 m_timeDependentPhi = GetVarTimeDependence("ShapeFunction", "Phi");
628
629 // If so, check if its velocity changes as well
630 m_timeDependentUp = GetVarTimeDependence("ShapeFunction", "Up");
631 switch (m_velocity.size())
632 {
633 case 2:
635 GetVarTimeDependence("ShapeFunction", "Vp");
636 break;
637 case 3:
639 GetVarTimeDependence("ShapeFunction", "Vp");
641 GetVarTimeDependence("ShapeFunction", "Wp");
642 break;
643 }
644 }
645}
646
647} // namespace Nektar
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
MultiRegions::ExpListSharedPtr m_pressure
Pointer to field holding pressure field.
NekDouble m_kinvis
Kinematic viscosity.
ExtrapolateSharedPtr m_extrapolation
Array< OneD, int > m_velocity
int which identifies which components of m_fields contains the velocity (u,v,w);
static std::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:223
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static std::string RegisterEnumValue(std::string pEnum, std::string pString, int pEnumValue)
Registers an enumeration value.
void DefineImplicitSolve(FuncPointerT func, ObjectPointerT obj)
void UpdateForcing(const Array< OneD, const Array< OneD, NekDouble > > &fields, NekDouble dt)
For a body with a velocity , the force applied to the fluid ensures that the IBC are met:
bool m_filePhi
Flag indicating that phi was defined in a file.
void v_GenerateSummary(SolverUtils::SummaryList &s) override
Generates the summary of the current simulation.
void SetCorrectionPressureBCs()
Updates the BCs for boundaries with Neumann or periodic BCs in the pressure:
void v_SolveUnsteadyStokesSystem(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, NekDouble time, NekDouble a_iixDt) override
Linear terms due to pressure and visosity are calculated here. After solving the velocity filed witho...
Array< OneD, Array< OneD, NekDouble > > m_up
Velocity of the immersed body(ies)
std::vector< std::string > m_velName
Vector storing the names of the components of \u_p.
NekDouble m_gamma0
Stiffly-stable scheme coefficient.
int m_forcesFilter
Position of "AeroForcesSPM" filter in 'm_session->GetFilters()'.
void v_InitObject(bool DeclareField=true) override
Initialisation object for EquationSystem.
bool m_timeDependentPhi
Flag that is true when phi depends on time.
static SolverUtils::EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
SmoothedProfileMethod(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Construct a new Smoothed Profile Method object.
void UpdatePhiUp(NekDouble time)
Calculates the values of the shape function.
SolverUtils::SessionFunctionSharedPtr m_phiEvaluator
Function that evaluates the values of \Phi.
void SolveCorrectedVelocity(Array< OneD, Array< OneD, NekDouble > > &Forcing, Array< OneD, Array< OneD, NekDouble > > &fields, NekDouble dt)
Corrects the velocity field so that the IBs are taken into account. Solves the explicit equation:
TiXmlElement * GetFunctionHdl(std::string functionName)
Returns a handle to the requested function. Returns NULL if it does not exist.
static std::string className
Name of class.
bool m_timeDependentUp
Flag signaling if depends on time.
MultiRegions::ExpListSharedPtr m_pressureP
Correction pressure field for SPM.
bool GetVarTimeDependence(std::string funcName, std::string attrName)
True if the function is timedependent, false otherwise.
Array< OneD, Array< OneD, NekDouble > > m_upPrev
Array< OneD, MultiRegions::ExpListSharedPtr > m_fs
Forcing function 'f_s'.
MultiRegions::ExpListSharedPtr m_phi
Shape function 'phi' as expansion list.
void SolveCorrectionPressure(const Array< OneD, NekDouble > &Forcing)
Solves the Poisson equation for the correction pressure :
void SetUpCorrectionPressure(const Array< OneD, const Array< OneD, NekDouble > > &fields, Array< OneD, Array< OneD, NekDouble > > &Forcing)
Sets the forcing term of the equation for the correction pressure :
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
enum HomogeneousType m_HomogeneousType
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
Get a SessionFunction by name.
Defines a forcing term to be explicitly applied.
Definition: Forcing.h:71
Base class for unsteady solvers.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
void v_GenerateSummary(SolverUtils::SummaryList &s) override
Print a summary of time stepping parameters.
Array< OneD, Array< OneD, NekDouble > > m_F
virtual void v_SolveUnsteadyStokesSystem(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time, const NekDouble a_iixDt)
void v_InitObject(bool DeclareField=true) override
Initialisation object for EquationSystem.
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
Definition: Blas.hpp:135
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:322
std::map< std::string, std::string > FieldMetaDataMap
Definition: FieldIO.h:50
std::shared_ptr< SessionReader > SessionReaderSharedPtr
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:87
std::vector< std::pair< std::string, std::string > > SummaryList
Definition: Misc.h:46
EquationSystemFactory & GetEquationSystemFactory()
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:47
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:174
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:430
StdRegions::ConstFactorMap factors
@ eSmoothedProfileMethod
double NekDouble
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 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
void Vvtvm(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvm (vector times vector minus vector): z = w*x - y
Definition: Vmath.hpp:381
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100
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