From fcf0ef42221b91153178e21e9243770a206e9dc3 Mon Sep 17 00:00:00 2001 From: Andreas Zach Date: Wed, 25 Mar 2026 00:51:18 +0100 Subject: [PATCH] Introduce a module for writeData.f90 --- QL-Balance/src/base/get_dql.f90 | 5 ++- QL-Balance/src/base/ramp_coil.f90 | 3 +- QL-Balance/src/base/time_evolution.f90 | 9 +++-- QL-Balance/src/base/writeData.f90 | 39 ++++++++++++------- .../src/stellarator/time_evol_stell.f90 | 3 +- 5 files changed, 38 insertions(+), 21 deletions(-) diff --git a/QL-Balance/src/base/get_dql.f90 b/QL-Balance/src/base/get_dql.f90 index a40c8e7d..2cc03985 100644 --- a/QL-Balance/src/base/get_dql.f90 +++ b/QL-Balance/src/base/get_dql.f90 @@ -26,6 +26,7 @@ subroutine get_dql use QLBalance_kinds, only: dp use PolyLagrangeInterpolation use logger_m, only: log_debug + use writeData_m, only: write_fields_currs_transp_coefs_to_h5, write_D_one_over_nu_to_h5 implicit none @@ -487,8 +488,8 @@ subroutine get_dql if (modulo(time_ind, save_prof_time_step) .eq. 0) then if (suppression_mode .eqv. .false.) then - call write_fields_currs_transp_coefs_to_h5 - call write_D_one_over_nu_to_h5 + call write_fields_currs_transp_coefs_to_h5(time_ind) + call write_D_one_over_nu_to_h5(time_ind) end if end if diff --git a/QL-Balance/src/base/ramp_coil.f90 b/QL-Balance/src/base/ramp_coil.f90 index 060e837d..7f0d28e3 100644 --- a/QL-Balance/src/base/ramp_coil.f90 +++ b/QL-Balance/src/base/ramp_coil.f90 @@ -390,6 +390,7 @@ subroutine check_linear_discr_pen_ratio br_stopping, discr_reached, br_abs, br_predicted, write_br_dqle22_time_data, write_kin_prof_data_to_disk use control_mod, only: suppression_mode use h5mod, only: write_reason_for_stop_to_h5 + use writeData_m, only: write_fields_currs_transp_coefs_to_h5 implicit none @@ -407,7 +408,7 @@ subroutine check_linear_discr_pen_ratio write(*,*) 'discrepancy to linearly predicted value of Br_abs_res > delta' if (modulo(time_ind, save_prof_time_step) .ne. 0) then if (suppression_mode .eqv. .false.) then - CALL write_fields_currs_transp_coefs_to_h5 + call write_fields_currs_transp_coefs_to_h5(time_ind) end if end if if (suppression_mode .eqv. .false.) then diff --git a/QL-Balance/src/base/time_evolution.f90 b/QL-Balance/src/base/time_evolution.f90 index 92265788..7adfcf0c 100644 --- a/QL-Balance/src/base/time_evolution.f90 +++ b/QL-Balance/src/base/time_evolution.f90 @@ -175,11 +175,11 @@ subroutine doStep(this) use baseparam_mod, only: factolmax, factolred use plasma_parameters, only: params, params_beg, params_begbeg, limit_temps_from_below use logger_m, only: log_debug - use recstep_mod, only: timstep_arr - use recstep_mod, only: tol + use recstep_mod, only: timstep_arr, tol use restart_mod, only: redostep use transp_coeffs_mod, only: rescale_transp_coeffs_by_ant_fac, & compute_antenna_factor_from_Ipar + use writeData_m, only: writefort9999 implicit none @@ -200,7 +200,7 @@ subroutine doStep(this) call message_Br_Dqle_values if (data_verbosity >= 2) then - call writefort9999 + call writefort9999(dqle11_prev, dqli11_prev) end if if (.true.) then @@ -446,6 +446,7 @@ subroutine write_br_dqle22_time_data subroutine check_linear_discr_pen_ratio + use writeData_m, only: write_fields_currs_transp_coefs_to_h5 implicit none @@ -463,7 +464,7 @@ subroutine check_linear_discr_pen_ratio write(*,*) 'discrepancy to linearly predicted value of Br_abs_res > delta' if (modulo(time_ind, save_prof_time_step) .ne. 0) then if (suppression_mode .eqv. .false.) then - CALL write_fields_currs_transp_coefs_to_h5 + call write_fields_currs_transp_coefs_to_h5(time_ind) end if end if if (suppression_mode .eqv. .false.) then diff --git a/QL-Balance/src/base/writeData.f90 b/QL-Balance/src/base/writeData.f90 index c9729ec3..5ae3d2f0 100644 --- a/QL-Balance/src/base/writeData.f90 +++ b/QL-Balance/src/base/writeData.f90 @@ -1,16 +1,29 @@ -subroutine write_fields_currs_transp_coefs_to_h5 +module writeData_m + use QLBalance_kinds, only: dp + + implicit none + private + + public :: write_fields_currs_transp_coefs_to_h5 + public :: writefort9999 + public :: writefort9999_stellarator + public :: write_D_one_over_nu_to_h5 + +contains + +subroutine write_fields_currs_transp_coefs_to_h5(time_ind) use grid_mod, only: npoib, dqle11, dqle12, dqle22, dqli11, dqli12, dqli22, & T_EM_phi_e, T_EM_phi_i - use baseparam_mod, only: e_charge, p_mass, c, e_mass, ev + use baseparam_mod, only: c use control_mod, only: ihdf5IO, data_verbosity, misalign_diffusion use logger_m, only: log_debug use wave_code_data use QLbalance_diag, only: iunit_diag - use time_evolution, only: time_ind use h5mod implicit none + integer, intent(in) :: time_ind integer :: ipoi character(len=1024) :: tempch @@ -145,21 +158,20 @@ subroutine write_dql_Br_Jp_profiles_to_hdf5(tempch) end subroutine -subroutine writefort9999 +subroutine writefort9999(dqle11_prev, dqli11_prev) use grid_mod, only: dqle11, dqli11, rb, rc, npoib use QLbalance_diag, only: timscal_dql, timscal_dqli, ind_dqle, ind_dqli - use time_evolution, only: dqle11_prev, dqli11_prev, determine_Dql_diagnostic use h5mod use logger_m, only: log_debug implicit none + real(dp), dimension(:), intent(in) :: dqle11_prev + real(dp), dimension(:), intent(in) :: dqli11_prev integer :: ipoi character(256) :: buf - call determine_Dql_diagnostic - write(buf, '(A,ES12.4,A,ES12.4)') 'timscal_dqle = ', & sngl(timscal_dql), ' timscal_dqli = ', sngl(timscal_dqli) call log_debug(trim(buf)) @@ -197,21 +209,20 @@ subroutine writefort9999 end subroutine -subroutine writefort9999_stellarator +subroutine writefort9999_stellarator(dqle11_prev, dqli11_prev) use grid_mod, only: dqle11, dqli11, rb, rc, npoib use QLbalance_diag, only: timscal_dql, timscal_dqli, ind_dqle, ind_dqli - use time_evolution_stellarator, only: dqle11_prev, dqli11_prev, determine_Dql_diagnostic use h5mod use logger_m, only: log_debug implicit none + real(dp), dimension(:), intent(in) :: dqle11_prev + real(dp), dimension(:), intent(in) :: dqli11_prev integer :: ipoi character(256) :: buf - call determine_Dql_diagnostic - write(buf, '(A,ES12.4,A,ES12.4)') 'timscal_dqle = ', & sngl(timscal_dql), ' timscal_dqli = ', sngl(timscal_dqli) call log_debug(trim(buf)) @@ -249,15 +260,15 @@ subroutine writefort9999_stellarator end subroutine -subroutine write_D_one_over_nu_to_h5 +subroutine write_D_one_over_nu_to_h5(time_ind) use grid_mod, only: Donue11, Donue12, Donue21, Donue22, & Donui11, Donui12, Donui21, Donui22 use h5mod - use time_evolution, only: time_ind implicit none + integer, intent(in) :: time_ind character(256) :: tempch tempch = "/"//trim(h5_mode_groupname)//"/LinearProfiles" write (tempch, '(A,"/",I0,"/")') trim(tempch), time_ind @@ -284,3 +295,5 @@ subroutine write_D_one_over_nu_to_h5 CALL h5_close(h5_id) CALL h5_deinit() end subroutine + +end module writeData_m diff --git a/QL-Balance/src/stellarator/time_evol_stell.f90 b/QL-Balance/src/stellarator/time_evol_stell.f90 index c233e653..3d34ea93 100644 --- a/QL-Balance/src/stellarator/time_evol_stell.f90 +++ b/QL-Balance/src/stellarator/time_evol_stell.f90 @@ -121,6 +121,7 @@ subroutine runTimeEvolution(this) use transp_coeffs_mod, only: rescale_transp_coeffs_by_ant_fac use recstep_mod, only: timstep_arr use grid_mod, only: Ipar + use writeData_m, only: writefort9999_stellarator implicit none @@ -148,7 +149,7 @@ subroutine runTimeEvolution(this) call message_Br_Dqle_values if (data_verbosity >= 2) then - call writefort9999_stellarator + call writefort9999_stellarator(dqle11_prev, dqli11_prev) end if if (.true.) then