Initialises a file for saving output of eigensolver to netcdf
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
character(len=*), | intent(in) | :: | fname | |||
logical, | intent(in) | :: | has_phi | |||
logical, | intent(in) | :: | has_apar | |||
logical, | intent(in) | :: | has_bpar | |||
type(EigNetcdfID), | intent(inout) | :: | IDs |
subroutine init_eigenfunc_file(fname, has_phi, has_apar, has_bpar, IDs)
use file_utils, only: error_unit
use theta_grid, only: ntgrid, theta
#ifdef NETCDF
use netcdf_utils, only: ensure_netcdf_dim_exists, ensure_netcdf_var_exists
#endif
implicit none
character(len=*), intent(in) :: fname
type(EigNetcdfID), intent(inout) :: IDs
logical, intent(in) :: has_phi, has_apar, has_bpar
#ifdef NETCDF
integer :: ierr
#endif
!Set nconv counter to 0
IDs%nconv_count=0
#ifdef NETCDF
!/Make file
ierr=nf90_create(fname, create_file_mode, IDs%ncid)
if(ierr/=NF90_NOERR) call netcdf_error(ierr,file=fname)
!/Define dimensions
call ensure_netcdf_dim_exists(IDs%ncid, "ri", 2, IDs%ri_dim_id)
call ensure_netcdf_dim_exists(IDs%ncid, "theta", 2*ntgrid+1, IDs%theta_dim_id)
call ensure_netcdf_dim_exists(IDs%ncid, "nconv", NF90_UNLIMITED, IDs%conv_dim_id)
!/Define variables
!Dimensions
call ensure_netcdf_var_exists(IDs%ncid, "theta", netcdf_real, &
IDs%theta_dim_id, IDs%theta_var_id)
call ensure_netcdf_var_exists(IDs%ncid, "conv", netcdf_real, &
IDs%conv_dim_id, IDs%conv_var_id)
!Fields
if(has_phi)then
call ensure_netcdf_var_exists(IDs%ncid, "phi", netcdf_real, &
[IDs%ri_dim_id, IDs%theta_dim_id, IDs%conv_dim_id], IDs%phi_var_id)
endif
if(has_apar)then
call ensure_netcdf_var_exists(IDs%ncid, "apar", netcdf_real, &
[IDs%ri_dim_id, IDs%theta_dim_id, IDs%conv_dim_id], IDs%apar_var_id)
endif
if(has_bpar)then
call ensure_netcdf_var_exists(IDs%ncid, "bpar", netcdf_real, &
[IDs%ri_dim_id, IDs%theta_dim_id, IDs%conv_dim_id], IDs%bpar_var_id)
endif
!Omega
call ensure_netcdf_var_exists(IDs%ncid, "omega", netcdf_real, &
[IDs%ri_dim_id, IDs%conv_dim_id], IDs%omega_var_id)
!End definitions
ierr=nf90_enddef(IDs%ncid)
if(ierr/=NF90_NOERR) call netcdf_error(ierr,file=fname)
!Now we can place some data in the file
ierr=nf90_put_var(IDs%ncid,IDs%theta_var_id,theta)
if(ierr/=NF90_NOERR) call netcdf_error(ierr,var="theta")
#endif
end subroutine init_eigenfunc_file