Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
IsentropicVortex.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File IsentropicVortex.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: Euler equations for isentropic vortex
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
37 
38 using namespace std;
39 
40 namespace Nektar
41 {
42  string IsentropicVortex::className =
44  "IsentropicVortex", IsentropicVortex::create,
45  "Euler equations for the isentropic vortex test case.");
46 
47  IsentropicVortex::IsentropicVortex(
49  : EulerCFE(pSession)
50  {
51  }
52 
53  /**
54  * @brief Destructor for EulerCFE class.
55  */
57  {
58  }
59 
60  /**
61  * @brief Print out a summary with some relevant information.
62  */
64  {
66  SolverUtils::AddSummaryItem(s, "Problem Type", "Isentropic Vortex");
67  }
68 
69  /**
70  * @brief Set the initial conditions.
71  */
73  NekDouble initialtime,
74  bool dumpInitialConditions,
75  const int domain)
76  {
78 
79  int nTotQuadPoints = GetTotPoints();
80  Array<OneD, NekDouble> x(nTotQuadPoints);
81  Array<OneD, NekDouble> y(nTotQuadPoints);
82  Array<OneD, NekDouble> z(nTotQuadPoints);
84 
85  m_fields[0]->GetCoords(x, y, z);
86 
87  for (int i = 0; i < m_spacedim + 2; ++i)
88  {
89  u[i] = Array<OneD, NekDouble>(nTotQuadPoints);
90  }
91 
92  EvaluateIsentropicVortex(x, y, z, u, initialtime);
93 
94  // Forward transform to fill the coefficient space
95  for(int i = 0; i < m_fields.num_elements(); ++i)
96  {
97  Vmath::Vcopy(nTotQuadPoints, u[i], 1, m_fields[i]->UpdatePhys(), 1);
98  m_fields[i]->SetPhysState(true);
99  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
100  m_fields[i]->UpdateCoeffs());
101  }
102 
103  if (dumpInitialConditions && m_checksteps)
104  {
105  // Dump initial conditions to file
107  }
108  }
109 
110  /**
111  * @brief Get the exact solutions for isentropic vortex and Ringleb
112  * flow problems.
113  */
115  unsigned int field,
116  Array<OneD, NekDouble> &outfield,
117  const NekDouble time)
118  {
119  int nTotQuadPoints = GetTotPoints();
120  Array<OneD, NekDouble> x(nTotQuadPoints);
121  Array<OneD, NekDouble> y(nTotQuadPoints);
122  Array<OneD, NekDouble> z(nTotQuadPoints);
124 
125  m_fields[0]->GetCoords(x, y, z);
126 
127  for (int i = 0; i < m_spacedim + 2; ++i)
128  {
129  u[i] = Array<OneD, NekDouble>(nTotQuadPoints);
130  }
131 
132  EvaluateIsentropicVortex(x, y, z, u, time);
133 
134  Vmath::Vcopy(nTotQuadPoints, u[field], 1, outfield, 1);
135  }
136 
138  const Array<OneD, NekDouble> &x,
139  const Array<OneD, NekDouble> &y,
140  const Array<OneD, NekDouble> &z,
142  NekDouble time,
143  const int o)
144  {
145  int nq = x.num_elements();
146 
147  // Flow parameters
148  const NekDouble x0 = 5.0;
149  const NekDouble y0 = 0.0;
150  const NekDouble beta = 5.0;
151  const NekDouble u0 = 1.0;
152  const NekDouble v0 = 0.5;
153  const NekDouble gamma = m_gamma;
154  NekDouble r, xbar, ybar, tmp;
155  NekDouble fac = 1.0/(16.0*gamma*M_PI*M_PI);
156 
157  // In 3D zero rhow field.
158  if (m_spacedim == 3)
159  {
160  Vmath::Zero(nq, &u[3][o], 1);
161  }
162 
163  // Fill storage
164  for (int i = 0; i < nq; ++i)
165  {
166  xbar = x[i] - u0*time - x0;
167  ybar = y[i] - v0*time - y0;
168  r = sqrt(xbar*xbar + ybar*ybar);
169  tmp = beta*exp(1-r*r);
170  u[0][i+o] = pow(1.0 - (gamma-1.0)*tmp*tmp*fac, 1.0/(gamma-1.0));
171  u[1][i+o] = u[0][i+o]*(u0 - tmp*ybar/(2*M_PI));
172  u[2][i+o] = u[0][i+o]*(v0 + tmp*xbar/(2*M_PI));
173  u[m_spacedim+1][i+o] = pow(u[0][i+o], gamma)/(gamma-1.0) +
174  0.5*(u[1][i+o]*u[1][i+o] + u[2][i+o]*u[2][i+o]) / u[0][i+o];
175  }
176  }
177 
178 }
void EvaluateIsentropicVortex(const Array< OneD, NekDouble > &x, const Array< OneD, NekDouble > &y, const Array< OneD, NekDouble > &z, Array< OneD, Array< OneD, NekDouble > > &u, NekDouble time, const int o=0)
Isentropic Vortex Test Case.
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.
virtual ~IsentropicVortex()
Destructor for EulerCFE class.
STL namespace.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
int m_checksteps
Number of steps between checkpoints.
SOLVER_UTILS_EXPORT int GetTotPoints()
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s)
Print a summary of time stepping parameters.
virtual void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Set the initial conditions.
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
EquationSystemFactory & GetEquationSystemFactory()
virtual void v_GenerateSummary(SolverUtils::SummaryList &s)
Print a summary of time stepping parameters.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
virtual void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time=0.0)
Get the exact solutions for isentropic vortex and Ringleb flow problems.
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:373
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215