Nektar++
Functions
Nektar::RiemannTests Namespace Reference

Functions

 BOOST_AUTO_TEST_CASE (RoeAlongXconstSolution)
 
 BOOST_AUTO_TEST_CASE (RoeAlongYconstSolution)
 
 BOOST_AUTO_TEST_CASE (RoeAlongZconstSolution)
 
 BOOST_AUTO_TEST_CASE (RoeAlongXdensityJump)
 
 BOOST_AUTO_TEST_CASE (RRR)
 

Function Documentation

◆ BOOST_AUTO_TEST_CASE() [1/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongXconstSolution  )

Definition at line 46 of file TestRiemann.cpp.

47  {
48  // LibUtilities::SessionReaderSharedPtr dummySession;
49  // SolverUtils::RiemannSolverSharedPtr riemannSolver;
50  // std::string riemannName = "Roe";
51  // riemannSolver = SolverUtils::GetRiemannSolverFactory()
52  // .CreateInstance(riemName, m_session);
53  // auto riemannSolver = RoeSolver();
54  auto riemannSolver = RoeSolverSIMD();
55  // Setting up parameters for Riemann solver
56  NekDouble gamma = 1.4;
57  riemannSolver.SetParam("gamma", [&gamma]()
58  -> NekDouble& {return gamma;});
59 
60  size_t spaceDim = 3;
61  size_t npts = 5; // so that avx spillover loop is engaged
62 
63  // Set up locations of velocity vector.
64  Array<OneD, Array<OneD, NekDouble>> vecLocs(1);
65  vecLocs[0] = Array<OneD, NekDouble>(spaceDim);
66  for (size_t i = 0; i < spaceDim; ++i)
67  {
68  vecLocs[0][i] = 1 + i;
69  }
70  riemannSolver.SetAuxVec("vecLocs", [&vecLocs]()
71  -> const Array<OneD, const Array<OneD, NekDouble>>&
72  {return vecLocs;});
73 
74  // setup normals
75  Array<OneD, Array<OneD, NekDouble>> normals(spaceDim);
76  for (size_t i = 0; i < spaceDim; ++i)
77  {
78  normals[i] = Array<OneD, NekDouble>(npts);
79  }
80  riemannSolver.SetVector("N", [&normals]()
81  -> const Array<OneD, const Array<OneD, NekDouble>>&
82  {return normals;});
83 
84  size_t nFields = spaceDim + 2;
85  Array<OneD, Array<OneD, NekDouble>> fwd(nFields), bwd(nFields),
86  flx(nFields), flxRef(nFields);
87  for (size_t i = 0; i < nFields; ++i)
88  {
89  fwd[i] = Array<OneD, NekDouble>(npts);
90  bwd[i] = Array<OneD, NekDouble>(npts);
91  flx[i] = Array<OneD, NekDouble>(npts);
92  flxRef[i] = Array<OneD, NekDouble>(npts);
93  }
94 
95  // up to now it is "boiler plate" code to set up the test
96  // below are the conditions for the test
97 
98  for (size_t i = 0; i < npts; ++i)
99  {
100  // density
101  NekDouble rho = 0.9;
102  fwd[0][i] = rho;
103  bwd[0][i] = rho;
104  // x-momentum
105  NekDouble rhou = rho*1.0;
106  fwd[1][i] = rhou;
107  bwd[1][i] = rhou;
108  // y-momentum
109  NekDouble rhov = rho*2.0;
110  fwd[2][i] = rhov;
111  bwd[2][i] = rhov;
112  // z-momentum
113  NekDouble rhow = rho*3.0;
114  fwd[3][i] = rhow;
115  bwd[3][i] = rhow;
116  // energy
117  NekDouble p = 1.0;
118  NekDouble rhoe = p / (gamma - 1.0);
119  NekDouble E = rhoe + 0.5*(rhou*rhou + rhov*rhov + rhow*rhow) / rho;
120  fwd[nFields-1][i] = E;
121  bwd[nFields-1][i] = E;
122  // set face normal along x
123  normals[0][i] = 1.0;
124  // Ref solution
125  flxRef[0][i] = rhou;
126  flxRef[1][i] = rhou*rhou/rho + p;
127  flxRef[2][i] = rhou*rhov/rho;
128  flxRef[3][i] = rhou*rhow/rho;
129  flxRef[nFields-1][i] = (E+p)*rhou/rho;
130  }
131 
132  riemannSolver.Solve(spaceDim, fwd, bwd, flx);
133 
134  // check fluxes
135  for (size_t i = 0; i < npts; ++i)
136  {
137  BOOST_CHECK_CLOSE(flxRef[0][i], flx[0][i], 1e-10);
138  BOOST_CHECK_CLOSE(flxRef[1][i], flx[1][i], 1e-10);
139  BOOST_CHECK_CLOSE(flxRef[2][i], flx[2][i], 1e-10);
140  BOOST_CHECK_CLOSE(flxRef[3][i], flx[3][i], 1e-10);
141  BOOST_CHECK_CLOSE(flxRef[4][i], flx[4][i], 1e-10);
142  }
143 
144  }
double NekDouble

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [2/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongXdensityJump  )

