fusion-energy / openmc-plasma-source

Creates a plasma source as an openmc.source object from input parameters that describe the plasma
MIT License
26 stars 11 forks source link

extra tests ideas #34

Open shimwell opened 2 years ago

shimwell commented 2 years ago

Test that point coordinates are generated within the ranges expected.

This is tricky for random samples as the results change

Allow for a wide tolerance to account for random distributions

Perhaps fixing the seed in openmc could help, openmc has a settings.seed value that defaults to 1

shimwell commented 2 years ago

testing the plotting routines ideas

check that plots result in file such as png files

check the axis ranges are correct

LiamPattinson commented 2 years ago

Hi Jonathan, for the first point, would I be right in thinking that you're interested in testing that TokamakSource is generating RZ coordinates correctly?

shimwell commented 2 years ago

That would be very nice, something that checks triangulation and elongation are correctly handled. What do you think @RemDelaporteMathurin ?

RemDelaporteMathurin commented 2 years ago

Hi @shimwell @LiamPattinson

So I think there are several things here and maybe we need seperated issues to make it clearer.

The TokamakSource class has a method called sample_sources() which creates a bunch of random (R, Z) coordinates based on the plasma geometry parameters (radii, elongation, triangularity, etc.). I think that's what @shimwell is refering to. https://github.com/fusion-energy/openmc-plasma-source/blob/3de9d21750b5e30ce3e38d5f8478f764cda48b8c/openmc_plasma_source/tokamak_source.py#L198-L215

One way to test this is to 1) test the method convert_a_alpha_to_R_Z() 2) for several sets of parameters (radii, elongation, etc.) test that the sampled coordinates are all in the shell produced by the plasma boundary equations : R = major_radius + a cos(alpha + triangularity sin(alpha) ) + shafranov_shift Z = elongation a sin(alpha)

https://github.com/fusion-energy/openmc-plasma-source/blob/3de9d21750b5e30ce3e38d5f8478f764cda48b8c/openmc_plasma_source/tokamak_source.py#L178

LiamPattinson commented 2 years ago

Thanks @RemDelaporteMathurin and @shimwell, that's what I assumed the issue was about. I've managed to get such a test working, though with a few caveats:

The first point is difficult to solve, as the plasma boundary is quite a complex shape to work with, but I might be able to constrain the test further with a bit more work. The second point should be quite easy to incorporate into my solution, as I'll just need to add a few more examples of tokamak geometries. Once I've done this I'll make a pull request.

RemDelaporteMathurin commented 2 years ago

Thanks for this. I guess that works but doesn't test the triangularity.

Perhaps we can focus on testing convert_a_alpha_to_R_Z. An "easy" way to do it is to create an inverse function convert_R_Z_to_a_alpha and check that convert_a_alpha_to_R_Z(convert_R_Z_to_a_alpha(R_test, Z_test)) == R_test, Z_test

So to create that inverse function, "simply" inverse the relations I gave above:

R = major_radius + a cos(alpha + triangularity sin(alpha) ) + shafranov_shift Z = elongation a sin(alpha)

Here's an idea of what it should look like: image

LiamPattinson commented 2 years ago

I've tried looking at an inverse function, but I believe it can only ever be approximate as the a-alpha to RZ relations aren't bijective (at least according to "Tokamak D-T neutron source models for different plasma physics confinement modes", Fausser et al.). I believe I have an approximate solution working, but it'll take some further testing to make sure it's robust.

RemDelaporteMathurin commented 2 years ago

It isn't bijective as is but if you constraint a to be positive than each R, Z couple corresponds to one unique (a, alpha) point doesn't it?

LiamPattinson commented 2 years ago

I believe so, but the inversion may need to be performed numerically rather than analytically.

RemDelaporteMathurin commented 2 years ago

Ah maybe. Well even an approximate is better than what we have currently I suppose :-) Thanks!