From 75ac548d0010c1bba08851690af7a458370edd68 Mon Sep 17 00:00:00 2001 From: Anchal Gupta Date: Wed, 27 Sep 2023 17:48:30 -0700 Subject: [PATCH 1/3] Added larger g file with unequal dims, fails test If we use sample/ITER_Lore_2296_00000/EQDSK/g002296.00200 then the geqdsk file has unequal dimensions for the r and z grid. For this test case, geqdsk_to_imas test fails. This shows that there might be a transpose issue in setting up of dd.equilibrium.time_slice[:].profiles_2d[:].psi Please verify this test and fix. --- test/runtests.jl | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index ff1818f..6f55c21 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -70,7 +70,7 @@ function define_default_sample_set() "/../sample/ITER_Lore_2296_00000/EQDSK/Baseline2008-li0.70.x4.mod2.eqdsk" return b2fgmtry, b2time, b2mn, gridspec, eqdsk end - +""" @testset "lightweight_utilities" begin # Gas unit converter flow_tls = 40.63 * Unitful.Torr * Unitful.L / Unitful.s @@ -252,14 +252,20 @@ end @test length(rho) > 10 @test maximum(rho) > 0 end +""" @testset "geqdsk_to_imas" begin sample_files = (splitdir(pathof(SD4SOLPS))[1] * "/../sample/") .* [ "g184833.03600", "geqdsk_iter_small_sample", ] + append!( + sample_files, + ["$(@__DIR__)/../sample/ITER_Lore_2296_00000/EQDSK/g002296.00200"], + ) tslice = 1 for sample_file ∈ sample_files + println(sample_file) dd = OMAS.dd() SD4SOLPS.geqdsk_to_imas(sample_file, dd; time_index=tslice) eqt = dd.equilibrium.time_slice[tslice] @@ -283,6 +289,7 @@ end p2 = eqt.profiles_2d[1] @test length(p2.grid.dim1) > 10 @test length(p2.grid.dim2) > 10 + println(size(p2.psi), (length(p2.grid.dim1), length(p2.grid.dim2))) @test size(p2.psi) == (length(p2.grid.dim1), length(p2.grid.dim2)) # derived From 349b6146125f9a4c66a69aefc591957d8ffa9820 Mon Sep 17 00:00:00 2001 From: David Eldon Date: Thu, 28 Sep 2023 08:26:28 -0700 Subject: [PATCH 2/3] Add additional tests to check equilibrium 2d profile dims --- test/runtests.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/test/runtests.jl b/test/runtests.jl index 6f55c21..45b2e40 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -288,6 +288,8 @@ end # 2d p2 = eqt.profiles_2d[1] @test length(p2.grid.dim1) > 10 + @test minimum(p2.grid.dim1) > 0 # R should be dim1 + @test minimum(p2.grid.dim2) < 0 @test length(p2.grid.dim2) > 10 println(size(p2.psi), (length(p2.grid.dim1), length(p2.grid.dim2))) @test size(p2.psi) == (length(p2.grid.dim1), length(p2.grid.dim2)) From 546f5cdd50a1b5ad4f3d253411e4920ae26b4aa3 Mon Sep 17 00:00:00 2001 From: David Eldon Date: Thu, 28 Sep 2023 09:08:19 -0700 Subject: [PATCH 3/3] Change transpose of psi in geqdsk_to_imas and add many tests --- src/SD4SOLPS.jl | 2 +- test/runtests.jl | 37 ++++++++++++++++++++++++++++--------- 2 files changed, 29 insertions(+), 10 deletions(-) diff --git a/src/SD4SOLPS.jl b/src/SD4SOLPS.jl index eed5d47..ebd9bfa 100644 --- a/src/SD4SOLPS.jl +++ b/src/SD4SOLPS.jl @@ -108,7 +108,7 @@ function geqdsk_to_imas(eqdsk_file, dd; time_index=1) p2 = eqt.profiles_2d[1] p2.grid.dim1 = collect(g.r) p2.grid.dim2 = collect(g.z) - p2.psi = Matrix(transpose(g.psirz)) # Not sure if transpose is correct + p2.psi = g.psirz # Not sure if transpose is correct # missing j_tor = pcurrt # Derived diff --git a/test/runtests.jl b/test/runtests.jl index 45b2e40..d26afcc 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -5,6 +5,7 @@ using EFIT: EFIT using Plots using Test using Unitful: Unitful +import Interpolations """ make_test_profile() @@ -285,21 +286,39 @@ end @test length(p1.pressure) == nprof @test length(p1.rho_tor_norm) == nprof + # boundary + r_bry = eqt.boundary.outline.r + z_bry = eqt.boundary.outline.z + @test length(r_bry) == length(z_bry) + # 2d + # Did the R and Z (dim1 and dim2) coordinates get written properly? p2 = eqt.profiles_2d[1] - @test length(p2.grid.dim1) > 10 - @test minimum(p2.grid.dim1) > 0 # R should be dim1 - @test minimum(p2.grid.dim2) < 0 - @test length(p2.grid.dim2) > 10 - println(size(p2.psi), (length(p2.grid.dim1), length(p2.grid.dim2))) - @test size(p2.psi) == (length(p2.grid.dim1), length(p2.grid.dim2)) + r_eq = p2.grid.dim1 + z_eq = p2.grid.dim2 + @test length(r_eq) > 10 + @test minimum(r_eq) > 0 # R should be dim1 + @test minimum(z_eq) < 0 + @test length(z_eq) > 10 + # Does psi match R and Z? + println(size(p2.psi), (length(r_eq), length(z_eq))) + @test size(p2.psi) == (length(r_eq), length(z_eq)) + # Does the psi grid look like it's transposed the right way? + # Many equilibrium reconstructions have equal numbers of cells in R and Z, + # so transposing the psi map incorrectly would not be detected by checking + # its array dimensions. It's also possible to incorrectly associate the psi + # map with R and Z (dim1 and dim2). So what recourse does that leave us? + # Well, the points on the boundary outline should be at psi_N = 1.0 to + # within some small tolerance. + psin2d = (p2.psi .- gq.psi_axis) ./ (gq.psi_boundary - gq.psi_axis) + tolerance = 2.0e-3 # It's not always a high res contour so cut some slack + psin_bry = Interpolations.LinearInterpolation((r_eq, z_eq), psin2d).(r_bry, z_bry) + @test maximum(psin_bry) < (1.0 + tolerance) + @test minimum(psin_bry) > (1.0 - tolerance) # derived @test gq.q_axis == p1.q[1] - # boundary - @test length(eqt.boundary.outline.r) == length(eqt.boundary.outline.z) - # wall limiter = dd.wall.description_2d[1].limiter @test length(limiter.unit[1].outline.r) > 10