Choice of ky definition in Numerics

[bpatel2107]

Currently we the definition of Numerics.ky is ambiguous. As of right now it is basically (effectively the choice made in GENE/GS2)

$k_y = \frac{n0}{L{ref}} \frac{\partial r_N}{\partial \psi_N}$

And we store this in pyro with units of [1 / rhoref], which technically seems fine and has dimensions of [1 / length], see below

The issue is that $\frac{\partial r_N}{\partial \psiN} = \frac{\partial r}{\partial \psi} B{ref} L_{ref}$ is normalised quantity with units of [Bref * Lref] so technically speaking ky should have units of [Bref * Lref / rhoref ] to ensure that $\frac{\partial r_N}{\partial \psi_N}$ transforms correctly, which feels very wrong but is technically right.

This all arises due ky using a dimensionless quantity which depends on the choice of normalisation.

This can cause issues specifically when wanting to transform to physical units as you wouldn't be able to convert something with units of [Bref * Lref / rhoref ] to $cm^{-1}$

We could avoid this by using the CGYRO defn where $k_y = \frac{n_0q}{r}$ or just bite the bullet and add in the "correct" dimensions in ky as is. We could also pivot entirely to using $n_0$ which is very safe but may not be user friendly as people tend to like working in ky rather than n.

Overall I am not really sure what to do so I would appreciate any thoughts on this.

[ZedThree]

What does IMAS do? I think it's probably best to have pyrokinetics store things in the "most sensible" manner, as long as we can convert to/from each code's idiosyncrasies.

The next question then is "what happens when we change definitions/units"? How would this affect existing work?

[bpatel2107]

IMAS store its own definition of $k\theta \rho{ref}$ which is basically the GKW input kthrho = $2\pi n0 \frac{\rho{ref}}{Lref} \sqrt{g^{\zeta\zeta}_{\theta=0}}$ which does not have this issue but is messy, but we do of course calculate this so its there...

In terms of previous runs there isn't zero impact as far as I can tell unless you use Bref != B0 or Bunit which is quite rare. Currently this is only possible in GS2 and before we hard coded in to force the use of B0. Using Bgeo has only been added recently supported (though I don't like it!)

A hacky fix now would be to sneak in the units of drho_dpsi when setting aky in the GS2 input file and leave as is. The problem then arises that if someone uses a non-standard B field, and then sets numerics.ky, the actual ky used in the simulation will differ slightly due to this sneak additional units.

[mrhardman]

In this discussion, how are B0 and Bgeo defined? I would guess Bgeo = I/Rgeo (and this is what I understood to be the usual GS2 normalisation).

[bpatel2107]

When we work with normalised units without physical reference values, it is generally important to capture how we transform from one unit to another i.e. knowing that $ cs = \sqrt{2} v{th}$ and then pint lets you convert from one to another unit to another

https://github.com/pyro-kinetics/pyrokinetics/blob/f085c348d7e5f1ff6206b5aebdd58496b76d292f/src/pyrokinetics/units.py#L169-L170

So when loading in the this file the only thing we need to know is how to map from B0 to Bgeo which we do via rmaj/rgeo

https://github.com/pyro-kinetics/pyrokinetics/blob/f085c348d7e5f1ff6206b5aebdd58496b76d292f/src/pyrokinetics/normalisation.py#L415-L416

Technically when we do have reference values B0 is defined as $B0 = I(\psi)/R{major}$ which is as you say and similarly. $B{geo} = I(\psi)/R{major}$.

The larger issue here is the units we give to numerics.ky which I think should be in units of [1 / rhoref] but due to the definition of $k_y = \frac{n0}{L{ref}} \frac{\partial r_{N}}{\partial \psi_N}$ means that we need to be much more careful when transforming.

Potentially the easiest solution would be to say that we force $k_y= \frac{n0}{L{ref}} \frac{\partial r{pyro}}{\partial \psi{pyro}}$ and then handle the extra units of dpsidr when reading/writing GS2/STELLA input files which is relatively straightforward to do.

[d7919]

I suppose you could say aky_gs2 = ky_gs2 rho_ref = ky_pyro rho_ref G where G is the geometrical factor? This allows ky_pyro to have whatever definition and normalisation you are most happy with?

[bpatel2107]

Yes that is effectively what I am suggesting to an extent. Will have a go it at some point this week

[bpatel2107]

So I think the solution is to have a technical definition of numerics.ky as $k_y = n_0 \frac{\partial r}{\partial \psi} B0$, and has units of $\rho{ref}$. Note all of these terms are physical units so is independent of choice of normalisation.

To map to GS2/STELLA ky which is $\frac{n0}{L{ref}}\frac{\partial\rho}{\partial \psi_N} = n0 \frac{\partial r}{\partial \psi} B{ref}$ we need a factor of $\frac{B_{ref}}{B_0}$.

To map to CGYRO/TGLF ky which is $n0 \frac{q}{r}$, using $B{unit} = \frac{q}{r} \frac{\partial r}{\partial \psi}$ we get $n0 \frac{\partial r}{\partial \psi} B{unit}$ we need a factor of $\frac{B_{unit}}{B_0}$

To map to GENE ky which is $\frac{n0}{L{ref}} C_y$ where $Cy = \frac{\partial r}{\partial \psi} B{ref}$ (same as GS2), we need a factor of $\frac{B_{ref}}{B_0}$

To map GKW is a whole thing but essentially we just need to throw in another factor of $\frac{B_{ref}}{B_0}$

So far only GS2/STELLA support Bref != B0 (CGYRO/TGLF always use Bunit and GENE/GKW always used B0) so they are the only code who can be impacted by this.

I've modified #324 to handle this now