Twistshift for field-aligned variables

[johnomotani]

When TwistShift = true, apply the shift only when communicating field-aligned variables. The twist-shift at the branch cut in the core region is needed for field-aligned variables, but incorrect for variables in x-z orthogonal coordinates.

When TwistShift = true it is necessary to load ShiftAngle from the gridfile or grid options, so it is useful to keep the option. FCI simulations and simulations in periodic slab geometry do not require a twist-shift, and so should be able to run without providing ShiftAngle.

I'm in favour of changing the default value of TwistShift from false to true. When TwistShift = true we will now either do the right thing, or raise an exception which tells the user how to fix the error (set TwistShift = false) and when it is safe to do so (for FCI, periodic slab or slab with only open field lines). If TwistShift = false there is a danger of a simulation running but giving incorrect results - then the best case is that the errors cause an instability that makes them obvious in the way that Fabio found. On the other hand changing the default may break existing simulations (e.g. many of the tests need to set TwistShift = false).

Fixes #1731.

[johnomotani]

These changes do raise the question of whether simulations using ParallelTransformIdentity should be run with field-aligned variables so that they do all get the twist-shift applied. I think this is not a common use case, but it would be nice to decide...

It might help to do that if the parallel derivative operators, etc. returned a field-aligned field when their original input was field-aligned. Looking at how the field-aligned variables are going in STORM, I think this could be a useful feature. While it is sometimes shorter to have the result of say Grad_par(f_aligned) be returned in the x-z orthogonal coordinates, it would sometimes be nice to get the result back without the transform from field-aligned to orthogonal; there's no performance penalty to doing fromFieldAligned(Grad_par(f_aligned)) and the checks will catch any places where the fromFieldAligned is missed by mistake.

Any thoughts?

[friva000]

Personally, I would prefer keeping the default TwistShift = false. Then, I would use an if in the communicate function: if the field is not field aligned, just communicate, if it is is field aligned and TwistShift = false, throw, if it is field aligned and TwistShift = true, communicate and apply the twist shift.

For the result of the operators in field aligned or in orthogonal coordinates: is it possible to introduce a flag in the operators, with default value returning orthogonal, but if changed not performing the fromFieldAligned and just returning the result in field aligned coordinates?

[johnomotani]

@friva000 yes, I like that! Will update...

[johnomotani]

I've rolled back to before changing the default of TwistShift.

After some discussion with Fabio about how the different ParallelTransforms should behave, I've come to the design I've just pushed: add a new twistShift(bool twist_shift_enabled) method of Field3D,
that then calls a twistShift method of ParallelTransform with both the value of the TwistShift flag, and the y-direction type (standard or aligned) of the Field3D. This enables us to have different behaviour depending on the ParallelTransform:

• ParallelTransformIdentity twist-shifts all fields if TwistShift=true, or none if TwistShift=false (this keeps the old behaviour of TwistShift).
• ShiftedMetric twist-shifts only field-aligned fields, throwing an exception if TwistShift=false and the field is field-aligned.
• FCI never twist-shifts as there cannot be field-aligned variables.

[ZedThree]

@friva000 @johnomotani Is there any chance there's some kind of simple test we can write to check we're handling twist-shifting correctly? For example, checking that something constant on a field-line, but with some shape perpendicular to it, is really periodic? Ideally as a unit test, but even as an integrated test would be good.

[johnomotani]

Testing the twist-shift operation requires a fully-initialised BoutMesh. I think that's hard to implement in a unit test, so I've made an integrated test.

The test uses ShiftedMetric so it has toFieldAligned/fromFieldAligned functions. It would be nice to have a test using ParallelTransformIdentity, especially as ShiftAngle is used to initialise zShift, which is used by to/fromFieldAligned. That test would need to construct field-aligned versions of poloidally-periodic functions for itself though, so would be more complicated to write.