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

Updated modelsolver to use the same induction equation as the updated DSL kernel

Browse Source
This commit is contained in:
jpekkila
2020-09-14 13:28:44 +03:00
parent 6d2be634d7
commit d56c7ef312
1 changed files with 7 additions and 5 deletions
Download Patch File Download Diff File Expand all files Collapse all files

12
src/utils/modelsolver.c
View File

@@ -702,17 +702,19 @@ momentum(const VectorData uu, const ScalarData lnrho
static inline Vector
induction(const VectorData uu, const VectorData aa)
{
Vector ind;
// Note: We do (-nabla^2 A + nabla(nabla dot A)) instead of (nabla x (nabla
// x A)) in order to avoid taking the first derivative twice (did the math,
// yes this actually works. See pg.28 in arXiv:astro-ph/0109497)
// u cross B - ETA * mu0 * (mu0^-1 * [- laplace A + grad div A ])
// u cross B - AC_eta * AC_mu0 * (AC_mu0^-1 * [- laplace A + grad div A ])
const Vector B = curl(aa);
const Vector grad_div = gradient_of_divergence(aa);
//MV: Due to gauge freedom we can reduce the gradient of scalar (divergence) from the equation
//const Vector grad_div = gradient_of_divergence(aa);
const Vector lap = laplace_vec(aa);
// Note, mu0 is cancelled out
ind = cross(vecvalue(uu), B) - getReal(AC_eta) * (grad_div - lap);
// Note, AC_mu0 is cancelled out
//MV: Due to gauge freedom we can reduce the gradient of scalar (divergence) from the equation
//const Vector ind = cross(value(uu), B) - getReal(AC_eta) * (grad_div - lap);
const Vector ind = cross(vecvalue(uu), B) + getReal(AC_eta) * lap;
return ind;
}