diff --git a/src/gray_project.f90 b/src/gray_project.f90 index f05ce1b..74e51c0 100644 --- a/src/gray_project.f90 +++ b/src/gray_project.f90 @@ -210,6 +210,7 @@ contains ! local variables real(wp_) :: P_tot, I_tot real(wp_) :: dVdrho_t, dAdrho_t + integer :: i_max ! ⟨ρ_t⟩, Δρ_t using dP/dV⋅dV as a PDF ! Note: the factor √8 is to match the full-width at max/e of a Gaussian @@ -218,7 +219,7 @@ contains rho_avg_P = sum(self%rho_t * dPdV*self%dV)/P_tot drho_avg_P = sqrt(8 * sum((self%rho_t - rho_avg_P)**2 * dPdV*self%dV)/P_tot) else - rho_avg_P = 0 + rho_avg_P = 1 drho_avg_P = 0 end if @@ -228,7 +229,7 @@ contains rho_avg_J = sum(self%rho_t * abs(J_phi)*self%dA)/I_tot drho_avg_J = sqrt(8* sum((self%rho_t - rho_avg_J)**2 * abs(J_phi)*self%dA)/I_tot) else - rho_avg_J = 0 + rho_avg_J = 1 drho_avg_J = 0 end if @@ -238,10 +239,12 @@ contains dPdV_peak = P_tot * 2/(sqrt(pi) * drho_avg_P * dVdrho_t) ! Compute numerically argmax, max and full-width at max/e of dP/dV - call profwidth(self%rho_t, dPdV, rho_max_P, dPdV_max, drho_max_P) + call profwidth(self%rho_t, dPdV, i_max, drho_max_P) + rho_max_P = self%rho_t(i_max) + dPdV_max = dPdV(i_max) else dPdV_peak = 0 - rho_max_P = 0 + rho_max_P = 1 dPdV_max = 0 drho_max_P = 0 end if @@ -250,46 +253,55 @@ contains ! Compute max of J_φ as if it were a Gaussian call equil%flux_average(equil%tor2pol(rho_avg_J), dAdrho_t=dAdrho_t, & ratio_astra_tor=ratio_Ja_max, ratio_jintrac_tor=ratio_Jb_max) - J_phi_peak = I_tot * 2/(sqrt(pi)*drho_avg_J*dAdrho_t) + J_phi_peak = sum(J_phi*self%dA) * 2/(sqrt(pi)*drho_avg_J*dAdrho_t) ! Compute numerically argmax, max and full-width at max/e of J_φ - call profwidth(self%rho_t, abs(J_phi), rho_max_J, J_phi_max, drho_max_J) + call profwidth(self%rho_t, abs(J_phi), i_max, drho_max_J) + rho_max_J = self%rho_t(i_max) + J_phi_max = J_phi(i_max) else J_phi_peak = 0 - rho_max_J = 0 + rho_max_J = 1 J_phi_max = 0 drho_max_J = 0 ratio_Ja_max = 0 ratio_Jb_max = 0 end if + + ! Signal that J_φ has alternating sign + if ((I_tot - abs(sum(J_phi*self%dA)))/ I_tot > 0.1) then + rho_avg_J = -rho_avg_J + drho_avg_J = -drho_avg_J + rho_max_J = -rho_max_J + drho_max_J = -drho_max_J + end if end subroutine statistics - subroutine profwidth(x, y, x_max, y_max, width) - ! Computes the argmax, max and full-width at max/e of the curve y(x) + subroutine profwidth(x, y, i_max, width) + ! Computes the argmax and full-width at max/e of the curve y(x) use const_and_precisions, only : emn1 use utils, only : locate_unord, linear_interp ! subroutine arguments real(wp_), intent(in) :: x(:), y(:) - real(wp_), intent(out) :: x_max, y_max, width + integer, intent(out) :: i_max + real(wp_), intent(out) :: width ! local variables - integer :: i, i_max, left, right + integer :: i, left, right integer :: nsols, sols(size(y)) real(wp_) :: x_left, x_right ! Find maximum i_max = maxloc(y, dim=1) - x_max = x(i_max) - y_max = y(i_max) ! Find width width = 0 ! Find solutions to y(x) = y_max⋅e⁻¹ - call locate_unord(y, y_max*emn1, sols, size(y), nsols) + call locate_unord(y, y(i_max)*emn1, sols, size(y), nsols) if (nsols < 2) return ! Find the two closest solutions to the maximum @@ -301,8 +313,8 @@ contains end do ! Interpolate linearly x(y) for better accuracy - x_left = linear_interp(y(left), x(left), y(left+1), x(left+1), y_max*emn1) - x_right = linear_interp(y(right), x(right), y(right+1), x(right+1), y_max*emn1) + x_left = linear_interp(y(left), x(left), y(left+1), x(left+1), y(i_max)*emn1) + x_right = linear_interp(y(right), x(right), y(right+1), x(right+1), y(i_max)*emn1) width = x_right - x_left end subroutine profwidth diff --git a/tests/03-TCV/outputs/summary.7.txt b/tests/03-TCV/outputs/summary.7.txt index 7f71639..f31435e 100644 --- a/tests/03-TCV/outputs/summary.7.txt +++ b/tests/03-TCV/outputs/summary.7.txt @@ -20,4 +20,4 @@ # COD iwarm ilarm imx ieccd: 2 5 -20 3 # COD ipec nrho istpr istpl: 1 1001 5 5 # I_cd P_abs J_φ_peak dPdV_peak ρ_max_J ρ_avg_J ρ_max_P ρ_avg_P Δρ_avg_J Δρ_avg_P ratio_Ja_max ratio_Jb_max s_max ψ χ index_rt J_φ_max dPdV_max Δρ_J Δρ_P P0 cpl1 cpl2 - 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.20000000E+002 -8.57133940E+001 -1.19957807E+000 2 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000 0.00000000E+000 0.00000000E+000 + 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.20000000E+002 -8.57133940E+001 -1.19957807E+000 2 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000 0.00000000E+000 0.00000000E+000 diff --git a/tests/04-JT60SA/outputs/summary.7.txt b/tests/04-JT60SA/outputs/summary.7.txt index 31afed8..8b41ed6 100644 --- a/tests/04-JT60SA/outputs/summary.7.txt +++ b/tests/04-JT60SA/outputs/summary.7.txt @@ -20,6 +20,6 @@ # COD iwarm ilarm imx ieccd: 2 4 -20 3 # COD ipec nrho istpr istpl: 1 2001 5 5 # I_cd P_abs J_φ_peak dPdV_peak ρ_max_J ρ_avg_J ρ_max_P ρ_avg_P Δρ_avg_J Δρ_avg_P ratio_Ja_max ratio_Jb_max s_max ψ χ index_rt J_φ_max dPdV_max Δρ_J Δρ_P P0 cpl1 cpl2 - 1.51692135E-009 4.54436052E-003 7.77985972E-007 1.71519960E-002 6.02647911E-001 6.01576172E-001 6.01680996E-001 6.01722246E-001 2.88901706E-003 1.76361881E-003 1.10064529E+000 1.03922461E+000 3.29500000E+002 2.14149518E+001 2.42883589E+001 1 1.06222518E-006 1.69653156E-002 9.94551841E-004 1.79307710E-003 1.00000000E+000 9.71180646E-001 2.88193543E-002 --4.77656743E-010 1.30591241E-003 3.37694971E-007 5.51787071E-003 8.95089111E-001 8.95558447E-001 8.96049366E-001 8.95746996E-001 1.53036611E-003 1.31588770E-003 1.17339615E+000 1.02259980E+000 7.25800000E+002 7.54874568E+000 -6.28239682E+000 3 2.60491790E-007 4.48770460E-003 2.11107756E-003 1.80616140E-003 9.66767251E-001 0.00000000E+000 0.00000000E+000 - 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 6.19300000E+002 -8.24512544E+001 6.28239682E+000 4 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 2.86883888E-002 0.00000000E+000 0.00000000E+000 + 1.51692135E-009 4.54436052E-003 7.77985972E-007 1.71519960E-002 -6.02647911E-001 -6.01576172E-001 6.01680996E-001 6.01722246E-001 -2.88901706E-003 1.76361881E-003 1.10064529E+000 1.03922461E+000 3.29500000E+002 2.14149518E+001 2.42883589E+001 1 1.06222518E-006 1.69653156E-002 -9.94551841E-004 1.79307710E-003 1.00000000E+000 9.71180646E-001 2.88193543E-002 +-4.77656743E-010 1.30591241E-003 3.37694971E-007 5.51787071E-003 -8.95089111E-001 -8.95558447E-001 8.96049366E-001 8.95746996E-001 -1.53036611E-003 1.31588770E-003 1.17339615E+000 1.02259980E+000 7.25800000E+002 7.54874568E+000 -6.28239682E+000 3 2.60491790E-007 4.48770460E-003 -2.11107756E-003 1.80616140E-003 9.66767251E-001 0.00000000E+000 0.00000000E+000 + 0.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 1.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 6.19300000E+002 -8.24512544E+001 6.28239682E+000 4 0.00000000E+000 0.00000000E+000 0.00000000E+000 0.00000000E+000 2.86883888E-002 0.00000000E+000 0.00000000E+000 diff --git a/tests/05-JT60SA-startup/outputs/summary.7.txt b/tests/05-JT60SA-startup/outputs/summary.7.txt index 4af8be0..3c2679c 100644 --- a/tests/05-JT60SA-startup/outputs/summary.7.txt +++ b/tests/05-JT60SA-startup/outputs/summary.7.txt @@ -20,10 +20,10 @@ # COD iwarm ilarm imx ieccd: 2 4 -20 3 # COD ipec nrho istpr istpl: 1 2001 5 5 # I_cd P_abs J_φ_peak dPdV_peak ρ_max_J ρ_avg_J ρ_max_P ρ_avg_P Δρ_avg_J Δρ_avg_P ratio_Ja_max ratio_Jb_max s_max ψ χ index_rt J_φ_max dPdV_max Δρ_J Δρ_P P0 cpl1 cpl2 - 1.09082396E-006 1.18757272E-001 -4.29115678E-004 -1.44252143E-001 7.67010753E-001 7.62216152E-001 7.62978352E-001 7.62550997E-001 1.39366449E-002 7.41198052E-003 1.02093556E+000 1.01041356E+000 4.33200000E+002 -8.31222486E+001 9.43560237E+000 2 1.67706595E-003 1.94713423E-001 2.35345649E-003 7.25633640E-003 1.00000000E+000 9.77799467E-001 2.22005331E-002 - 9.94865477E-008 7.31906936E-004 -3.50540708E-006 -9.36832811E-004 7.75802182E-001 7.78029924E-001 7.77841628E-001 7.78840037E-001 1.19589063E-002 6.94332053E-003 1.02092468E+000 1.01040817E+000 7.24600000E+002 6.74902073E+000 -5.32965597E+000 5 9.57658383E-006 1.27596579E-003 4.07144595E-003 6.92180492E-003 8.61678670E-001 2.66131402E-003 9.97338686E-001 - 4.49555314E-008 9.48537732E-004 -4.07568278E-006 -1.36699445E-003 8.11726185E-001 8.07325441E-001 8.08218554E-001 8.07968922E-001 1.14139213E-002 6.03393306E-003 1.02090168E+000 1.01039679E+000 9.26400000E+002 -8.32509793E+001 5.32965597E+000 6 1.22712358E-005 1.90616572E-003 2.30445993E-003 5.96877041E-003 1.95640583E-002 4.34298009E-002 9.56570199E-001 --7.83315679E-010 4.40706477E-006 -2.07836538E-008 -4.45175802E-006 6.87427661E-001 6.90422730E-001 6.89916553E-001 6.91288993E-001 1.50339943E-002 9.50373557E-003 1.02096106E+000 1.01042618E+000 9.50100000E+002 7.31932505E+000 7.96902175E+000 11 5.29735208E-008 5.76098737E-006 5.08608352E-003 9.41593883E-003 2.29124969E-003 0.00000000E+000 0.00000000E+000 - 4.51312962E-006 1.99291903E-001 -6.59730755E-004 -2.11481112E-001 6.98678817E-001 6.93617676E-001 6.94287713E-001 6.93536568E-001 1.63845675E-002 9.02658927E-003 1.02096109E+000 1.01042619E+000 9.49800000E+002 -8.26806750E+001 -7.96902175E+000 12 3.34567968E-003 2.68462628E-001 2.17779776E-003 8.76562927E-003 8.58655513E-001 0.00000000E+000 0.00000000E+000 + 1.09082396E-006 1.18757272E-001 -4.29115678E-004 -1.44252143E-001 -7.67010753E-001 -7.62216152E-001 7.62978352E-001 7.62550997E-001 -1.39366449E-002 7.41198052E-003 1.02093556E+000 1.01041356E+000 4.33200000E+002 -8.31222486E+001 9.43560237E+000 2 1.67706595E-003 1.94713423E-001 -2.35345649E-003 7.25633640E-003 1.00000000E+000 9.77799467E-001 2.22005331E-002 + 9.94865477E-008 7.31906936E-004 -3.50540708E-006 -9.36832811E-004 -7.75802182E-001 -7.78029924E-001 7.77841628E-001 7.78840037E-001 -1.19589063E-002 6.94332053E-003 1.02092468E+000 1.01040817E+000 7.24600000E+002 6.74902073E+000 -5.32965597E+000 5 9.57658383E-006 1.27596579E-003 -4.07144595E-003 6.92180492E-003 8.61678670E-001 2.66131402E-003 9.97338686E-001 + 4.49555314E-008 9.48537732E-004 -4.07568278E-006 -1.36699445E-003 -8.11726185E-001 -8.07325441E-001 8.08218554E-001 8.07968922E-001 -1.14139213E-002 6.03393306E-003 1.02090168E+000 1.01039679E+000 9.26400000E+002 -8.32509793E+001 5.32965597E+000 6 1.22712358E-005 1.90616572E-003 -2.30445993E-003 5.96877041E-003 1.95640583E-002 4.34298009E-002 9.56570199E-001 +-7.83315679E-010 4.40706477E-006 -2.07836538E-008 -4.45175802E-006 -6.87427661E-001 -6.90422730E-001 6.89916553E-001 6.91288993E-001 -1.50339943E-002 9.50373557E-003 1.02096106E+000 1.01042618E+000 9.50100000E+002 7.31932505E+000 7.96902175E+000 11 5.29735208E-008 5.76098737E-006 -5.08608352E-003 9.41593883E-003 2.29124969E-003 0.00000000E+000 0.00000000E+000 + 4.51312962E-006 1.99291903E-001 -6.59730755E-004 -2.11481112E-001 -6.98678817E-001 -6.93617676E-001 6.94287713E-001 6.93536568E-001 -1.63845675E-002 9.02658927E-003 1.02096109E+000 1.01042619E+000 9.49800000E+002 -8.26806750E+001 -7.96902175E+000 12 3.34567968E-003 2.68462628E-001 -2.17779776E-003 8.76562927E-003 8.58655513E-001 0.00000000E+000 0.00000000E+000 2.92468740E-010 3.00247904E-006 -3.87384950E-008 -3.94324559E-004 6.11080841E-001 6.11093245E-001 6.11080841E-001 6.11081132E-001 1.44628454E-003 8.00518460E-005 1.02089926E+000 1.01039560E+000 1.18250000E+003 6.51056937E+000 -6.98374197E+000 13 9.60260701E-008 5.52072301E-005 7.31377568E-004 7.23821974E-004 8.08468354E-004 0.00000000E+000 0.00000000E+000 1.41767167E-007 4.17030110E-003 -1.55290680E-005 -6.14842229E-001 6.18552626E-001 6.18598876E-001 6.18552626E-001 6.18553161E-001 1.73225661E-003 7.06369029E-005 1.02091159E+000 1.01040170E+000 1.17790000E+003 -8.34894306E+001 6.98374197E+000 14 4.33974961E-005 7.59608065E-002 7.63590695E-004 7.29963129E-004 1.78070522E-002 0.00000000E+000 0.00000000E+000 diff --git a/tests/06-ITER-startup/outputs/mixed/summary.7.txt b/tests/06-ITER-startup/outputs/mixed/summary.7.txt index 8364b20..af231e6 100644 --- a/tests/06-ITER-startup/outputs/mixed/summary.7.txt +++ b/tests/06-ITER-startup/outputs/mixed/summary.7.txt @@ -20,5 +20,5 @@ # COD iwarm ilarm imx ieccd: 2 4 -20 3 # COD ipec nrho istpr istpl: 1 501 5 5 # I_cd P_abs J_φ_peak dPdV_peak ρ_max_J ρ_avg_J ρ_max_P ρ_avg_P Δρ_avg_J Δρ_avg_P ratio_Ja_max ratio_Jb_max s_max ψ χ index_rt J_φ_max dPdV_max Δρ_J Δρ_P P0 cpl1 cpl2 --6.84761868E-007 6.70958966E-003 8.41854855E-006 5.31882140E-004 3.31457142E-001 3.48277613E-001 3.50892072E-001 3.50807962E-001 6.39977138E-002 4.07899012E-002 1.02361209E+000 1.00001579E+000 8.00000000E+002 3.36364721E+001 -3.41970443E+001 1 1.07442271E-005 5.26793259E-004 3.16299561E-002 4.14324774E-002 5.71190361E-001 0.00000000E+000 0.00000000E+000 - 1.67484070E-005 1.64840181E-001 1.99720749E-004 1.28981575E-002 3.70133716E-001 3.56067546E-001 3.53038507E-001 3.53670278E-001 6.44956402E-002 4.10084319E-002 1.02377382E+000 1.00001668E+000 8.00000000E+002 -5.63635279E+001 3.41970443E+001 2 2.44482189E-004 1.27979295E-002 2.96804467E-002 4.15437790E-002 4.28809639E-001 0.00000000E+000 0.00000000E+000 +-6.84761868E-007 6.70958966E-003 8.41854855E-006 5.31882140E-004 -3.31457142E-001 -3.48277613E-001 3.50892072E-001 3.50807962E-001 -6.39977138E-002 4.07899012E-002 1.02361209E+000 1.00001579E+000 8.00000000E+002 3.36364721E+001 -3.41970443E+001 1 1.07442271E-005 5.26793259E-004 -3.16299561E-002 4.14324774E-002 5.71190361E-001 0.00000000E+000 0.00000000E+000 + 1.67484070E-005 1.64840181E-001 1.99720749E-004 1.28981575E-002 -3.70133716E-001 -3.56067546E-001 3.53038507E-001 3.53670278E-001 -6.44956402E-002 4.10084319E-002 1.02377382E+000 1.00001668E+000 8.00000000E+002 -5.63635279E+001 3.41970443E+001 2 2.44482189E-004 1.27979295E-002 -2.96804467E-002 4.15437790E-002 4.28809639E-001 0.00000000E+000 0.00000000E+000 diff --git a/tests/06-ITER-startup/outputs/x-mode/summary.7.txt b/tests/06-ITER-startup/outputs/x-mode/summary.7.txt index 8cbfc87..093cad3 100644 --- a/tests/06-ITER-startup/outputs/x-mode/summary.7.txt +++ b/tests/06-ITER-startup/outputs/x-mode/summary.7.txt @@ -20,4 +20,4 @@ # COD iwarm ilarm imx ieccd: 2 4 -20 3 # COD ipec nrho istpr istpl: 1 501 5 5 # I_cd P_abs J_φ_peak dPdV_peak ρ_max_J ρ_avg_J ρ_max_P ρ_avg_P Δρ_avg_J Δρ_avg_P ratio_Ja_max ratio_Jb_max s_max ψ χ index_rt J_φ_max dPdV_max Δρ_J Δρ_P P0 cpl1 cpl2 - 3.90425628E-005 3.84383056E-001 4.66399221E-004 3.00793382E-002 3.70133716E-001 3.55972396E-001 3.53038507E-001 3.53578490E-001 6.45085061E-002 4.10148056E-002 1.02377231E+000 1.00001666E+000 8.00000000E+002 -5.63635279E+001 3.41970443E+001 2 5.70581310E-004 2.98529841E-002 2.96939387E-002 4.15444773E-002 1.00000000E+000 0.00000000E+000 0.00000000E+000 +-3.90425628E-005 3.84383056E-001 4.66399221E-004 3.00793382E-002 -3.70133716E-001 -3.55972396E-001 3.53038507E-001 3.53578490E-001 -6.45085061E-002 4.10148056E-002 1.02377231E+000 1.00001666E+000 8.00000000E+002 -5.63635279E+001 3.41970443E+001 2 5.70581310E-004 2.98529841E-002 -2.96939387E-002 4.15444773E-002 1.00000000E+000 0.00000000E+000 0.00000000E+000