Write various moments to netcdf
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | file_id |
NetCDF ID of the file to write to |
||
integer, | intent(in) | :: | nout |
Current output timestep |
||
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | ntot |
Total perturbed density |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | density |
Non-adiabatic part of the perturbed density |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | upar |
Parallel velocity |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | tpar |
Parallel temperature |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | tperp |
Perpendicular temperature |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | qparflux |
Parallel heat flux |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | pperpj1 |
Pressure part of the particle flux |
|
complex, | intent(in), | dimension (-ntgrid:,:,:,:) | :: | qpperpj1 |
Pressure part of the heat flux |
|
logical, | intent(in), | optional | :: | ob_midplane |
If true, write out moments at the outboard midplane only, otherwise write them out along the entire flux surface |
subroutine nc_write_moments (file_id, nout, ntot, density, upar, tpar, tperp, qparflux, pperpj1, qpperpj1, ob_midplane)
use theta_grid, only: ntgrid
implicit none
!> NetCDF ID of the file to write to
integer, intent(in) :: file_id
!> Total perturbed density
complex, dimension (-ntgrid:,:,:,:), intent (in) :: ntot
!> Non-adiabatic part of the perturbed density
complex, dimension (-ntgrid:,:,:,:), intent (in) :: density
!> Parallel velocity
complex, dimension (-ntgrid:,:,:,:), intent (in) :: upar
!> Parallel temperature
complex, dimension (-ntgrid:,:,:,:), intent (in) :: tpar
!> Perpendicular temperature
complex, dimension (-ntgrid:,:,:,:), intent (in) :: tperp
!> Parallel heat flux
complex, dimension (-ntgrid:,:,:,:), intent (in) :: qparflux
!> Pressure part of the particle flux
complex, dimension (-ntgrid:,:,:,:), intent (in) :: pperpj1
!> Pressure part of the heat flux
complex, dimension (-ntgrid:,:,:,:), intent (in) :: qpperpj1
!> Current output timestep
integer, intent (in) :: nout
!> If true, write out moments at the outboard midplane only,
!> otherwise write them out along the entire flux surface
logical, optional, intent(in) :: ob_midplane
#ifdef NETCDF
if (present(ob_midplane)) then
if (ob_midplane ) then
call netcdf_write_complex(file_id, "ntot0", ntot(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Total perturbed density (theta=0) over time")
call netcdf_write_complex(file_id, "density0", density(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Nonadiabatic piece of perturbed density (theta=0) over time")
call netcdf_write_complex(file_id, "upar0", upar(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Upar (theta=0) over time")
call netcdf_write_complex(file_id, "tpar0", tpar(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Tpar (theta=0) over time")
call netcdf_write_complex(file_id, "tperp0", tperp(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Tperp (theta=0) over time")
call netcdf_write_complex(file_id, "qparflux0", qparflux(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Qparflux (theta=0) over time")
call netcdf_write_complex(file_id, "pperpj10", pperpj1(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Pperpj1 (theta=0) over time")
call netcdf_write_complex(file_id, "qpperpj10", qpperpj1(0, :, :, :), mom_dims, start=starts(5, nout), &
long_name="Qpperpj1 (theta=0) over time")
return
endif
endif
call netcdf_write_complex(file_id, 'ntot_t', ntot, mom_t_dims, start=starts(6, nout), &
long_name="Total perturbed density over time")
call netcdf_write_complex(file_id, 'density_t', density, mom_t_dims, start=starts(6, nout), &
long_name="Nonadiabatic piece of perturbed density over time")
call netcdf_write_complex(file_id, 'upar_t', upar, mom_t_dims, start=starts(6, nout), &
long_name="Upar over time")
call netcdf_write_complex(file_id, 'tpar_t', tpar, mom_t_dims, start=starts(6, nout), &
long_name="Tpar over time")
call netcdf_write_complex(file_id, 'tperp_t', tperp, mom_t_dims, start=starts(6, nout), &
long_name="Tperp over time")
call netcdf_write_complex(file_id, 'qparflux_t', qparflux, mom_t_dims, start=starts(6, nout), &
long_name="Qparflux over time")
call netcdf_write_complex(file_id, 'pperpj1_t', pperpj1, mom_t_dims, start=starts(6, nout), &
long_name="Pperpj1 over time")
call netcdf_write_complex(file_id, 'qpperpj1_t', qpperpj1, mom_t_dims, start=starts(6, nout), &
long_name="Qpperpj1 over time")
# endif
end subroutine nc_write_moments