Definition at line 336 of file TestRiemann.cpp.

337  {
338  // LibUtilities::SessionReaderSharedPtr dummySession;
339  // SolverUtils::RiemannSolverSharedPtr riemannSolver;
340  // std::string riemannName = "Roe";
341  // riemannSolver = SolverUtils::GetRiemannSolverFactory()
342  // .CreateInstance(riemName, m_session);
343  // auto riemannSolver = RoeSolver();
344  auto riemannSolver = RoeSolverSIMD();
345  // Setting up parameters for Riemann solver
346  NekDouble gamma = 1.4;
347  riemannSolver.SetParam("gamma", [&gamma]()
348  -> NekDouble& {return gamma;});
349 
350  size_t spaceDim = 3;
351  size_t npts = 11; // so that avx spillover loop is engaged
352 
353  // Set up locations of velocity vector.
354  Array<OneD, Array<OneD, NekDouble>> vecLocs(1);
355  vecLocs[0] = Array<OneD, NekDouble>(spaceDim);
356  for (size_t i = 0; i < spaceDim; ++i)
357  {
358  vecLocs[0][i] = 1 + i;
359  }
360  riemannSolver.SetAuxVec("vecLocs", [&vecLocs]()
361  -> const Array<OneD, const Array<OneD, NekDouble>>&
362  {return vecLocs;});
363 
364  // setup normals
365  Array<OneD, Array<OneD, NekDouble>> normals(spaceDim);
366  for (size_t i = 0; i < spaceDim; ++i)
367  {
368  normals[i] = Array<OneD, NekDouble>(npts);
369  }
370  riemannSolver.SetVector("N", [&normals]()
371  -> const Array<OneD, const Array<OneD, NekDouble>>&
372  {return normals;});
373 
374  size_t nFields = spaceDim + 2;
375  Array<OneD, Array<OneD, NekDouble>> fwd(nFields), bwd(nFields),
376  flx(nFields), flxRef(nFields);
377  for (size_t i = 0; i < nFields; ++i)
378  {
379  fwd[i] = Array<OneD, NekDouble>(npts);
380  bwd[i] = Array<OneD, NekDouble>(npts);
381  flx[i] = Array<OneD, NekDouble>(npts);
382  flxRef[i] = Array<OneD, NekDouble>(npts);
383  }
384 
385  // up to now it is "boiler plate" code to set up the test
386  // below are the conditions for the test
387 
388  for (size_t i = 0; i < npts; ++i)
389  {
390  // density
391  NekDouble rhoL = 1.0;
392  NekDouble rhoR = 2*rhoL;
393  fwd[0][i] = rhoL;
394  bwd[0][i] = rhoR;
395  // x-momentum
396  NekDouble rhou = rhoL*1.0;
397  fwd[1][i] = rhou;
398  bwd[1][i] = rhou;
399  // y-momentum
400  NekDouble rhov = rhoL*2.0;
401  fwd[2][i] = rhov;
402  bwd[2][i] = rhov;
403  // z-momentum
404  NekDouble rhow = rhoL*3.0;
405  fwd[3][i] = rhow;
406  bwd[3][i] = rhow;
407  // energy
408  NekDouble p = 1.0;
409  NekDouble rhoe = p / (gamma - 1.0);
410  NekDouble EL = rhoe + 0.5*(rhou*rhou + rhov*rhov + rhow*rhow) / rhoL;
411  NekDouble ER = rhoe + 0.5*(rhou*rhou + rhov*rhov + rhow*rhow) / rhoR;
412  fwd[nFields-1][i] = EL;
413  bwd[nFields-1][i] = ER;
414  // set face normal along x
415  normals[0][i] = 1.0;
416  // Ref solution (from point solve)
417  flxRef[0][i] = 0.87858599768171342;
418  flxRef[1][i] = 2.0449028304431223;
419  flxRef[2][i] = 1.8282946712594808;
420  flxRef[3][i] = 2.7424420068892208;
421  flxRef[nFields-1][i] = 9.8154698039903128;
422  }
423 
424  riemannSolver.Solve(spaceDim, fwd, bwd, flx);
425 
426  // check fluxes
427  for (size_t i = 0; i < npts; ++i)
428  {
429  BOOST_CHECK_CLOSE(flxRef[0][i], flx[0][i], 1e-10);
430  BOOST_CHECK_CLOSE(flxRef[1][i], flx[1][i], 1e-10);
431  BOOST_CHECK_CLOSE(flxRef[2][i], flx[2][i], 1e-10);
432  BOOST_CHECK_CLOSE(flxRef[3][i], flx[3][i], 1e-10);
433  BOOST_CHECK_CLOSE(flxRef[4][i], flx[4][i], 1e-10);
434  }
435 
436  }

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [3/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongYconstSolution  )

