| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | #ifdef USE_FORTPLOT | ||
| 2 | module diag_meiss | ||
| 3 | !> Diagnostic routines for field_can_meiss.f90 | ||
| 4 | !> Provides visualization and analysis tools for canonical coordinate transformations | ||
| 5 | |||
| 6 | use, intrinsic :: iso_fortran_env, only: dp => real64 | ||
| 7 | use field_can_meiss, only: rh_can, n_r, xmin, xmax | ||
| 8 | use fortplot, only: figure, plot, savefig, xlabel, ylabel, title | ||
| 9 | |||
| 10 | implicit none | ||
| 11 | private | ||
| 12 | |||
| 13 | public :: plot_rh_can_vs_rc | ||
| 14 | |||
| 15 | contains | ||
| 16 | |||
| 17 | ✗ | subroutine plot_rh_can_vs_rc(i_th_in, i_phi_in, filename) | |
| 18 | !> Plot rh_can over r_c for fixed i_th and i_phi indices | ||
| 19 | !> | ||
| 20 | !> @param i_th_in Theta grid index (default: 1) | ||
| 21 | !> @param i_phi_in Phi grid index (default: 1) | ||
| 22 | !> @param filename Output filename (default: "diag_meiss.pdf") | ||
| 23 | |||
| 24 | integer, intent(in), optional :: i_th_in, i_phi_in | ||
| 25 | character(len=*), intent(in), optional :: filename | ||
| 26 | |||
| 27 | ! Local variables | ||
| 28 | integer :: i_th, i_phi | ||
| 29 | integer :: i, n_points | ||
| 30 | ✗ | real(dp), dimension(:), allocatable :: r_c_array, dz1_vals, dz2_vals | |
| 31 | real(dp), dimension(2) :: z, dz | ||
| 32 | character(len=100) :: output_file, output_file1, output_file2 | ||
| 33 | character(len=50) :: plot_title | ||
| 34 | |||
| 35 | ! Set defaults | ||
| 36 | ✗ | i_th = 1 | |
| 37 | ✗ | i_phi = 1 | |
| 38 | ✗ | if (present(i_th_in)) i_th = i_th_in | |
| 39 | ✗ | if (present(i_phi_in)) i_phi = i_phi_in | |
| 40 | |||
| 41 | ✗ | output_file = "diag_meiss.pdf" | |
| 42 | ✗ | if (present(filename)) output_file = trim(filename) | |
| 43 | |||
| 44 | ! Create separate filenames for the two plots | ||
| 45 | ✗ | output_file1 = trim(output_file(1:len_trim(output_file) - 4))//"_dz1.pdf" | |
| 46 | ✗ | output_file2 = trim(output_file(1:len_trim(output_file) - 4))//"_dz2.pdf" | |
| 47 | |||
| 48 | ! Create r_c array - use same grid as field_can_meiss | ||
| 49 | ✗ | n_points = n_r | |
| 50 | ✗ | allocate (r_c_array(n_points)) | |
| 51 | ✗ | allocate (dz1_vals(n_points)) | |
| 52 | ✗ | allocate (dz2_vals(n_points)) | |
| 53 | |||
| 54 | ! Fill r_c array | ||
| 55 | ✗ | do i = 1, n_points | |
| 56 | ✗ | r_c_array(i) = xmin(1) + (xmax(1) - xmin(1))*real(i - 1, dp)/real(n_points - 1, dp) | |
| 57 | end do | ||
| 58 | |||
| 59 | ! Initialize z (lam_phi=0, chi_gauge=0 for diagnostic purposes) | ||
| 60 | ✗ | z = [0.0_dp, 0.0_dp] | |
| 61 | |||
| 62 | ! Compute both components of rh_can | ||
| 63 | ✗ | do i = 1, n_points | |
| 64 | ✗ | call rh_can(r_c_array(i), z, dz, i_th, i_phi) | |
| 65 | ✗ | dz1_vals(i) = dz(1) ! dz(1) = -hr/hp | |
| 66 | ✗ | dz2_vals(i) = dz(2) ! dz(2) = Ar + Ap*dz(1) | |
| 67 | end do | ||
| 68 | |||
| 69 | ! Create two separate plots for better visualization | ||
| 70 | |||
| 71 | ! First plot: dz(1) = -hr/hp | ||
| 72 | ✗ | call figure() | |
| 73 | ✗ | call plot(r_c_array, dz1_vals, label="dz(1) = -hr/hp", linestyle="b-") | |
| 74 | ✗ | call xlabel("r_c") | |
| 75 | ✗ | call ylabel("dz(1) = -hr/hp") | |
| 76 | ✗ | write (plot_title, '(A,I0,A,I0,A)') "dz(1) vs r_c (i_th=", i_th, ", i_phi=", i_phi, ")" | |
| 77 | ✗ | call title(trim(plot_title)) | |
| 78 | |||
| 79 | ! Save first plot | ||
| 80 | ✗ | call savefig(trim(output_file1)) | |
| 81 | |||
| 82 | ! Second plot: dz(2) = Ar + Ap*dz(1) | ||
| 83 | ✗ | call figure() | |
| 84 | ✗ | call plot(r_c_array, dz2_vals, label="dz(2) = Ar + Ap*dz(1)", linestyle="r-") | |
| 85 | ✗ | call xlabel("r_c") | |
| 86 | ✗ | call ylabel("dz(2) = Ar + Ap*dz(1)") | |
| 87 | ✗ | write (plot_title, '(A,I0,A,I0,A)') "dz(2) vs r_c (i_th=", i_th, ", i_phi=", i_phi, ")" | |
| 88 | ✗ | call title(trim(plot_title)) | |
| 89 | |||
| 90 | ! Save second plot | ||
| 91 | ✗ | call savefig(trim(output_file2)) | |
| 92 | |||
| 93 | ! Cleanup | ||
| 94 | ✗ | deallocate (r_c_array, dz1_vals, dz2_vals) | |
| 95 | |||
| 96 | ✗ | print *, "Diagnostic plots saved to: ", trim(output_file1), " and ", trim(output_file2) | |
| 97 | |||
| 98 | ✗ | end subroutine plot_rh_can_vs_rc | |
| 99 | |||
| 100 | end module diag_meiss | ||
| 101 | #else | ||
| 102 | module diag_meiss | ||
| 103 | !> Stub module when fortplot is not available (e.g., nvfortran builds) | ||
| 104 | implicit none | ||
| 105 | private | ||
| 106 | public :: plot_rh_can_vs_rc | ||
| 107 | contains | ||
| 108 | subroutine plot_rh_can_vs_rc(i_th_in, i_phi_in, filename) | ||
| 109 | integer, intent(in), optional :: i_th_in, i_phi_in | ||
| 110 | character(len=*), intent(in), optional :: filename | ||
| 111 | print *, "Warning: plot_rh_can_vs_rc requires fortplot (disabled for this build)" | ||
| 112 | end subroutine | ||
| 113 | end module diag_meiss | ||
| 114 | #endif | ||
| 115 |