chalmersplasmatheory / DREAM

The Disruption Runaway Electron Analysis Model
MIT License
24 stars 5 forks source link

Bug found in RunawayFluid.cpp in partial derivative of Sauter conductivity #212

Open peterhalldestam opened 1 year ago

peterhalldestam commented 1 year ago

Bug in $\texttt{RunawayFluid::evaluatePartialContributionSauterConductivity}$.

When calculating the partial derivate of the Sauter conductivity $\sigma$ with respect to $Z_{\rm eff}$ (lines 809-811), there seem to be in the numerical derivative. As of now, it is something along the lines of:

$$ \frac{\partial\sigma}{\partial Z{\rm eff}} \approx \frac{ \sigma(Z{\rm eff}+h) - \sigma(Z_{\rm eff}+h) }{ h} , $$

which is always zero!

I'm unsure if it supposed to be a central difference or a forward difference. If the former, to minimise the number of required function evaluations (of $\sigma$), I would suggest using the already evaluated sigma ($\texttt{electricConductivity[ir]}$, that is) in the lines directly after, where the additional contribution from the $\ln{\Lambda}$ derivative is evaluated.

peterhalldestam commented 1 year ago

I noticed now that there are terms in the formula for $\sigma$ (as well as the bootstrap current in the Redl-Sauter paper) that include $\sqrt{Z\text{eff}-1}$. This would cause problems for the Newton solver when taking $\partial\sigma/\partial Z\text{eff}$ for whenever $Z\text{eff}=1$ (i.e. purely hydrogenic plasmas). Really, this will be a problem whenever $Z\text{eff}-h<1$ in the numerical derivative.

Therefore, I'd suggest using a forward difference here, i.e. $$\frac{\partial\sigma}{\partial Z\text{eff}} \approx\frac{\sigma(Z\text{eff}+h)-\sigma(Z_\text{eff})}{h}$$