[Inverse solver] Issue creating single null equilibria

[kpentland]

Prior issues (now closed) raised a problem in the inverse solver that needs to be investigated around generating single null equilibria (think JET/ITER-like plasmas).

Following some playing around with the "Anamak" tokamak example (#61), I tried to generate a single null plasma:

While the inverse solve appears to have worked, the solver itself is suggesting non-convergence:

This should be investigated and is probably related to the fact that when I pass this "converged looking" but "unconverged really" equilibrium into a static forward solve, the forward solve then returns the following:

My example can be found here:

example07a - Anamak_single_null.ipynb

[ubaidali]

I've played around with it a bit and found the following:

• Setting max_iter_per_update=10/20 produces a very divergent solution,
• Setting force_up_down_symmetric=False produces a very divergent solution,
• Setting Picard_handover=1e-4 with force_up_down_symmetric=False produces a much better solution - still doesn't show convergence in the inverse solve but the subsequent forward solve produces a very similar solution.

This then led me to work out that the the NK shouldn't be used in the initial stages and that the cost function has a local minimum that the linear update is unable to get out of, so i then tried the following settings and got a converged solution...

GSStaticSolver.inverse_solve(eq=eq,
profiles=profiles,
constrain=constrain,
target_relative_tolerance=1e-5,
target_relative_psit_update=1e-3,
verbose=True,
l2_reg=np.array([1e-12]*N_actives),
# max_iter_per_update=20,
Picard_handover=1e-4,
full_jacobian_handover=[2e-2,1e-2],
force_up_down_symmetric=False,
)

[kpentland]

This is great, thanks @ubaidali! We should dig into this some more...

There are still times when the constrains look like they're almost there but no quite and yet the solver still says the solution "converged". I'll have a think...

[ubaidali]

I think what it comes down to is that formally the solution is not differentiable (when the problem is written in the strong form in finite difference) when the plasma has a x-point. So Newton or Newton-Krylov will not show full quadratic convergence because the solution is not locally quadratic and small changes in current lead to topological changes.... in the limited or double-null case the discontinuity in phase space is less severe so NK is able to converge closer and closer... this same issue might also come up in virtual circuits - not sure if you have seen inaccuracies/lack-of-consistency in the FD gradients in diverted plasmas...