Definition at line 147 of file TestRiemann.cpp.

148  {
149  // LibUtilities::SessionReaderSharedPtr dummySession;
150  // SolverUtils::RiemannSolverSharedPtr riemannSolver;
151  // std::string riemannName = "Roe";
152  // riemannSolver = SolverUtils::GetRiemannSolverFactory()
153  // .CreateInstance(riemName, m_session);
154  // auto riemannSolver = RoeSolver();
155  auto riemannSolver = RoeSolverSIMD();
156  // Setting up parameters for Riemann solver
157  NekDouble gamma = 1.4;
158  riemannSolver.SetParam("gamma", [&gamma]()
159  -> NekDouble& {return gamma;});
160 
161  size_t spaceDim = 3;
162  size_t npts = 5;
163 
164  // Set up locations of velocity vector.
165  Array<OneD, Array<OneD, NekDouble>> vecLocs(1);
166  vecLocs[0] = Array<OneD, NekDouble>(spaceDim);
167  for (size_t i = 0; i < spaceDim; ++i)
168  {
169  vecLocs[0][i] = 1 + i;
170  }
171  riemannSolver.SetAuxVec("vecLocs", [&vecLocs]()
172  -> const Array<OneD, const Array<OneD, NekDouble>>&
173  {return vecLocs;});
174 
175  // setup normals
176  Array<OneD, Array<OneD, NekDouble>> normals(spaceDim);
177  for (size_t i = 0; i < spaceDim; ++i)
178  {
179  normals[i] = Array<OneD, NekDouble>(npts);
180  }
181  riemannSolver.SetVector("N", [&normals]()
182  -> const Array<OneD, const Array<OneD, NekDouble>>&
183  {return normals;});
184 
185  size_t nFields = spaceDim + 2;
186  Array<OneD, Array<OneD, NekDouble>> fwd(nFields), bwd(nFields),
187  flx(nFields), flxRef(nFields);
188  for (size_t i = 0; i < nFields; ++i)
189  {
190  fwd[i] = Array<OneD, NekDouble>(npts);
191  bwd[i] = Array<OneD, NekDouble>(npts);
192  flx[i] = Array<OneD, NekDouble>(npts);
193  flxRef[i] = Array<OneD, NekDouble>(npts);
194  }
195 
196  // up to now it is "boiler plate" code to set up the test
197  // below are the conditions for the test
198 
199  // density
200  NekDouble rho = 0.9;
201  fwd[0][0] = rho;
202  bwd[0][0] = rho;
203  // x-momentum
204  NekDouble rhou = rho*1.0;
205  fwd[1][0] = rhou;
206  bwd[1][0] = rhou;
207  // y-momentum
208  NekDouble rhov = rho*2.0;
209  fwd[2][0] = rhov;
210  bwd[2][0] = rhov;
211  // z-momentum
212  NekDouble rhow = rho*3.0;
213  fwd[3][0] = rhow;
214  bwd[3][0] = rhow;
215  // energy
216  NekDouble p = 1.0;
217  NekDouble rhoe = p / (gamma - 1.0);
218  NekDouble E = rhoe + 0.5*(rhou*rhou + rhov*rhov + rhow*rhow) / rho;
219  fwd[nFields-1][0] = E;
220  bwd[nFields-1][0] = E;
221  // set face normal along y
222  normals[1][0] = 1.0;
223  // Ref solution
224  flxRef[0][0] = rhov;
225  flxRef[1][0] = rhov*rhou/rho;
226  flxRef[2][0] = rhov*rhov/rho + p;
227  flxRef[3][0] = rhov*rhow/rho;
228  flxRef[nFields-1][0] = (E+p)*rhov/rho;
229 
230  riemannSolver.Solve(spaceDim, fwd, bwd, flx);
231 
232  // check fluxes
233  BOOST_CHECK_CLOSE(flxRef[0][0], flx[0][0], 1e-10);
234  BOOST_CHECK_CLOSE(flxRef[1][0], flx[1][0], 1e-10);
235  BOOST_CHECK_CLOSE(flxRef[2][0], flx[2][0], 1e-10);
236  BOOST_CHECK_CLOSE(flxRef[3][0], flx[3][0], 1e-10);
237  BOOST_CHECK_CLOSE(flxRef[4][0], flx[4][0], 1e-10);
238 
239  }

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [4/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongZconstSolution  )

