cwpearson/astaroth
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Added a comment about helical forcing

Browse Source
This commit is contained in:
jpekkila
2019-08-06 19:47:03 +03:00
parent 6b53eb31ef
commit a6fca069a7
1 changed files with 7 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files

11
acc/mhd_solver/stencil_process.sps
View File

@@ -202,9 +202,6 @@ heat_transfer(in Vector uu, in Scalar lnrho, in Scalar tt)
#if LFORCING #if LFORCING
Vector Vector
simple_vortex_forcing(Vector a, Vector b, Scalar magnitude) simple_vortex_forcing(Vector a, Vector b, Scalar magnitude)
{ {
@@ -222,7 +219,13 @@ simple_outward_flow_forcing(Vector a, Vector b, Scalar magnitude)
Vector Vector
helical_forcing(Scalar magnitude, Vector k_force, Vector xx, Vector ff_re, Vector ff_im, Scalar phi) helical_forcing(Scalar magnitude, Vector k_force, Vector xx, Vector ff_re, Vector ff_im, Scalar phi)
{ {
// JP: This looks wrong:
// 1) Should it be dsx * nx instead of dsx * ny?
// 2) Should you also use globalGrid.n instead of the local n?
// 3) Also final point: can we do this with vectors/quaternions instead?
// Tringonometric functions are much more expensive and inaccurate/
// Fun related article:
// https://randomascii.wordpress.com/2014/10/09/intel-underestimates-error-bounds-by-1-3-quintillion/
xx.x = xx.x*(2.0*M_PI/(dsx*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min)))); xx.x = xx.x*(2.0*M_PI/(dsx*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min))));
xx.y = xx.y*(2.0*M_PI/(dsy*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min)))); xx.y = xx.y*(2.0*M_PI/(dsy*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min))));
xx.z = xx.z*(2.0*M_PI/(dsz*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min)))); xx.z = xx.z*(2.0*M_PI/(dsz*(DCONST_INT(AC_ny_max) - DCONST_INT(AC_ny_min))));