From dc037fcd29ad7c8635241bf9544f4010ffd5ae2f Mon Sep 17 00:00:00 2001 From: StuartBenjamin Date: Mon, 12 Oct 2026 16:39:15 -0700 Subject: [PATCH 1/2] Equilibrium - BUGFIX - Write gse.h5 and gsei.h5 when diagnose_src is set equilibrium_gse! referred to an undefined `flux` for the psi and theta grids, so diagnose_src = true stopped with an UndefVarError before writing gse.h5 and gsei.h5. Use equil.rzphi_xs / rzphi_ys. New testset checks the three files and their grids. Co-Authored-By: Claude Opus 5.5 --- src/Equilibrium/Equilibrium.jl | 6 +++--- test/runtests_equil.jl | 18 ++++++++++++++++++ 2 files changed, 21 insertions(+), 3 deletions(-) diff --git a/src/Equilibrium/Equilibrium.jl b/src/Equilibrium/Equilibrium.jl index 9bb0bd3ae..65a9307fe 100644 --- a/src/Equilibrium/Equilibrium.jl +++ b/src/Equilibrium/Equilibrium.jl @@ -625,8 +625,8 @@ function equilibrium_gse!(equil::PlasmaEquilibrium) gse_data = Array{Float32,3}(undef, mpsi + 1, mtheta + 1, 7) for ipsi in 0:mpsi for itheta in 0:mtheta - gse_data[ipsi+1, itheta+1, 1] = Float32(flux.ys[itheta+1]) - gse_data[ipsi+1, itheta+1, 2] = Float32(flux.xs[ipsi+1]) + gse_data[ipsi+1, itheta+1, 1] = Float32(equil.rzphi_ys[itheta+1]) + gse_data[ipsi+1, itheta+1, 2] = Float32(equil.rzphi_xs[ipsi+1]) gse_data[ipsi+1, itheta+1, 3] = Float32(flux_fs[ipsi+1, itheta+1, 1]) gse_data[ipsi+1, itheta+1, 4] = Float32(flux_fs[ipsi+1, itheta+1, 2]) gse_data[ipsi+1, itheta+1, 5] = Float32(source[ipsi+1, itheta+1]) @@ -639,7 +639,7 @@ function equilibrium_gse!(equil::PlasmaEquilibrium) # Write integrated error criterion h5open(joinpath(dirname(equil.config.eq_filename), "gsei.h5"), "w") do file - file["xs"] = Float32.(flux.xs) + file["xs"] = Float32.(equil.rzphi_xs) file["term"] = Float32.(term) file["totali"] = Float32.(totali) file["errori"] = Float32.(errori) diff --git a/test/runtests_equil.jl b/test/runtests_equil.jl index 2d33380c6..26b1249cf 100644 --- a/test/runtests_equil.jl +++ b/test/runtests_equil.jl @@ -1,3 +1,5 @@ +using HDF5 + @testset "Equilibrium Unit Tests" begin @@ -594,4 +596,20 @@ @test occursin("sol.retcode != ReturnCode.Success", src) @test occursin("isapprox(sol.t[end], 2π", src) end + + @testset "equilibrium_gse! writes gse.h5 with diagnose_src" begin + dir = mktempdir() + cp(joinpath(data_dir, "EQDSK_COCOS_02"), joinpath(dir, "EQDSK_COCOS_02")) + config = GeneralizedPerturbedEquilibrium.Equilibrium.EquilibriumConfig(; + eq_filename=joinpath(dir, "EQDSK_COCOS_02"), eq_type="efit", jac_type="boozer", + grid_type="ldp", psilow=0.01, psihigh=0.994) + eq = GeneralizedPerturbedEquilibrium.Equilibrium.setup_equilibrium(config) + eq.params.diagnose_src = true + GeneralizedPerturbedEquilibrium.Equilibrium.equilibrium_gse!(eq) + @test all(isfile(joinpath(dir, f)) for f in ("gsec.h5", "gse.h5", "gsei.h5")) + gse = HDF5.h5read(joinpath(dir, "gse.h5"), "gse_data") + @test gse[:, 1, 2] == Float32.(eq.rzphi_xs) + @test gse[1, :, 1] == Float32.(eq.rzphi_ys) + @test HDF5.h5read(joinpath(dir, "gsei.h5"), "xs") == Float32.(eq.rzphi_xs) + end end From 5dfa356478566904f5f709af6c6e6ece1be624d5 Mon Sep 17 00:00:00 2001 From: StuartBenjamin Date: Tue, 13 Oct 2026 01:30:24 -0700 Subject: [PATCH 2/2] =?UTF-8?q?Equilibrium=20-=20MINOR=20-=20Trim=20the=20?= =?UTF-8?q?gse.h5=20test=20to=20the=20file=20check=20and=20the=20=CF=88=20?= =?UTF-8?q?column?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-Authored-By: Claude Opus 5.5 --- test/runtests_equil.jl | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/test/runtests_equil.jl b/test/runtests_equil.jl index 26b1249cf..6468eaa18 100644 --- a/test/runtests_equil.jl +++ b/test/runtests_equil.jl @@ -1,6 +1,5 @@ using HDF5 - @testset "Equilibrium Unit Tests" begin # --- Directory Configuration --- @@ -607,9 +606,7 @@ using HDF5 eq.params.diagnose_src = true GeneralizedPerturbedEquilibrium.Equilibrium.equilibrium_gse!(eq) @test all(isfile(joinpath(dir, f)) for f in ("gsec.h5", "gse.h5", "gsei.h5")) - gse = HDF5.h5read(joinpath(dir, "gse.h5"), "gse_data") - @test gse[:, 1, 2] == Float32.(eq.rzphi_xs) - @test gse[1, :, 1] == Float32.(eq.rzphi_ys) - @test HDF5.h5read(joinpath(dir, "gsei.h5"), "xs") == Float32.(eq.rzphi_xs) + # Column 2 is ψ, as plot_gse_by_theta reads it + @test h5read(joinpath(dir, "gse.h5"), "gse_data")[:, 1, 2] == Float32.(eq.rzphi_xs) end end