diagnostics_braginskii_s.f90 Source File


Contents


Source Code

submodule (diagnostics_braginskii_m) diagnostics_braginskii_s
    implicit none

contains

    module subroutine init_diagnostics(self, comm_handler, filepath_lineout, equi, mesh_cano, &
                                mesh_stag, polars_cano, polars_stag, tau, isnaps, &
                                ne, te, ti, pot, upar, jpar, apar, apar_fluct, &
                                neutrals_dens, neutrals_parmom, neutrals_pressure)
        class(diagnostics_braginskii_t), intent(inout) :: self
        type(comm_handler_t), intent(in) :: comm_handler
        character(len=*), intent(in) :: filepath_lineout
        class(equilibrium_t), intent(in) :: equi
        type(mesh_cart_t), intent(in) :: mesh_cano
        type(mesh_cart_t), intent(in) :: mesh_stag
        type(polars_t), intent(in) :: polars_cano
        type(polars_t), intent(in) :: polars_stag
        real(GP), intent(in) :: tau
        integer, intent(in) :: isnaps
        type(variable_t), intent(in) :: ne
        type(variable_t), intent(in) :: te
        type(variable_t), intent(in) :: ti
        type(variable_t), intent(in) :: pot
        type(variable_t), intent(in) :: upar
        type(variable_t), intent(in) :: jpar
        type(variable_t), intent(in) :: apar
        type(variable_t), intent(in) :: apar_fluct
        type(variable_t), intent(in) :: neutrals_dens
        type(variable_t), intent(in) :: neutrals_parmom
        type(variable_t), intent(in) :: neutrals_pressure

        self%idiag = 0

        if ( ( zonal_diagnostics_on .or. lineout_diagnostics_on .or. field_diagnostics_on ) &
              .and. .not. equi%is_axisymmetric() ) then
            call handle_error('For non-axisymmetric equilibria only scalar diagnostics are available', &
                               GRILLIX_ERR_OTHER, __LINE__, __FILE__)
        endif

        if ( scalar_diagnostics_on .or. zonal_diagnostics_on &
             .or. lineout_diagnostics_on .or. field_diagnostics_on) then

            if (is_rank_info_writer) then
                write(get_stdout(),*)'---------------------------------------------'
                write(get_stdout(),*)'Initiating diagnostics ----------------------'
            endif

            if (scalar_diagnostics_on) call self%scalar%init_diagnostics_group(comm_handler, 'scalar')
            if (zonal_diagnostics_on) call self%zonal%init_diagnostics_group(comm_handler, 'zonal')
            if (lineout_diagnostics_on) call self%lineout%init_lineout(comm_handler, 'lineout', filepath_lineout)
            ! Field diagnostics are not initialized here, but during the processing step 

            ! Initiate multistep storage (cano and stag) and fill with initial state
            call self%multistep_storage_cano%init_storage(ndim=mesh_cano%get_n_points(), &
                                                          nstorage=2, &
                                                          nfields=6)
            call self%multistep_storage_cano%shift_storage( &
                                            [ne%vals, te%vals, ti%vals, pot%vals, & 
                                             neutrals_dens%vals, neutrals_pressure%vals])

            call self%multistep_storage_stag%init_storage(ndim=mesh_stag%get_n_points(), &
                                                          nstorage=2, &
                                                          nfields=5)
            call self%multistep_storage_stag%shift_storage( &
                                            [upar%vals, jpar%vals, apar%vals, &
                                             apar_fluct%vals, neutrals_parmom%vals])

            call self%clock%init(major_interval=tau_snaps, &
                                 minor_interval=tau_diags, &
                                 major_tally_offset=isnaps)

        endif

    end subroutine

    module subroutine drive_diagnostics(self, comm_handler, equi, equi_on_cano, equi_on_stag, mesh_cano, mesh_stag, map, &
                                        penalisation_cano, penalisation_stag, polars_cano, polars_stag, &
                                        parflux_utils_cano, parflux_utils_stag, perp_bnd_flux_cano, perp_bnd_flux_stag, &
                                        opsinplane_cano, opsinplane_stag, boundaries, boundaries_neutrals, &
                                        cf_buffer_cano, cf_buffer_stag, cf_diss_cano, cf_diss_stag, isnaps, tau, &
                                        ne, te, ti, pot, vort, upar, jpar, apar, apar_fluct, gstress, &
                                        neutrals_dens, neutrals_parmom, neutrals_pressure, &
                                        src_ne, src_te, src_ti, src_upar, src_vort, &
                                        sources_external)
        class(diagnostics_braginskii_t), intent(inout) :: self
        type(comm_handler_t), intent(in) :: comm_handler
        class(equilibrium_t), intent(inout) :: equi
        class(equilibrium_storage_t), intent(in) :: equi_on_cano
        class(equilibrium_storage_t), intent(in) :: equi_on_stag
        type(mesh_cart_t), intent(in) :: mesh_cano
        type(mesh_cart_t), intent(in) :: mesh_stag
        type(parallel_map_t), intent(in) :: map
        type(penalisation_t), intent(in) :: penalisation_cano
        type(penalisation_t), intent(in) :: penalisation_stag
        type(polars_t), intent(in) :: polars_cano
        type(polars_t), intent(in) :: polars_stag
        type(parallel_target_flux_t), intent(in) :: parflux_utils_cano
        type(parallel_target_flux_t), intent(in) :: parflux_utils_stag
        type(perp_bnd_flux_t), intent(in) :: perp_bnd_flux_cano
        type(perp_bnd_flux_t), intent(in) :: perp_bnd_flux_stag
        type(inplane_operators_t), intent(inout) :: opsinplane_cano
        type(inplane_operators_t), intent(inout) :: opsinplane_stag
        type(boundaries_braginskii_t), intent(in) :: boundaries
        type(boundaries_neutrals_t), intent(in) :: boundaries_neutrals
        real(GP), dimension(mesh_cano%get_n_points()), intent(in) :: cf_buffer_cano
        real(GP), dimension(mesh_stag%get_n_points()), intent(in) :: cf_buffer_stag
        real(GP), dimension(mesh_cano%get_n_points()), intent(in) :: cf_diss_cano
        real(GP), dimension(mesh_stag%get_n_points()), intent(in) :: cf_diss_stag
        integer, intent(in) :: isnaps
        real(GP), intent(in) :: tau
        type(variable_t), intent(inout) :: ne 
        type(variable_t), intent(inout) :: te
        type(variable_t), intent(inout) :: ti
        type(variable_t), intent(inout) :: pot
        type(variable_t), intent(inout) :: vort
        type(variable_t), intent(inout) :: upar
        type(variable_t), intent(inout) :: jpar
        type(variable_t), intent(inout) :: apar
        type(variable_t), intent(inout) :: apar_fluct
        type(gyroviscosity_t), intent(inout) :: gstress
        type(variable_t), intent(inout) :: neutrals_dens
        type(variable_t), intent(inout) :: neutrals_parmom
        type(variable_t), intent(inout) :: neutrals_pressure
        real(GP), dimension(mesh_cano%get_n_points_inner()), intent(in) :: src_ne
        real(GP), dimension(mesh_cano%get_n_points_inner()), intent(in) :: src_te
        real(GP), dimension(mesh_cano%get_n_points_inner()), intent(in) :: src_ti
        real(GP), dimension(mesh_stag%get_n_points_inner()), intent(in) :: src_upar
        real(GP), dimension(mesh_cano%get_n_points_inner()), intent(in) :: src_vort
        type(sources_external_t), intent(inout) :: sources_external

        logical :: process_diags
        ! Flag whether to run diagnostics routines this timestep
        logical :: start_new_file
        ! If true, a new file will be written out

        type(diagnostics_packet_t), allocatable :: diag_packet

        process_diags = .false. 

        if ( scalar_diagnostics_on .or. zonal_diagnostics_on &
            .or. lineout_diagnostics_on .or. field_diagnostics_on ) then
            ! Advance multistep storage at every call
            call self%multistep_storage_cano%shift_storage( &
                        [ne%vals, te%vals, ti%vals, pot%vals, &
                        neutrals_dens%vals, neutrals_pressure%vals])
            call self%multistep_storage_stag%shift_storage( &
                        [upar%vals, jpar%vals, apar%vals, &
                        apar_fluct%vals, neutrals_parmom%vals])
            
            ! Drive clock
            call self%clock%drive(tau, dtau)
            process_diags = self%clock%on_minor_checkpoint()
            start_new_file = self%clock%on_major_checkpoint()
            self%idiag = self%clock%get_minor_tally_bwd()
        endif
        
        if ( analyse_snapshots_on ) then
            process_diags = .true.
            start_new_file = .true.
            self%idiag  = 1
        endif

        ! Only perform computation on timesteps where demanded
        if ( .not. (process_diags) ) then
            return
        endif

        ! otherwise, perform diagnostics

        allocate(diag_packet)

        call diag_packet%compute(comm_handler, equi, equi_on_cano, equi_on_stag, mesh_cano, mesh_stag, &
                            map, penalisation_cano, penalisation_stag, polars_cano, polars_stag, &
                            parflux_utils_cano, parflux_utils_stag, perp_bnd_flux_cano, perp_bnd_flux_stag, &
                            opsinplane_cano, opsinplane_stag, boundaries, boundaries_neutrals, &
                            cf_buffer_cano, cf_buffer_stag, cf_diss_cano, cf_diss_stag, tau, &
                            ne, te, ti, pot, vort, upar, jpar, apar, apar_fluct, gstress, &
                            neutrals_dens, neutrals_parmom, neutrals_pressure, sources_external, &
                            src_upar, src_vort, self%multistep_storage_cano, self%multistep_storage_stag)

        if ( scalar_diagnostics_on ) then
            call self%scalar%process_to_file(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, self%idiag, start_new_file)
        endif
        if ( zonal_diagnostics_on ) then
            call self%zonal%process_to_file(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, &
                                src_ne, src_te, src_ti, src_upar, src_vort, &
                                isnaps, self%idiag, start_new_file)
        endif
        if ( lineout_diagnostics_on ) then
            call self%lineout%process_to_file(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, self%idiag, start_new_file)
        endif
        if ( field_diagnostics_on ) then
            call self%field%process_to_file(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, self%idiag, start_new_file)
        endif

        deallocate(diag_packet)

    end subroutine

end submodule