module diagnostics_lineout_m !! Module for computing lineout diagnostics for the BRAGINSKII model use MPI use netcdf use comm_handler_m, only : comm_handler_t use screen_io_m, only : get_stdout, get_stderr use parallelisation_setup_m, only : is_rank_info_writer use error_handling_grillix_m, only: handle_error_netcdf, handle_error, error_info_t use status_codes_grillix_m, only : GRILLIX_ERR_OTHER use diagnostic_variable_m, only : diagnostic_variable_t use variable_m, only : variable_t use precision_grillix_m, only : GP, MPI_GP use perf_m, only : perf_start, perf_stop use sources_external_m, only : sources_external_t 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 diagnostics_helpers_m, only : interpolate_lineout use diagnostics_group_m, only : diagnostics_group_t use diagnostics_packet_m, only : diagnostics_packet_t implicit none type, private :: lineout_coords_t !! Helper type to store variable length arrays of coordinates and names real(GP), dimension(:), allocatable, private :: x !! x-coordinates of sample points real(GP), dimension(:), allocatable, private :: y !! y-coordinates of sample points end type type, private :: params_diag_lineout_t !! Parameters related to defining lineouts character(len=:), allocatable, private :: lineout_select !! Type of lineout character(len=:), allocatable, private :: lineout_path !! Path to file where lineout params are specified integer, private :: n_lineouts !! Number of lineouts integer, dimension(:), allocatable, private :: n_points !! Number of points per lineout to interpolate integer, dimension(:), allocatable, private :: int_order !! Interpolation order for computing lineout values type(lineout_coords_t), dimension(:), allocatable, private :: coords !! X/Y coordinates for each lineout character(len=256), dimension(:), allocatable, private :: lineout_names !! Name of lineout contains procedure, public :: read_params => read_params_lineout procedure, public :: write_params => write_params_lineout end type interface module subroutine read_params_lineout(self, filepath) !! Reads parameters related to lineouts class(params_diag_lineout_t), intent(inout) :: self !! Instance of the type character(len=*), intent(in) :: filepath !! Filepath to read from end subroutine module subroutine write_params_lineout(self, filepath) !! Writes parameters related to lineouts class(params_diag_lineout_t), intent(in) :: self !! Instance of the type character(len=*), intent(in), optional :: filepath !! If present writes to file, otherwise to stdout end subroutine end interface type, public, extends(diagnostics_group_t) :: & diagnostics_lineout_t !! Class containing lineout diagnostics for the BRAGINSKII model type(params_diag_lineout_t), private :: params !! Parameters for lineout diagnostics contains procedure, public :: init_lineout procedure, public :: process_to_file procedure, public :: project_lineout end type contains subroutine init_lineout(self, comm_handler, groupname, filepath) !! Initialize lineout diagnostic class(diagnostics_lineout_t), intent(inout) :: self !! Instance of class type(comm_handler_t), intent(in) :: comm_handler !! MPI communication handler character(len=*), intent(in) :: groupname !! Directory to write snapsfiles into character(len=*), intent(in) :: filepath !! Filepath to read lineouts parameters from ! Read parameters call self%params%read_params(filepath) ! Display parameters if (is_rank_info_writer) then call self%params%write_params() endif ! Initiate base module call self%init_diagnostics_group(comm_handler, groupname) end subroutine 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_lineout_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 :: ind !! Reference index of last initialized diagnostic integer :: i_lineout, k call perf_start('diagnostics_braginskii_lineout') call self%allocate_diags(n_diags = 1 + 3*self%params%n_lineouts & ! for time and lineout coords + 7 * self%params%n_lineouts) ! << ADJUST WHEN ADDING DIAGNOSTICS ! Time call self%init_next_diagnostic(1, 'dim_scalar', 'tau', ind) self%diags(ind)%vals(1) = tau do i_lineout = 1, self%params%n_lineouts ! Lineout coordinates call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='x_'//self%params%lineout_names(i_lineout), & ind_out=ind ) self%diags(ind)%vals = self%params%coords(i_lineout)%x call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='y_'//self%params%lineout_names(i_lineout), & ind_out=ind ) self%diags(ind)%vals = self%params%coords(i_lineout)%y ! Add rho of lineout call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='rho_'//self%params%lineout_names(i_lineout), & ind_out=ind ) do k = 1, self%params%n_points(i_lineout) self%diags(ind)%vals(k) = equi%rho(self%params%coords(i_lineout)%x(k), & self%params%coords(i_lineout)%y(k), & mesh_cano%get_phi()) enddo ! Actual diagnostic variables ! Plasma variables call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='ne_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & ne%vals, mode='TORAVG') call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='te_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & te%vals, mode='TORAVG') call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='ti_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & ti%vals, mode='TORAVG') ! Parallel heat fluxes call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='qpar_cond_e_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & diag_packet%qpar_cond_e, mode='TORAVG') call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='qpar_cond_i_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & diag_packet%qpar_cond_i, mode='TORAVG') call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='qpar_conv_e_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & diag_packet%qpar_conv_e, mode='TORAVG') call self%init_next_diagnostic( & ndim=self%params%n_points(i_lineout), & dimname='dim_'//self%params%lineout_names(i_lineout), & varname='qpar_conv_i_toravg_'//self%params%lineout_names(i_lineout), & ind_out=ind ) call self%project_lineout(comm_handler, mesh_cano, i_lineout, ind, & diag_packet%qpar_conv_i, mode='TORAVG') enddo call self%write_diagnostics(tau, isnaps, idiag, start_new_file) deallocate(self%diags) call perf_stop('diagnostics_braginskii_lineout') end subroutine subroutine project_lineout(self, comm_handler, mesh, i_lineout, ind, u, mode) !! Project input onto lineout diagnostic dimension nrho class(diagnostics_lineout_t), intent(inout) :: self !! Instance of class type(comm_handler_t), intent(in) :: comm_handler !! MPI communication handler type(mesh_cart_t), intent(in) :: mesh !! Mesh within poloidal plane integer, intent(in) :: i_lineout !! Work index of lineout list integer, intent(in) :: ind !! Work index in diagnostics array real(GP), dimension(:), intent(in) :: u !! Input array character(len=*), intent(in) :: mode !! Projection mode, select from : !! 'TORAVG' - toroidal average at each lineout point integer :: l, comm, nplanes, ierr comm = comm_handler%get_comm_planes() nplanes = comm_handler%get_nplanes() select case(mode) case('TORAVG') call interpolate_lineout(mesh, self%params%n_points(i_lineout), & self%params%int_order(i_lineout), & self%params%coords(i_lineout)%x, & self%params%coords(i_lineout)%y, & u_in=u, u_out=self%diags(ind)%vals) ! Take average over planes call MPI_allreduce(MPI_IN_PLACE, self%diags(ind)%vals, & self%params%n_points(i_lineout), MPI_GP, & MPI_SUM, comm, ierr) !$omp parallel do private(l) do l = 1, self%params%n_points(i_lineout) self%diags(ind)%vals(l) = self%diags(ind)%vals(l) / nplanes enddo !$omp end parallel do case default call handle_error('Invalid mode for lineout projection', & GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end subroutine end module