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",
58  [&gamma]() -> 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(
71  "vecLocs",
72  [&vecLocs]() -> const Array<OneD, const Array<OneD, NekDouble>> & {
73  return vecLocs;
74  });
75 
76  // setup normals
77  Array<OneD, Array<OneD, NekDouble>> normals(spaceDim);
78  for (size_t i = 0; i < spaceDim; ++i)
79  {
80  normals[i] = Array<OneD, NekDouble>(npts);
81  }
82  riemannSolver.SetVector(
83  "N", [&normals]() -> const Array<OneD, const Array<OneD, NekDouble>> & {
84  return normals;
85  });
86 
87  size_t nFields = spaceDim + 2;
88  Array<OneD, Array<OneD, NekDouble>> fwd(nFields), bwd(nFields),
89  flx(nFields), flxRef(nFields);
90  for (size_t i = 0; i < nFields; ++i)
91  {
92  fwd[i] = Array<OneD, NekDouble>(npts);
93  bwd[i] = Array<OneD, NekDouble>(npts);
94  flx[i] = Array<OneD, NekDouble>(npts);
95  flxRef[i] = Array<OneD, NekDouble>(npts);
96  }
97 
98  // up to now it is "boiler plate" code to set up the test
99  // below are the conditions for the test
100 
101  for (size_t i = 0; i < npts; ++i)
102  {
103  // density
104  NekDouble rho = 0.9;
105  fwd[0][i] = rho;
106  bwd[0][i] = rho;
107  // x-momentum
108  NekDouble rhou = rho * 1.0;
109  fwd[1][i] = rhou;
110  bwd[1][i] = rhou;
111  // y-momentum
112  NekDouble rhov = rho * 2.0;
113  fwd[2][i] = rhov;
114  bwd[2][i] = rhov;
115  // z-momentum
116  NekDouble rhow = rho * 3.0;
117  fwd[3][i] = rhow;
118  bwd[3][i] = rhow;
119  // energy
120  NekDouble p = 1.0;
121  NekDouble rhoe = p / (gamma - 1.0);
122  NekDouble E =
123  rhoe + 0.5 * (rhou * rhou + rhov * rhov + rhow * rhow) / rho;
124  fwd[nFields - 1][i] = E;
125  bwd[nFields - 1][i] = E;
126  // set face normal along x
127  normals[0][i] = 1.0;
128  // Ref solution
129  flxRef[0][i] = rhou;
130  flxRef[1][i] = rhou * rhou / rho + p;
131  flxRef[2][i] = rhou * rhov / rho;
132  flxRef[3][i] = rhou * rhow / rho;
133  flxRef[nFields - 1][i] = (E + p) * rhou / rho;
134  }
135 
136  riemannSolver.Solve(spaceDim, fwd, bwd, flx);
137 
138  // check fluxes
139  for (size_t i = 0; i < npts; ++i)
140  {
141  BOOST_CHECK_CLOSE(flxRef[0][i], flx[0][i], 1e-10);
142  BOOST_CHECK_CLOSE(flxRef[1][i], flx[1][i], 1e-10);
143  BOOST_CHECK_CLOSE(flxRef[2][i], flx[2][i], 1e-10);
144  BOOST_CHECK_CLOSE(flxRef[3][i], flx[3][i], 1e-10);
145  BOOST_CHECK_CLOSE(flxRef[4][i], flx[4][i], 1e-10);
146  }
147 }
double NekDouble

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [2/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongXdensityJump  )

Definition at line 341 of file TestRiemann.cpp.

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

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [3/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongYconstSolution  )

Definition at line 149 of file TestRiemann.cpp.

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

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [4/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RoeAlongZconstSolution  )

Definition at line 245 of file TestRiemann.cpp.

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

References CellMLToNektar.cellml_metadata::p.

◆ BOOST_AUTO_TEST_CASE() [5/5]

Nektar::RiemannTests::BOOST_AUTO_TEST_CASE ( RRR  )

Definition at line 43 of file UnitTestRiemann.cpp.

44 {
45 
46  BOOST_CHECK_CLOSE(1., 1., 1e-10);
47 }