!> A module for writing out moments of the distribution function, !! for example density, temperature, parallel flow module diagnostics_moments implicit none private public :: write_moments, write_full_moments_notgc contains !> FIXME : Add documentation subroutine write_moments(gnostics) use theta_grid, only: ntgrid use kt_grids, only: naky, ntheta0 use species, only: nspec use dist_fn, only: getmoms use fields_arrays, only: phinew, bparnew use diagnostics_config, only: diagnostics_type implicit none type(diagnostics_type), intent(in) :: gnostics complex, dimension (-ntgrid:ntgrid,ntheta0,naky,nspec) :: ntot, density, & upar, tpar, tperp, qparflux, pperpj1, qpperpj1 call getmoms (phinew, bparnew, ntot, density, upar, tpar, tperp, qparflux, pperpj1, qpperpj1) call write_standard_moment_properties(gnostics, & 'ntot', 'The total perturbed species density', 'n_ref', ntot, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'density', 'The non-adiabatic part of the perturbed species density', 'n_ref', density, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'upar', 'The perturbed parallel flow', 'v_ths', upar, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'tpar', 'The perturbed parallel temperature', 'T_ref', tpar, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'tperp', 'The perturbed perpendicular temperature', 'T_ref', tperp, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'qparflux', 'The parallel heat flux', 'n_ref T_ref v_thr', qparflux, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'pperpj1', 'A modified perpendicular pressure which give the particle flux from bpar',& 'n_ref T_ref / q_ref', pperpj1, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'qpperpj1', 'A modified perpendicular pressure * energy which give the heat flux from bpar',& 'n_ref T_ref^2 / q_ref', qpperpj1, gnostics%distributed) end subroutine write_moments !> FIXME : Add documentation subroutine write_full_moments_notgc(gnostics) use theta_grid, only: ntgrid use kt_grids, only: naky, ntheta0 use species, only: nspec use dist_fn, only: getmoms_notgc use fields_arrays, only: phinew, bparnew use diagnostics_config, only: diagnostics_type implicit none type(diagnostics_type), intent(in) :: gnostics complex, dimension (-ntgrid:ntgrid,ntheta0,naky,nspec) :: ntot, density, upar, tpar, tperp call getmoms_notgc(phinew, bparnew, density,upar,tpar,tperp,ntot) call write_standard_moment_properties(gnostics, & 'ntot_notgc', 'The total perturbed species density & & in non-guiding centre coordinates', 'n_r', ntot, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'density_notgc', 'The non-adiabatic part of the perturbed species density & & in non-guiding centre coordinates', 'n_r', density, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'upar_notgc', 'The perturbed parallel flow & & in non-guiding centre coordinates', 'v_ths', upar, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'tpar_notgc', 'The perturbed parallel temperature & & in non-guiding centre coordinates', 'T_r', tpar, gnostics%distributed) call write_standard_moment_properties(gnostics, & 'tperp_notgc', 'The perturbed perpendicular temperature & & in non-guiding centre coordinates', 'T_r', tperp, gnostics%distributed) end subroutine write_full_moments_notgc !> FIXME : Add documentation subroutine write_standard_moment_properties(gnostics, moment_name, moment_description, & moment_units, moment_value, distributed) use volume_averages, only: average_theta, average_kx, average_ky use fields_parallelization, only: field_k_local use mp, only: sum_allreduce use theta_grid, only: ntgrid use kt_grids, only: ntheta0, naky use species, only: nspec use diagnostics_config, only: diagnostics_type use gs2_io, only: starts, mom_dims, mom_t_dims, fluxk_dims, netcdf_write_complex, complex_dim, & kx_dim, ky_dim, species_dim, time_dim, theta_dim use neasyf, only: neasyf_write implicit none type(diagnostics_type), intent(in) :: gnostics character(*), intent(in) :: moment_name, moment_description, moment_units complex, dimension(-ntgrid:,:,:,:), intent(inout) :: moment_value logical, intent(in) :: distributed real, dimension(ntheta0, naky, nspec) :: moment2_by_mode real, dimension(naky,nspec) :: moment2_by_ky real, dimension(ntheta0,nspec) :: moment2_by_kx complex, dimension(ntheta0, naky, nspec) :: moment_by_mode complex, dimension(ntheta0, naky, nspec) :: moment_igomega_by_mode complex, dimension(ntheta0, nspec) :: moment_flx_surfavg logical :: write_moment_by_time integer :: it !, ik if (.not. gnostics%writing) return call average_theta(moment_value, moment_value, moment2_by_mode, distributed) call average_kx(moment2_by_mode, moment2_by_ky, distributed) call neasyf_write(gnostics%file_id, moment_name//"2_by_ky", moment2_by_ky, & dim_names=[ky_dim, species_dim, time_dim], start=starts(3, gnostics%nout), & long_name=moment_description//" squared and averaged over theta and kx", & units="("//moment_units//")^2") call average_ky(moment2_by_mode, moment2_by_kx, distributed) call neasyf_write(gnostics%file_id, moment_name//"2_by_kx", moment2_by_kx, & dim_names=[kx_dim, species_dim, time_dim], start=starts(3, gnostics%nout), & long_name=moment_description//" squared and averaged over theta and ky", & units="("//moment_units//")^2") call neasyf_write(gnostics%file_id, moment_name//"2", sum(moment2_by_kx, 1), & dim_names=[species_dim, time_dim], start=starts(2, gnostics%nout), & long_name=moment_description//" squared and averaged over theta, kx and ky", & units="("//moment_units//")^2") call average_theta(moment_value, moment_by_mode, gnostics%distributed) moment_igomega_by_mode(:,:,:) = moment_value(gnostics%igomega,:,:,:) ! moment_by_mode could be distributed so we have to be careful here moment_flx_surfavg(:,:) = 0.0 do it = 1,ntheta0 if ((.not. gnostics%distributed).or.field_k_local(it, 1))& moment_flx_surfavg(it, :) = moment_by_mode(it, 1, :) end do if (gnostics%distributed) call sum_allreduce(moment_flx_surfavg) call netcdf_write_complex(gnostics%file_id, moment_name//"_flxsurf_avg", moment_flx_surfavg, & dim_names=[complex_dim, kx_dim, species_dim, time_dim], start=starts(4, gnostics%nout), & long_name=moment_description//" flux surface averaged over theta, at ky=0 (actally ik==1)", & units=moment_units) call netcdf_write_complex(gnostics%file_id, moment_name, moment_value, & dim_names=[complex_dim, theta_dim, kx_dim, ky_dim, species_dim], & long_name=moment_description, units=moment_units) call netcdf_write_complex( & gnostics%file_id, moment_name//"_igomega_by_mode", moment_igomega_by_mode, & dim_names=mom_dims, start=starts(5, gnostics%nout), & long_name=moment_description//" at ig=igomega" , units=moment_units) call netcdf_write_complex(gnostics%file_id, moment_name//"_by_mode", moment_by_mode, & dim_names=mom_dims, start=starts(5, gnostics%nout), & long_name=moment_description//" averaged over theta" , & units=moment_units) call neasyf_write(gnostics%file_id, moment_name//"2_by_mode", moment2_by_mode, & dim_names=fluxk_dims, start=starts(4, gnostics%nout), & long_name=moment_description//" squared and averaged over theta" , & units=moment_units) write_moment_by_time = .false. if (moment_name .eq. 'ntot' .and. gnostics%write_ntot_over_time ) write_moment_by_time = .true. if (moment_name .eq. 'density' .and. gnostics%write_density_over_time) write_moment_by_time = .true. if (moment_name .eq. 'upar' .and. gnostics%write_upar_over_time) write_moment_by_time = .true. if (moment_name .eq. 'tperp' .and. gnostics%write_tperp_over_time) write_moment_by_time = .true. if (write_moment_by_time) then call netcdf_write_complex(gnostics%file_id, moment_name//"_t", moment_value, & dim_names=mom_t_dims, start=starts(6, gnostics%nout), & long_name=moment_description//": the whole moment" , & units=moment_units) end if end subroutine write_standard_moment_properties end module diagnostics_moments