Nektar++
SubSteppingExtrapolate.h
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: SubSteppingExtrapolate.h
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: Abstract base class for SubSteppingExtrapolate.
32//
33///////////////////////////////////////////////////////////////////////////////
34
35#ifndef NEKTAR_SOLVERS_SUBSTEPPINGEXTRAPOLATE_H
36#define NEKTAR_SOLVERS_SUBSTEPPINGEXTRAPOLATE_H
37
46
47namespace Nektar
48{
49//--------
50// Substepping
51// --------
52
53class SubSteppingExtrapolate;
54
55typedef std::shared_ptr<SubSteppingExtrapolate> SubSteppingExtrapolateSharedPtr;
56
58{
59public:
61
62 /// Creates an instance of this class
67 const SolverUtils::AdvectionSharedPtr &advObject)
68 {
71 pSession, pFields, pPressure, pVel, advObject);
72 return p;
73 }
74
75 /// Name of class
76 static std::string className;
77
78protected:
82
85
90
94 const Array<OneD, int> pVel,
95 const SolverUtils::AdvectionSharedPtr advObject);
96
97 ~SubSteppingExtrapolate() override = default;
98
100 const Array<OneD, const Array<OneD, NekDouble>> &fields,
101 const Array<OneD, const Array<OneD, NekDouble>> &N,
102 NekDouble kinvis) override;
103
105 const LibUtilities::TimeIntegrationSchemeSharedPtr &IntegrationScheme)
106 override;
107
108 void v_SubStepSaveFields(int nstep) override;
109
111 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
112 NekDouble Aii_Dt, NekDouble kinvis) override;
113
114 void v_AccelerationBDF(Array<OneD, Array<OneD, NekDouble>> &array) override;
115
116 void v_SubStepAdvance(int nstep, NekDouble time) override;
117
118 void v_MountHOPBCs(int HBCdata, NekDouble kinvis, Array<OneD, NekDouble> &Q,
120
121 std::string v_GetSubStepName(void) override;
122
123 void SubStepAdvection(
124 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
125 Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble time);
126
128 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
129 Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble time);
130
134
136 const Array<OneD, NekDouble> &Vn,
137 const Array<OneD, const Array<OneD, NekDouble>> &physfield,
138 Array<OneD, Array<OneD, NekDouble>> &outarray);
139
141};
142} // namespace Nektar
143
144#endif
Binds a set of functions for use by time integration schemes.
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
An abstract base class encapsulating the concept of advection of a vector field.
Definition: Advection.h:81
~SubSteppingExtrapolate() override=default
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
void v_MountHOPBCs(int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection) override
LibUtilities::TimeIntegrationSchemeSharedPtr m_subStepIntegrationScheme
static std::string className
Name of class.
void v_SubStepSaveFields(int nstep) override
static ExtrapolateSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, MultiRegions::ExpListSharedPtr &pPressure, const Array< OneD, int > &pVel, const SolverUtils::AdvectionSharedPtr &advObject)
Creates an instance of this class.
void SubStepAdvection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void SubStepProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
void v_SubSteppingTimeIntegration(const LibUtilities::TimeIntegrationSchemeSharedPtr &IntegrationScheme) override
void v_SubStepAdvance(int nstep, NekDouble time) override
std::string v_GetSubStepName(void) override
Array< OneD, Array< OneD, NekDouble > > m_previousVnFields
void v_AccelerationBDF(Array< OneD, Array< OneD, NekDouble > > &array) override
void v_SubStepSetPressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_Dt, NekDouble kinvis) override
SubSteppingExtrapolate(const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
void AddAdvectionPenaltyFlux(const Array< OneD, NekDouble > &Vn, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
void SubStepExtrapolateField(NekDouble toff, Array< OneD, Array< OneD, NekDouble > > &ExtVel, Array< OneD, NekDouble > &ExtVn)
Array< OneD, Array< OneD, NekDouble > > m_previousVelFields
void v_EvaluatePressureBCs(const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis) override
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< TimeIntegrationScheme > TimeIntegrationSchemeSharedPtr
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< Advection > AdvectionSharedPtr
A shared pointer to an Advection object.
Definition: Advection.h:54
std::shared_ptr< SubSteppingExtrapolate > SubSteppingExtrapolateSharedPtr
std::shared_ptr< Extrapolate > ExtrapolateSharedPtr
Definition: Extrapolate.h:60
double NekDouble