boutproject / BOUT-dev

BOUT++: Plasma fluid finite-difference simulation code in curvilinear coordinate systems
http://boutproject.github.io/
GNU Lesser General Public License v3.0
182 stars 95 forks source link

Handling of negative poloidal field, left-hand coordinates #1617

Open bendudson opened 5 years ago

bendudson commented 5 years ago

The coordinate system and operators were mostly derived while assuming that the coordinate system is right handed, so the Jacobian J is positive. Sorting this out becomes important in tokamak simulations with reversed current, where Bp < 0.

In the definition of the coordinate system:https://bout-dev.readthedocs.io/en/latest/user_docs/coordinates.html#id1 there is a sigma_Btheta variable, containing the sign of the poloidal field. This needs to be carefully propagated through the calculation. If the Jacobian changes sign, since it is hthe / Bp, then it matters whether sqrt(g_22) or J * Bxy are used. Currently sqrt(g_22) is used in the Grad_par operators, for example: https://github.com/boutproject/BOUT-dev/blob/master/src/mesh/coordinates.cxx#L670

Changing sign of the Z -> toroidal angle mapping may also affect the ShiftZ and TwistShift settings, probably requiring changes in Hypnotoad.

A note for Hypnotoad: most EFIT input files have poloidal flux increasing with minor radius, regardless of which direction the poloidal field is in. I think the sign is supposed to be inferred from the sign of the plasma current.

Thanks to Haruki Seto-san and Ben Zhu for work on this.

friva000 commented 5 years ago

Would a possible solution be to change also the sign of y, such that the (x,y,z) coordinate system is always right handed (i.e., y=sigma_Btheta*theta)? This would ensure that J is always positive and that b=e_y/\sqrt{g_yy}.

mikekryjak commented 1 year ago

@bendudson This has come up again in my meshing of ST-40. What workaround do you recommend? I have talked to @johnomotani and he suggested I could use the reverse_current option in Hypnotoad, and that this may be OK even though the signs of Bp and Bt will be incorrect if we are running without drifts.

bendudson commented 1 year ago

Hey @mikekryjak ,

Using reverse_current is probably a good workaround for now, at least as a first pass. If the plasma is up-down symmetric then it should look similar if the current is reversed and the machine is turned upside down. There is a paper somewhere on symmetries in MHD that would be useful here..

The correct thing to do would be to handle both current directions. John and I went through the coordinate system calculation (again) and I'm 99% sure this is now correct: https://bout-dev.readthedocs.io/en/latest/user_docs/coordinates.html#right-handed-field-aligned-coordinates . The calculations are however very fiddly and it's very easy to accidentally miss out a minus sign. The other part of this is ensuring that all the differential operators in the code correctly handle both left and right-handed coordinate systems.

mikekryjak commented 1 year ago

Thanks for the summary @bendudson. Any idea what pitfalls to expect when using reverse_current apart from issues when drifts are enabled?