Definition at line 241 of file TestRiemann.cpp.

242  {
243  // LibUtilities::SessionReaderSharedPtr dummySession;
244  // SolverUtils::RiemannSolverSharedPtr riemannSolver;
245  // std::string riemannName = "Roe";
246  // riemannSolver = SolverUtils::GetRiemannSolverFactory()
247  // .CreateInstance(riemName, m_session);
248  // auto riemannSolver = RoeSolver();
249  auto riemannSolver = RoeSolverSIMD();
250  // Setting up parameters for Riemann solver
251  NekDouble gamma = 1.4;
252  riemannSolver.SetParam("gamma", [&gamma]()
253  -> NekDouble& {return gamma;});
254 
255  size_t spaceDim = 3;
256  size_t npts = 1;
257 
258  // Set up locations of velocity vector.
259  Array<OneD, Array<OneD, NekDouble>> vecLocs(1);
260  vecLocs[0] = Array<OneD, NekDouble>(spaceDim);
261  for (size_t i = 0; i < spaceDim; ++i)
262  {
263  vecLocs[0][i] = 1 + i;
264  }
265  riemannSolver.SetAuxVec("vecLocs", [&vecLocs]()
266  -> const Array<OneD, const Array<OneD, NekDouble>>&
267  {return vecLocs;});
268 
269  // setup normals
270  Array<OneD, Array<OneD, NekDouble>> normals(spaceDim);
271  for (size_t i = 0; i < spaceDim; ++i)
272  {
273  normals[i] = Array<OneD, NekDouble>(npts);
274  }
275  riemannSolver.SetVector("N", [&normals]()
276  -> const Array<OneD, const Array<OneD, NekDouble>>&
277  {return normals;});
278 
279  size_t nFields = spaceDim + 2;
280  Array<OneD, Array<OneD, NekDouble>> fwd(nFields), bwd(nFields),
281  flx(nFields), flxRef(nFields);
282  for (size_t i = 0; i < nFields; ++i)
283  {
284  fwd[i] = Array<OneD, NekDouble>(npts);
285  bwd[i] = Array<OneD, NekDouble>(npts);
286  flx[i] = Array<OneD, NekDouble>(npts);
287  flxRef[i] = Array<OneD, NekDouble>(npts);
288  }
289 
290  // up to now it is "boiler plate" code to set up the test
291  // below are the conditions for the test
292 
293  // density
294  NekDouble rho = 0.9;
295  fwd[0][0] = rho;
296  bwd[0][0] = rho;
297  // x-momentum
298  NekDouble rhou = rho*1.0;
299  fwd[1][0] = rhou;
300  bwd[1][0] = rhou;
301  // y-momentum
302  NekDouble rhov = rho*2.0;
303  fwd[2][0] = rhov;
304  bwd[2][0] = rhov;
305  // z-momentum
306  NekDouble rhow = rho*3.0;
307  fwd[3][0] = rhow;
308  bwd[3][0] = rhow;
309  // energy
310  NekDouble p = 1.0;
311  NekDouble rhoe = p / (gamma - 1.0);
312  NekDouble E = rhoe + 0.5*(rhou*rhou + rhov*rhov + rhow*rhow) / rho;
313  fwd[nFields-1][0] = E;
314  bwd[nFields-1][0] = E;
315  // set face normal along y
316  normals[2][0] = 1.0;
317  // Ref solution
318  flxRef[0][0] = rhow;
319  flxRef[1][0] = rhow*rhou/rho;
320  flxRef[2][0] = rhow*rhov/rho;
321  flxRef[3][0] = rhow*rhow/rho + p;
322  flxRef[nFields-1][0] = (E+p)*rhow/rho;
323 
324  riemannSolver.Solve(spaceDim, fwd, bwd, flx);
325 
326  // check fluxes
327  BOOST_CHECK_CLOSE(flxRef[0][0], flx[0][0], 1e-10);
328  BOOST_CHECK_CLOSE(flxRef[1][0], flx[1][0], 1e-10);
329  BOOST_CHECK_CLOSE(flxRef[2][0], flx[2][0], 1e-10);
330  BOOST_CHECK_CLOSE(flxRef[3][0], flx[3][0], 1e-10);
331  BOOST_CHECK_CLOSE(flxRef[4][0], flx[4][0], 1e-10);
332 
333  }

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [5/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RRR  )

Definition at line 12 of file UnitTestRiemann.cpp.

13  {
14 
15  BOOST_CHECK_CLOSE(1., 1., 1e-10);
16  }