module diagnostics_field_m !! Module for writing out full field diagnostics for the BRAGINSKII model use netcdf use MPI use comm_handler_m, only : comm_handler_t use variable_m, only : variable_t use screen_io_m, only : get_stdout use precision_grillix_m, only : GP use perf_m, only : perf_start, perf_stop use error_handling_grillix_m, only: handle_error_netcdf, handle_error, error_info_t use status_codes_grillix_m, only : GRILLIX_ERR_OTHER use mesh_cart_m, only: mesh_cart_t use parallel_map_m, only : parallel_map_t use inplane_operators_m, only : inplane_operators_t use equilibrium_m, only: equilibrium_t use penalisation_m, only : penalisation_t use polars_m, only : polars_t use equilibrium_storage_m, only : equilibrium_storage_t use parallel_target_flux_m, only : parallel_target_flux_t use snapshot_m, only : snapshot_t, char_arr_t use diagnostics_packet_m, only : diagnostics_packet_t use neutrals_module_m, only: diff_co, diff_co_presmooth implicit none type, public :: diagnostics_field_t !! Datatype for full-field diagnostics contains procedure, public :: process_to_file end type contains subroutine process_to_file(self, diag_packet, comm_handler, equi, equi_on_cano, equi_on_stag, mesh_cano, & mesh_stag, map, penalisation_cano, penalisation_stag, polars_cano, polars_stag, & tau, ne, te, ti, pot, vort, upar, jpar, apar, apar_fluct, & neutrals_dens, neutrals_parmom, neutrals_pressure, & isnaps, idiag, start_new_file) !! Main diagnostics computation routine class(diagnostics_field_t), intent(inout) :: self !! Instance of class type(diagnostics_packet_t), intent(inout) :: diag_packet !! Diagnostic packet type(comm_handler_t), intent(in) :: comm_handler !! MPI communication handler class(equilibrium_t), intent(in) :: equi !! Magnetic equilibrium class(equilibrium_storage_t), intent(in) :: equi_on_cano !! Equilibrium storage on canonical plane enabling faster performance at certain locations class(equilibrium_storage_t), intent(in) :: equi_on_stag !! Equilibrium storage on staggered plane enabling faster performance at certain locations type(mesh_cart_t), intent(in) :: mesh_cano !! Mesh (canonical) within poloidal plane type(mesh_cart_t), intent(in) :: mesh_stag !! Mesh (staggered) within poloidal plane type(parallel_map_t), intent(in) :: map !! Parallel map type(penalisation_t), intent(in) :: penalisation_cano !! Penalisation (canonical) type(penalisation_t), intent(in) :: penalisation_stag !! Penalisation (staggered) type(polars_t), intent(in) :: polars_cano !! Polar grid and operators (canonical) type(polars_t), intent(in) :: polars_stag !! Polar grid and operators (staggered) real(GP), intent(in) :: tau !! Time type(variable_t), intent(in) :: ne !! Electron density type(variable_t), intent(in) :: te !! Electron temperature type(variable_t), intent(in) :: ti !! Ion temperature type(variable_t), intent(in) :: pot !! Electrostatic potential type(variable_t), intent(in) :: vort !! Generalised vorticity type(variable_t), intent(in) :: upar !! Parallel ion velocity type(variable_t), intent(in) :: jpar !! Parallel current type(variable_t), intent(in) :: apar !! Parallel electromagnetic potential type(variable_t), intent(in) :: apar_fluct !! Parallel electromagnetic potential fluctation type(variable_t), intent(in) :: neutrals_dens !! Neutrals density type(variable_t), intent(in) :: neutrals_parmom !! Neutrals parallel momentum type(variable_t), intent(in) :: neutrals_pressure !! Neutrals pressure integer, intent(in) :: isnaps !! Braginskii model snapshot number integer, intent(in) :: idiag !! Diagnostic snapshot number logical, intent(in) :: start_new_file !! Flag whether to write to new file integer :: ndiags type(char_arr_t), allocatable, dimension(:) :: diag_names type(variable_t), allocatable, dimension(:) :: diag_vars type(snapshot_t) :: diag_snapshot ndiags = 13 ! << IMPORTANT: MUST ADJUST WHEN ADDING DIAGNOSTICS !!! call perf_start('diagnostics_braginskii_field') allocate(diag_names(ndiags)) allocate(diag_vars(ndiags)) diag_names(1)%str = 'src_iz' call diag_vars(1)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(1)%vals = diag_packet%src_iz_ne diag_names(2)%str = 'src_rec' call diag_vars(2)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(2)%vals = diag_packet%src_rec_ne diag_names(3)%str = 'src_neut_ti' call diag_vars(3)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(3)%vals = diag_packet%src_neut_ti diag_names(4)%str = 'src_neut_pressure' call diag_vars(4)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(4)%vals = diag_packet%src_neut_pressure diag_names(5)%str = 'src_neut_rcy' call diag_vars(5)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(5)%vals = diag_packet%src_neut_rcy diag_names(6)%str = 'P_rad_rec' call diag_vars(6)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(6)%vals = diag_packet%p_rad_rec diag_names(7)%str = 'P_rad_imp' call diag_vars(7)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(7)%vals = diag_packet%p_rad_imp diag_names(8)%str = 'qpar_cond_e' call diag_vars(8)%init(comm_handler, mesh_cano%get_n_points(), staggered=.true.) diag_vars(8)%vals = diag_packet%qpar_cond_e diag_names(9)%str = 'qpar_cond_i' call diag_vars(9)%init(comm_handler, mesh_cano%get_n_points(), staggered=.true.) diag_vars(9)%vals = diag_packet%qpar_cond_i diag_names(10)%str = 'qpar_conv_e' call diag_vars(10)%init(comm_handler, mesh_cano%get_n_points(), staggered=.true.) diag_vars(10)%vals = diag_packet%qpar_conv_e diag_names(11)%str = 'qpar_conv_i' call diag_vars(11)%init(comm_handler, mesh_cano%get_n_points(), staggered=.true.) diag_vars(11)%vals = diag_packet%qpar_conv_i diag_names(12)%str = 'diff_co_presmooth' call diag_vars(12)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(12)%vals = diff_co_presmooth diag_names(13)%str = 'diff_co_postsmooth' call diag_vars(13)%init(comm_handler, mesh_cano%get_n_points(), staggered=.false.) diag_vars(13)%vals = diff_co%vals call diag_snapshot%init(comm_handler, mesh_cano, mesh_cano, 'diagnostics_field', & ndiags, diag_names) call diag_snapshot%write_snaps(comm_handler, isnaps, tau, diag_vars, & isubsnaps=idiag, allow_overwrite=.true.) deallocate(diag_names) deallocate(diag_vars) call perf_stop('diagnostics_braginskii_field') end subroutine end module