From 752d6ed469f9ab69759d3af409a62d6c6e73b0c4 Mon Sep 17 00:00:00 2001 From: Orso Meneghini Date: Sun, 17 Dec 2023 17:43:14 -0800 Subject: [PATCH] No find_psi_2nd_separatrix number of x-point < 1 --- src/physics/fluxsurfaces.jl | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/src/physics/fluxsurfaces.jl b/src/physics/fluxsurfaces.jl index 3293a66c..cd4f88ee 100644 --- a/src/physics/fluxsurfaces.jl +++ b/src/physics/fluxsurfaces.jl @@ -183,7 +183,7 @@ end find_psi_2nd_separatrix(eqt::IMAS.equilibrium__time_slice, PSI_interpolant::Interpolations.AbstractInterpolation) Returns psi of the second magentic separatrix. This relies on the fact that find_x_points! saves the x points in such a way -that the last one is the null with Z opposite to the first x point which is the closest in psi to the lcfs +that the last one is the null with Z opposite to the first x point which is the closest in psi to the lcfs. """ function find_psi_2nd_separatrix(eqt::IMAS.equilibrium__time_slice, PSI_interpolant::Interpolations.AbstractInterpolation) psi2nd = PSI_interpolant.(eqt.boundary.x_point[end].r, eqt.boundary.x_point[end].z) @@ -731,12 +731,14 @@ function flux_surfaces(eqt::equilibrium__time_slice{T}, B0::T, R0::T; upsample_f find_strike_points!(eqt) # secondary separatrix - psi2nd = find_psi_2nd_separatrix(eqt, PSI_interpolant) - tmp, _ = flux_surface(r, z, PSI, eqt.profiles_1d.psi, eqt.global_quantities.magnetic_axis.r, eqt.global_quantities.magnetic_axis.z, psi2nd, :encircling) - if !isempty(tmp) - (pr2nd, pz2nd) = tmp[1] - eqt.boundary_secondary_separatrix.outline.r = pr2nd - eqt.boundary_secondary_separatrix.outline.z = pz2nd + if length(eqt.boundary.x_point) > 1 + psi2nd = find_psi_2nd_separatrix(eqt, PSI_interpolant) + tmp, _ = flux_surface(r, z, PSI, eqt.profiles_1d.psi, eqt.global_quantities.magnetic_axis.r, eqt.global_quantities.magnetic_axis.z, psi2nd, :encircling) + if !isempty(tmp) + (pr2nd, pz2nd) = tmp[1] + eqt.boundary_secondary_separatrix.outline.r = pr2nd + eqt.boundary_secondary_separatrix.outline.z = pz2nd + end end return eqt