Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NavierStokesCFE.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File NavierStokesCFE.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: Navier Stokes equations in conservative variables
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
37 
38 using namespace std;
39 
40 namespace Nektar
41 {
42  string NavierStokesCFE::className =
44  "NavierStokesCFE", NavierStokesCFE::create,
45  "NavierStokes equations in conservative variables.");
46 
47  NavierStokesCFE::NavierStokesCFE(
49  : CompressibleFlowSystem(pSession)
50  {
51  }
52 
54  {
56 
57  if(m_session->DefinesSolverInfo("PROBLEMTYPE"))
58  {
59 
60  std::string ProblemTypeStr = m_session->GetSolverInfo("PROBLEMTYPE");
61  for(int i = 0; i < (int) SIZE_ProblemType; ++i)
62  {
63  if(NoCaseStringCompare(ProblemTypeMap[i],ProblemTypeStr) == 0)
64  {
66  break;
67  }
68  }
69  }
70  else
71  {
73  }
74 
76  {
79  }
80  else
81  {
82  ASSERTL0(false, "Implicit CFE not set up.");
83  }
84  }
85 
87  {
88 
89  }
90 
92  {
94  SolverUtils::AddSummaryItem(s, "Problem Type",
96  }
97 
99  NekDouble initialtime,
100  bool dumpInitialConditions,
101  const int domain)
102  {
103  EquationSystem::v_SetInitialConditions(initialtime, false);
104 
105  // insert white noise in initial condition
106  NekDouble Noise;
107  int phystot = m_fields[0]->GetTotPoints();
108  Array<OneD, NekDouble> noise(phystot);
109 
110  m_session->LoadParameter("Noise", Noise,0.0);
111  int m_nConvectiveFields = m_fields.num_elements();
112 
113  if (Noise > 0.0)
114  {
115  for (int i = 0; i < m_nConvectiveFields; i++)
116  {
117  Vmath::FillWhiteNoise(phystot, Noise, noise, 1,
118  m_comm->GetColumnComm()->GetRank()+1);
119  Vmath::Vadd(phystot, m_fields[i]->GetPhys(), 1,
120  noise, 1, m_fields[i]->UpdatePhys(), 1);
121  m_fields[i]->FwdTrans_IterPerExp(m_fields[i]->GetPhys(),
122  m_fields[i]->UpdateCoeffs());
123  }
124  }
125 
127 
128  if (dumpInitialConditions)
129  {
130  // Dump initial conditions to file
132  }
133  }
134 
136  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
137  Array<OneD, Array<OneD, NekDouble> > &outarray,
138  const NekDouble time)
139  {
140  int i;
141  int nvariables = inarray.num_elements();
142  int npoints = GetNpoints();
143 
144 
146  Array<OneD, Array<OneD, NekDouble> > outarrayAdv(nvariables);
147  Array<OneD, Array<OneD, NekDouble> > outarrayDiff(nvariables);
148 
149  Array<OneD, Array<OneD, NekDouble> > inarrayTemp(nvariables-1);
150  Array<OneD, Array<OneD, NekDouble> > inarrayDiff(nvariables-1);
151 
152  for (i = 0; i < nvariables; ++i)
153  {
154  outarrayAdv[i] = Array<OneD, NekDouble>(npoints, 0.0);
155  outarrayDiff[i] = Array<OneD, NekDouble>(npoints, 0.0);
156  }
157 
158  for (i = 0; i < nvariables-1; ++i)
159  {
160  inarrayTemp[i] = Array<OneD, NekDouble>(npoints, 0.0);
161  inarrayDiff[i] = Array<OneD, NekDouble>(npoints, 0.0);
162  }
163 
164  // Advection term in physical rhs form
165  m_advection->Advect(nvariables, m_fields, advVel, inarray,
166  outarrayAdv, time);
167 
168  // Extract pressure and temperature
169  Array<OneD, NekDouble > pressure (npoints, 0.0);
170  Array<OneD, NekDouble > temperature(npoints, 0.0);
171  GetPressure(inarray, pressure);
172  GetTemperature(inarray, pressure, temperature);
173 
174  // Extract velocities
175  for (i = 1; i < nvariables-1; ++i)
176  {
177  Vmath::Vdiv(npoints,
178  inarray[i], 1,
179  inarray[0], 1,
180  inarrayTemp[i-1], 1);
181  }
182 
183  // Copy velocities into new inarrayDiff
184  for (i = 0; i < nvariables-2; ++i)
185  {
186  Vmath::Vcopy(npoints, inarrayTemp[i], 1, inarrayDiff[i], 1);
187  }
188 
189  // Copy temperature into new inarrayDiffusion
190  Vmath::Vcopy(npoints,
191  temperature, 1,
192  inarrayDiff[nvariables-2], 1);
193 
194  // Diffusion term in physical rhs form
195  m_diffusion->Diffuse(nvariables, m_fields, inarrayDiff, outarrayDiff);
196 
197  for (i = 0; i < nvariables; ++i)
198  {
199  Vmath::Vsub(npoints,
200  outarrayDiff[i], 1,
201  outarrayAdv[i], 1,
202  outarray[i], 1);
203  }
204 
205  // Add sponge layer if defined in the session file
206  std::vector<SolverUtils::ForcingSharedPtr>::const_iterator x;
207  for (x = m_forcing.begin(); x != m_forcing.end(); ++x)
208  {
209  (*x)->Apply(m_fields, inarray, outarray, time);
210  }
211  }
212 
214  const Array<OneD, const Array<OneD, NekDouble> > &inarray,
215  Array<OneD, Array<OneD, NekDouble> > &outarray,
216  const NekDouble time)
217  {
218  int i;
219  int nvariables = inarray.num_elements();
220 
221  switch(m_projectionType)
222  {
224  {
225  // Just copy over array
226  int npoints = GetNpoints();
227 
228  for(i = 0; i < nvariables; ++i)
229  {
230  Vmath::Vcopy(npoints, inarray[i], 1, outarray[i], 1);
231  }
232  SetBoundaryConditions(outarray, time);
233  break;
234  }
237  {
238  ASSERTL0(false, "No Continuous Galerkin for full compressible "
239  "Navier-Stokes equations");
240  break;
241  }
242  default:
243  ASSERTL0(false, "Unknown projection scheme");
244  break;
245  }
246  }
247 
249  Array<OneD, Array<OneD, NekDouble> > &inarray,
250  NekDouble time)
251  {
252  std::string varName;
253  int cnt = 0;
254 
255  // loop over Boundary Regions
256  for (int n = 0; n < m_fields[0]->GetBndConditions().num_elements(); ++n)
257  {
258  std::string type = m_fields[0]->GetBndConditions()[n]->GetUserDefined();
259  SetCommonBC(type,n,time,cnt,inarray);
260  }
261  }
262 }
virtual void v_InitObject()
Initialization object for CompressibleFlowSystem class.
void DoOdeProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
virtual void v_InitObject()
Initialization object for CompressibleFlowSystem class.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
SOLVER_UTILS_EXPORT int NoCaseStringCompare(const std::string &s1, const std::string &s2)
Perform a case-insensitive string comparison.
std::vector< std::pair< std::string, std::string > > SummaryList
Definition: Misc.h:47
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
virtual void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
STL namespace.
void GetTemperature(const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure, Array< OneD, NekDouble > &temperature)
void Vdiv(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.cpp:227
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
LibUtilities::CommSharedPtr m_comm
Communicator.
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
virtual void v_GenerateSummary(SolverUtils::SummaryList &s)
Print a summary of time stepping parameters.
void SetBoundaryConditions(Array< OneD, Array< OneD, NekDouble > > &physarray, NekDouble time)
void SetCommonBC(const std::string &userDefStr, const int n, const NekDouble time, int &cnt, Array< OneD, Array< OneD, NekDouble > > &inarray)
Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem ...
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:50
int m_spacedim
Spatial dimension (>= expansion dim).
double NekDouble
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Length of enum list.
Definition: EulerADCFE.h:47
EquationSystemFactory & GetEquationSystemFactory()
SolverUtils::AdvectionSharedPtr m_advection
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.cpp:329
SolverUtils::DiffusionSharedPtr m_diffusion
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
virtual void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
void FillWhiteNoise(int n, const T eps, T *x, const int incx, int outseed)
Fills a vector with white noise.
Definition: Vmath.cpp:138
SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr GetPressure()
Get pressure field if available.
ProblemType
Definition: EulerADCFE.h:44
virtual void v_GenerateSummary(SolverUtils::SummaryList &s)
Print a summary of time stepping parameters.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
std::vector< SolverUtils::ForcingSharedPtr > m_forcing
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.cpp:285
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215
const char *const ProblemTypeMap[]
Definition: EulerADCFE.h:50