MODULE iostart_pem !----------------------------------------------------------------------- ! NAME ! iostart_pem ! ! DESCRIPTION ! Read and write specific netcdf files for the PEM (start and restart) ! Provides interfaces for field and variable I/O operations ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! Inspired by iostart from the PCM. ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! MODULE VARIABLES ! ---------------- PRIVATE INTEGER,SAVE :: nid_start ! NetCDF file identifier for startfi.nc file INTEGER,SAVE :: nid_restart ! NetCDF file identifier for restartfi.nc file ! restartfi.nc file dimension identifiers: (see open_restartphy()) INTEGER,SAVE :: idim1 ! "index" dimension INTEGER,SAVE :: idim2 ! "physical_points" dimension INTEGER,SAVE :: idim3 ! "subsurface_layers" dimension INTEGER,SAVE :: idim4 ! "nlayer_plus_1" dimension INTEGER,SAVE :: idim5 ! "number_of_advected_fields" dimension INTEGER,SAVE :: idim6 ! "nlayer" dimension INTEGER,SAVE :: idim7 ! "Time" dimension INTEGER,SAVE :: idim8 ! "nslope" dimension INTEGER,SAVE :: idim9 ! "inter slope" dimension INTEGER,SAVE :: idim10 ! "nb_str_max" dimension INTEGER,SAVE :: timeindex ! current time index (for time-dependent fields) INTEGER,PARAMETER :: length=100 ! size of tab_cntrl array INTERFACE get_field MODULE PROCEDURE Get_field_r1,Get_field_r2,Get_field_r3 END INTERFACE get_field INTERFACE get_var MODULE PROCEDURE get_var_r0,Get_var_r1,Get_var_r2,Get_var_r3 END INTERFACE get_var INTERFACE put_field MODULE PROCEDURE put_field_r1,put_field_r2,put_field_r3 END INTERFACE put_field INTERFACE put_var MODULE PROCEDURE put_var_r0,put_var_r1,put_var_r2,put_var_r3 END INTERFACE put_var PUBLIC nid_start, length PUBLIC get_field,get_var,put_field,put_var PUBLIC inquire_dimension, inquire_dimension_length PUBLIC inquire_field, inquire_field_ndims PUBLIC open_startphy,close_startphy,open_restartphy,close_restartphy CONTAINS !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !======================================================================= SUBROUTINE open_startphy(filename) !----------------------------------------------------------------------- ! NAME ! open_startphy ! ! DESCRIPTION ! Open the startphy netCDF file for reading. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_OPEN, NF90_NOERR, NF90_NOWRITE, nf90_strerror USE mod_phys_lmdz_para, only: is_master, bcast ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*) :: filename ! LOCAL VARIABLES ! --------------- INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr = NF90_OPEN (filename, NF90_NOWRITE, nid_start) IF (ierr.NE.NF90_NOERR) THEN write(*,*)'open_startphy: problem opening file '//trim(filename) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_startphy","Cannot open file",1) ENDIF ENDIF CALL bcast(nid_start) ! tell all procs about nid_start END SUBROUTINE open_startphy !======================================================================= !======================================================================= SUBROUTINE close_startphy !----------------------------------------------------------------------- ! NAME ! close_startphy ! ! DESCRIPTION ! Close the startphy netCDF file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_CLOSE USE mod_phys_lmdz_para, only: is_master ! DECLARATION ! ----------- implicit none ! LOCAL VARIABLES ! --------------- INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr = NF90_CLOSE (nid_start) ENDIF END SUBROUTINE close_startphy !======================================================================= !======================================================================= FUNCTION inquire_field(Field_name) !----------------------------------------------------------------------- ! NAME ! inquire_field ! ! DESCRIPTION ! Check if a given field is present in the input file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_INQ_VARID, NF90_NOERR USE mod_phys_lmdz_para, only: is_master, bcast ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: Field_name ! LOCAL VARIABLES ! --------------- LOGICAL :: inquire_field INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr=NF90_INQ_VARID(nid_start,Field_name,varid) IF (ierr==NF90_NOERR) THEN Inquire_field=.TRUE. ELSE Inquire_field=.FALSE. ENDIF ENDIF CALL bcast(inquire_field) END FUNCTION inquire_field !======================================================================= !======================================================================= FUNCTION inquire_field_ndims(Field_name) !----------------------------------------------------------------------- ! NAME ! inquire_field_ndims ! ! DESCRIPTION ! Give the number of dimensions of "Field_name" stored in the input file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: nf90_inq_varid, nf90_inquire_variable, & NF90_NOERR, nf90_strerror USE mod_phys_lmdz_para, only: is_master, bcast ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: Field_name ! LOCAL VARIABLES ! --------------- INTEGER :: inquire_field_ndims INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr=nf90_inq_varid(nid_start,Field_name,varid) ierr=nf90_inquire_variable(nid_start,varid,& ndims=inquire_field_ndims) IF (ierr.NE.NF90_NOERR) THEN write(*,*)'inquire_field_ndims: problem obtaining ndims of '& //trim(field_name) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("inquire_field_ndims","Failed to get ndims",1) ENDIF ENDIF CALL bcast(inquire_field_ndims) END FUNCTION inquire_field_ndims !======================================================================= !======================================================================= FUNCTION inquire_dimension(Field_name) !----------------------------------------------------------------------- ! NAME ! inquire_dimension ! ! DESCRIPTION ! Check if a given dimension is present in the input file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: nf90_inq_dimid, NF90_NOERR USE mod_phys_lmdz_para, only: is_master, bcast ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: Field_name ! LOCAL VARIABLES ! --------------- LOGICAL :: inquire_dimension INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr=NF90_INQ_DIMID(nid_start,Field_name,varid) IF (ierr==NF90_NOERR) THEN Inquire_dimension=.TRUE. ELSE Inquire_dimension=.FALSE. ENDIF ENDIF CALL bcast(inquire_dimension) END FUNCTION inquire_dimension !======================================================================= !======================================================================= FUNCTION inquire_dimension_length(Field_name) !----------------------------------------------------------------------- ! NAME ! inquire_dimension_length ! ! DESCRIPTION ! Give the length of the "Field_name" dimension stored in the input file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: nf90_inquire_dimension, nf90_inq_dimid, & NF90_NOERR, nf90_strerror USE mod_phys_lmdz_para, only: is_master, bcast ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: Field_name ! LOCAL VARIABLES ! --------------- INTEGER :: inquire_dimension_length INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr=nf90_inq_dimid(nid_start,Field_name,varid) ierr=nf90_inquire_dimension(nid_start,varid,& len=inquire_dimension_length) IF (ierr.NE.NF90_NOERR) THEN write(*,*)'inquire_field_length: problem obtaining length of '& //trim(field_name) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("inquire_field_ndims","Failed to get length",1) ENDIF ENDIF CALL bcast(inquire_dimension_length) END FUNCTION inquire_dimension_length !======================================================================= !======================================================================= SUBROUTINE Get_Field_r1(field_name,field,found,timeindex) !----------------------------------------------------------------------- ! NAME ! Get_Field_r1 ! ! DESCRIPTION ! Read a surface field (1D) from the start file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use mod_grid_phy_lmdz, only: klon_glo ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name REAL, INTENT(INOUT) :: field(:) LOGICAL, INTENT(OUT), OPTIONAL :: found INTEGER, INTENT(IN), OPTIONAL :: timeindex ! LOCAL VARIABLES ! --------------- integer :: edges(4), corners(4) ! CODE ! ---- edges(1)=klon_glo edges(2:4)=1 corners(1:4)=1 if (PRESENT(timeindex)) then corners(2)=timeindex endif IF (PRESENT(found)) THEN CALL Get_field_rgen(field_name,field,1,corners,edges,found) ELSE CALL Get_field_rgen(field_name,field,1,corners,edges) ENDIF END SUBROUTINE Get_Field_r1 !======================================================================= !======================================================================= SUBROUTINE Get_Field_r2(field_name,field,found,timeindex) !----------------------------------------------------------------------- ! NAME ! Get_Field_r2 ! ! DESCRIPTION ! Read a "3D" horizontal-vertical field (2D) from the start file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use mod_grid_phy_lmdz, only: klon_glo ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name REAL, INTENT(INOUT) :: field(:,:) LOGICAL, INTENT(OUT), OPTIONAL :: found INTEGER, INTENT(IN), OPTIONAL :: timeindex ! LOCAL VARIABLES ! --------------- integer :: edges(4), corners(4) ! CODE ! ---- edges(1)=klon_glo edges(2)=size(field,2) edges(3:4)=1 corners(1:4)=1 if (PRESENT(timeindex)) then corners(3)=timeindex endif IF (PRESENT(found)) THEN CALL Get_field_rgen(field_name,field,size(field,2),& corners,edges,found) ELSE CALL Get_field_rgen(field_name,field,size(field,2),& corners,edges) ENDIF END SUBROUTINE Get_Field_r2 !======================================================================= !======================================================================= SUBROUTINE Get_Field_r3(field_name,field,found,timeindex) !----------------------------------------------------------------------- ! NAME ! Get_Field_r3 ! ! DESCRIPTION ! Read a "4D" field surf/alt/?? (3D) from the start file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use mod_grid_phy_lmdz, only: klon_glo ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name REAL, INTENT(INOUT) :: field(:,:,:) LOGICAL, INTENT(OUT), OPTIONAL :: found INTEGER, INTENT(IN), OPTIONAL :: timeindex ! LOCAL VARIABLES ! --------------- integer :: edges(4), corners(4) ! CODE ! ---- edges(1)=klon_glo edges(2)=size(field,2) edges(3)=size(field,3) edges(4)=1 corners(1:4)=1 if (PRESENT(timeindex)) then corners(4)=timeindex endif IF (PRESENT(found)) THEN CALL Get_field_rgen(field_name,field,size(field,2)*size(field,3),& corners,edges,found) ELSE CALL Get_field_rgen(field_name,field,size(field,2)*size(field,3),& corners,edges) ENDIF END SUBROUTINE Get_Field_r3 !======================================================================= !======================================================================= SUBROUTINE Get_field_rgen(field_name,field,field_size, & corners,edges,found) !----------------------------------------------------------------------- ! NAME ! Get_field_rgen ! ! DESCRIPTION ! Generic subroutine to read fields from start file with scatter to processes. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf USE dimphy USE mod_grid_phy_lmdz USE mod_phys_lmdz_para ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name INTEGER, INTENT(IN) :: field_size REAL, INTENT(OUT) :: field(klon,field_size) INTEGER, INTENT(IN) :: corners(4) INTEGER, INTENT(IN) :: edges(4) LOGICAL, INTENT(OUT), OPTIONAL :: found ! LOCAL VARIABLES ! --------------- REAL :: field_glo(klon_glo,field_size) LOGICAL :: tmp_found INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr=NF90_INQ_VARID(nid_start,field_name,varid) IF (ierr==NF90_NOERR) THEN CALL body(field_glo) tmp_found=.TRUE. ELSE tmp_found=.FALSE. ENDIF ENDIF CALL bcast(tmp_found) IF (tmp_found) THEN CALL scatter(field_glo,field) ENDIF IF (PRESENT(found)) THEN found=tmp_found ELSE IF (.NOT. tmp_found) THEN write(*,*) 'get_field_rgen: Field <'//field_name//'> not found' CALL abort_physic("get_field_rgen","Failed to get field",1) ENDIF ENDIF CONTAINS !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !======================================================================= SUBROUTINE body(field_glo) !----------------------------------------------------------------------- ! NAME ! body ! ! DESCRIPTION ! Read the requested field from the start file into the global array. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! Helper contained in Get_field_rgen. !----------------------------------------------------------------------- ! DECLARATION ! ----------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- REAL, INTENT(OUT) :: field_glo(klon_glo*field_size) ! Buffer for global field ! CODE ! ---- ierr=NF90_GET_VAR(nid_start,varid,field_glo,corners,edges) IF (ierr/=NF90_NOERR) THEN ! La variable exist dans le fichier mais la lecture a echouee. write(*,*) 'get_field_rgen: Failed reading <'//field_name//'>' ! IF (field_name=='CLWCON' .OR. field_name=='RNEBCON' .OR. field_name=='RATQS') THEN ! ! Essaye de lire le variable sur surface uniqument, comme fait avant ! field_glo(:)=0. ! ierr=NF90_GET_VAR(nid_start,varid,field_glo(1:klon_glo)) ! IF (ierr/=NF90_NOERR) THEN ! write(*,*) 'phyetat0: Lecture echouee aussi en 2D pour <'//field_name//'>' ! CALL stop_clean ! ELSE ! write(*,*) 'phyetat0: La variable <'//field_name//'> lu sur surface seulement'!, selon ancien format, le reste mis a zero' ! END IF ! ELSE CALL abort_physic("get_field_rgen","Failed to read field",1) ! ENDIF ENDIF END SUBROUTINE body END SUBROUTINE Get_field_rgen !======================================================================= !======================================================================= SUBROUTINE get_var_r0(var_name,var,found) !----------------------------------------------------------------------- ! NAME ! get_var_r0 ! ! DESCRIPTION ! Get a scalar variable from input file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name REAL, INTENT(INOUT) :: var LOGICAL, INTENT(OUT), OPTIONAL :: found ! LOCAL VARIABLES ! --------------- REAL :: varout(1) ! CODE ! ---- IF (PRESENT(found)) THEN CALL Get_var_rgen(var_name,varout,size(varout),found) ELSE CALL Get_var_rgen(var_name,varout,size(varout)) ENDIF var=varout(1) END SUBROUTINE get_var_r0 !======================================================================= !======================================================================= SUBROUTINE get_var_r1(var_name,var,found) !----------------------------------------------------------------------- ! NAME ! get_var_r1 ! ! DESCRIPTION ! Get a vector (1D) variable from input file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name REAL, INTENT(INOUT) :: var(:) LOGICAL, INTENT(OUT), OPTIONAL :: found ! CODE ! ---- IF (PRESENT(found)) THEN CALL Get_var_rgen(var_name,var,size(var),found) ELSE CALL Get_var_rgen(var_name,var,size(var)) ENDIF END SUBROUTINE get_var_r1 !======================================================================= !======================================================================= SUBROUTINE get_var_r2(var_name,var,found) !----------------------------------------------------------------------- ! NAME ! get_var_r2 ! ! DESCRIPTION ! Get a 2D field variable from input file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name REAL, INTENT(OUT) :: var(:,:) LOGICAL, INTENT(OUT), OPTIONAL :: found ! CODE ! ---- IF (PRESENT(found)) THEN CALL Get_var_rgen(var_name,var,size(var),found) ELSE CALL Get_var_rgen(var_name,var,size(var)) ENDIF END SUBROUTINE get_var_r2 !======================================================================= !======================================================================= SUBROUTINE get_var_r3(var_name,var,found) !----------------------------------------------------------------------- ! NAME ! get_var_r3 ! ! DESCRIPTION ! Get a 3D field variable from input file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name REAL, INTENT(INOUT) :: var(:,:,:) LOGICAL, INTENT(OUT), OPTIONAL :: found ! CODE ! ---- IF (PRESENT(found)) THEN CALL Get_var_rgen(var_name,var,size(var),found) ELSE CALL Get_var_rgen(var_name,var,size(var)) ENDIF END SUBROUTINE get_var_r3 !======================================================================= !======================================================================= SUBROUTINE Get_var_rgen(var_name,var,var_size,found) !----------------------------------------------------------------------- ! NAME ! Get_var_rgen ! ! DESCRIPTION ! Generic subroutine to read variables from input file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf USE dimphy USE mod_grid_phy_lmdz USE mod_phys_lmdz_para ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name INTEGER, INTENT(IN) :: var_size REAL, INTENT(OUT) :: var(var_size) LOGICAL, INTENT(OUT), OPTIONAL :: found ! LOCAL VARIABLES ! --------------- LOGICAL :: tmp_found INTEGER :: varid INTEGER :: ierr ! CODE ! ---- IF (is_mpi_root .AND. is_omp_root) THEN ierr=NF90_INQ_VARID(nid_start,var_name,varid) IF (ierr==NF90_NOERR) THEN ierr=NF90_GET_VAR(nid_start,varid,var) IF (ierr/=NF90_NOERR) THEN write(*,*) 'phyetat0: Failed loading <'//trim(var_name)//'>' CALL abort_physic("get_var_rgen","Failed to read variable",1) ENDIF tmp_found=.TRUE. ELSE tmp_found=.FALSE. ENDIF ENDIF CALL bcast(tmp_found) IF (tmp_found) THEN CALL bcast(var) ENDIF IF (PRESENT(found)) THEN found=tmp_found ELSE IF (.NOT. tmp_found) THEN write(*,*) 'phyetat0: Variable <'//trim(var_name)//'> not found' CALL abort_physic("get_var_rgen","Failed to read variable",1) ENDIF ENDIF END SUBROUTINE Get_var_rgen !======================================================================= !======================================================================= SUBROUTINE open_restartphy(filename) !----------------------------------------------------------------------- ! NAME ! open_restartphy ! ! DESCRIPTION ! Open (or create) the restart netCDF file for writing. ! Define dimensions and global attributes on first call. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_CREATE, NF90_CLOBBER, NF90_64BIT_OFFSET, & NF90_NOERR, nf90_strerror, & NF90_PUT_ATT, NF90_GLOBAL, NF90_DEF_DIM, & NF90_UNLIMITED, NF90_ENDDEF, & NF90_WRITE, NF90_OPEN USE mod_phys_lmdz_para, only: is_master USE mod_grid_phy_lmdz, only: klon_glo USE dimphy, only: klev, klevp1 USE tracer_mod, only: nqmx USE soil, only: nsoilmx_PEM USE comslope_mod, only: nslope use layered_deposits, only: nb_str_max ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: filename ! LOCAL VARIABLES ! --------------- INTEGER :: ierr LOGICAL,SAVE :: already_created=.false. ! CODE ! ---- IF (is_master) THEN if (.not.already_created) then ! At the very first call, create the file ierr=NF90_CREATE(filename,IOR(NF90_CLOBBER,NF90_64BIT_OFFSET), & nid_restart) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem creating file '//trim(filename) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed creating file",1) ENDIF already_created=.true. else ! Just open the file ierr=NF90_OPEN(filename,NF90_WRITE,nid_restart) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem opening file '//trim(filename) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed opening file",1) ENDIF return endif ! of if (.not.already_created) ierr=NF90_PUT_ATT(nid_restart,NF90_GLOBAL,"title",& "Physics start file") IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem writing title ' write(*,*)trim(nf90_strerror(ierr)) ENDIF ierr=NF90_DEF_DIM(nid_restart,"index",length,idim1) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining index dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining index",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"physical_points",klon_glo,idim2) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining physical_points dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining physical_points",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"subsurface_layers",nsoilmx_PEM,idim3) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining subsurface_layers dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining subsurface_layers",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"nlayer_plus_1",klevp1,idim4) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining nlayer_plus_1 dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining nlayer_plus_1",1) ENDIF if (nqmx.ne.0) then ierr=NF90_DEF_DIM(nid_restart,"number_of_advected_fields",nqmx,idim5) else ! pretend nqmx=1 because 0 size implies unlimited dimension for NetCDF ierr=NF90_DEF_DIM(nid_restart,"number_of_advected_fields",1,idim5) endif IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining number_of_advected_fields dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining number_of_advected_fields",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"nlayer",klev,idim6) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining nlayer dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining nlayer",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"Time",NF90_UNLIMITED,idim7) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem defining Time dimension ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed defining Time",1) ENDIF ierr=NF90_DEF_DIM(nid_restart,"nslope",nslope,idim8) IF (ierr/=NF90_NOERR) THEN write(*,*)'phyredem: problem defining nslope dimension' write(*,*)trim(nf90_strerror(ierr)) CALL ABORT ENDIF ierr=NF90_DEF_DIM(nid_restart,"inter slope",nslope+1,idim9) IF (ierr/=NF90_NOERR) THEN write(*,*)'phyredem: problem defining inter slope dimension' write(*,*)trim(nf90_strerror(ierr)) CALL ABORT ENDIF ierr=NF90_DEF_DIM(nid_restart,"nb_str_max",nb_str_max,idim10) IF (ierr/=NF90_NOERR) THEN write(*,*)'phyredem: problem defining nb_str_max dimension' write(*,*)trim(nf90_strerror(ierr)) CALL ABORT ENDIF ierr=NF90_ENDDEF(nid_restart) IF (ierr/=NF90_NOERR) THEN write(*,*)'open_restartphy: problem ending definition mode ' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("open_restartphy","Failed ending definition mode",1) ENDIF ENDIF END SUBROUTINE open_restartphy !======================================================================= !======================================================================= SUBROUTINE close_restartphy !----------------------------------------------------------------------- ! NAME ! close_restartphy ! ! DESCRIPTION ! Close the restart netCDF file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_CLOSE USE mod_phys_lmdz_para, only: is_master ! DECLARATION ! ----------- implicit none ! LOCAL VARIABLES ! --------------- INTEGER :: ierr ! CODE ! ---- IF (is_master) THEN ierr = NF90_CLOSE (nid_restart) ENDIF END SUBROUTINE close_restartphy !======================================================================= !======================================================================= SUBROUTINE put_field_r1(field_name,title,field,time) !----------------------------------------------------------------------- ! NAME ! put_field_r1 ! ! DESCRIPTION ! Write a surface field to the restart file. ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: field(:) REAL, INTENT(IN), OPTIONAL :: time ! CODE ! ---- IF (present(time)) THEN ! if timeindex is present, it is a time-dependent variable CALL put_field_rgen(field_name,title,field,1,time) ELSE CALL put_field_rgen(field_name,title,field,1) ENDIF END SUBROUTINE put_field_r1 !======================================================================= !======================================================================= SUBROUTINE put_field_r2(field_name,title,field,time) !----------------------------------------------------------------------- ! NAME ! put_field_r2 ! ! DESCRIPTION ! Write a "3D" horizontal-vertical field (2D) to the restart file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: field(:,:) REAL, INTENT(IN), OPTIONAL :: time ! CODE ! ---- IF (present(time)) THEN ! if timeindex is present, it is a time-dependent variable CALL put_field_rgen(field_name,title,field,size(field,2),time) ELSE CALL put_field_rgen(field_name,title,field,size(field,2)) ENDIF END SUBROUTINE put_field_r2 !======================================================================= !======================================================================= SUBROUTINE put_field_r3(field_name,title,field,time) !----------------------------------------------------------------------- ! NAME ! put_field_r3 ! ! DESCRIPTION ! Write a "4D" field surf/alt/?? (3D) to the restart file. ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: field_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: field(:,:,:) REAL, INTENT(IN), OPTIONAL :: time ! CODE ! ---- IF (present(time)) THEN ! if timeindex is present, it is a time-dependent variable CALL put_field_rgen(field_name,title,field,size(field,2)*size(field,3),& time) ELSE CALL put_field_rgen(field_name,title,field,size(field,2)*size(field,3)) ENDIF END SUBROUTINE put_field_r3 !======================================================================= !======================================================================= SUBROUTINE put_field_rgen(field_name,title,field,field_size,time) !----------------------------------------------------------------------- ! NAME ! put_field_rgen ! ! DESCRIPTION ! Generic subroutine to write fields to restart file with gather from processes ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf USE dimphy USE soil, only: nsoilmx_PEM USE comslope_mod, ONLY: nslope USE mod_grid_phy_lmdz USE mod_phys_lmdz_para ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*),INTENT(IN) :: field_name CHARACTER(LEN=*),INTENT(IN) :: title INTEGER,INTENT(IN) :: field_size REAL,INTENT(IN) :: field(klon,field_size) REAL,OPTIONAL,INTENT(IN) :: time ! LOCAL VARIABLES ! --------------- REAL :: field_glo(klon_glo,field_size) INTEGER :: ierr INTEGER :: nvarid ! CODE ! ---- CALL gather(field,field_glo) IF (is_master) THEN IF (index(field_name,"stratif") == 0) then IF (field_size==1) THEN ! input is a 1D "surface field" array if (.not.present(time)) then ! for a time-independent field ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,timeindex/)) endif ! of if (.not.present(timeindex)) ELSE IF (field_size==klev) THEN ! input is a 2D "atmospheric field" array if (.not.present(time)) then ! for a time-independent field ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim6/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim6/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim6,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim6,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,1,timeindex/)) endif ! of if (.not.present(time)) ELSE IF (field_size==klevp1) THEN ! input is a 2D "interlayer atmospheric field" array if (.not.present(time)) then ! for a time-independent field ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim4/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim4/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr = NF90_ENDDEF(nid_restart) ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim4,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim4,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,1,timeindex/)) endif ! of if (.not.present(timeindex)) ELSE IF (field_size==nsoilmx_PEM) THEN ! input is a 2D "subsurface field" array if (.not.present(time)) then ! for a time-independent field ierr = NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim3/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim3/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr = NF90_ENDDEF(nid_restart) ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim3,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim3,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define "//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,1,timeindex/)) endif ! of if (.not.present(time)) ELSE IF (field_size==nslope) THEN ! input is a 2D "subsurface field" array if (.not.present(time)) then ! for a time-independent field ierr = NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim8/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim8/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define"//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr = NF90_ENDDEF(nid_restart) ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim8,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim8,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define"//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,1,timeindex/)) endif ! of if (.not.present(time)) ELSE write(*,*) "Error phyredem(put_field_rgen) : wrong dimension for ",trim(field_name) write(*,*) " field_size =",field_size CALL abort_physic("put_field_rgen","wrong field dimension",1) ENDIF ELSE ! input is a 2D "stratification" array if (.not.present(time)) then ! for a time-independent field ierr = NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim8/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim8/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define"//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr = NF90_ENDDEF(nid_restart) ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo) else ! check if the variable has already been defined: ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid) if (ierr/=NF90_NOERR) then ! variable not found, define it ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim10,idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim10,idim7/),nvarid) #endif if (ierr.ne.NF90_NOERR) then write(*,*)"put_field_rgen error: failed to define"//trim(field_name) write(*,*)trim(nf90_strerror(ierr)) endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) endif ! Write the variable ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,& start=(/1,1,timeindex/)) endif ! of if (.not.present(time)) ENDIF ! Check the writting of field to file went OK if (ierr.ne.NF90_NOERR) then write(*,*) " Error phyredem(put_field_rgen) : failed writing ",trim(field_name) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("put_field_rgen","Failed writing field",1) endif ENDIF ! of IF (is_master) END SUBROUTINE put_field_rgen !======================================================================= !======================================================================= SUBROUTINE put_var_r0(var_name,title,var) !----------------------------------------------------------------------- ! NAME ! put_var_r0 ! ! DESCRIPTION ! Write a scalar variable to the restart file ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: var ! LOCAL VARIABLES ! --------------- REAL :: varin(1) ! CODE ! ---- varin(1)=var CALL put_var_rgen(var_name,title,varin,size(varin)) END SUBROUTINE put_var_r0 !======================================================================= !======================================================================= SUBROUTINE put_var_r1(var_name,title,var) !----------------------------------------------------------------------- ! NAME ! put_var_r1 ! ! DESCRIPTION ! Write a vector (1D) variable to the restart file ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: var(:) ! CODE ! ---- CALL put_var_rgen(var_name,title,var,size(var)) END SUBROUTINE put_var_r1 !======================================================================= !======================================================================= SUBROUTINE put_var_r2(var_name,title,var) !----------------------------------------------------------------------- ! NAME ! put_var_r2 ! ! DESCRIPTION ! Write a 3D field variable to the restart file ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: var(:,:) ! CODE ! ---- CALL put_var_rgen(var_name,title,var,size(var)) END SUBROUTINE put_var_r2 !======================================================================= !======================================================================= SUBROUTINE put_var_r3(var_name,title,var) !----------------------------------------------------------------------- ! NAME ! put_var_r3 ! ! DESCRIPTION ! Write a 3D field variable to the restart file ! ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name CHARACTER(LEN=*), INTENT(IN) :: title REAL, INTENT(IN) :: var(:,:,:) ! CODE ! ---- CALL put_var_rgen(var_name,title,var,size(var)) END SUBROUTINE put_var_r3 !======================================================================= !======================================================================= SUBROUTINE put_var_rgen(var_name,title,var,var_size) !----------------------------------------------------------------------- ! NAME ! put_var_rgen ! ! DESCRIPTION ! Generic subroutine to write variables to the restart file ! AUTHORS & DATE ! L. Lange ! JB Clement, 2023-2025 ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ USE netcdf, only: NF90_REDEF, NF90_DEF_VAR, NF90_ENDDEF, NF90_PUT_VAR, & NF90_FLOAT, NF90_DOUBLE, & NF90_PUT_ATT, NF90_NOERR, nf90_strerror, & nf90_inq_dimid, nf90_inquire_dimension, NF90_INQ_VARID USE soil, only: nsoilmx_PEM USE comslope_mod, only: nslope USE mod_phys_lmdz_para, only: is_master ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- CHARACTER(LEN=*), INTENT(IN) :: var_name CHARACTER(LEN=*), INTENT(IN) :: title INTEGER, INTENT(IN) :: var_size REAL, INTENT(IN) :: var(var_size) ! LOCAL VARIABLES ! --------------- INTEGER :: ierr INTEGER :: nvarid INTEGER :: idim1d LOGICAL, SAVE :: firsttime = .true. ! CODE ! ---- IF (is_master) THEN IF (var_name=="Time") THEN ! Very specific case of "Time" variable if (firsttime) then ! Create the "Time variable" ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_DOUBLE,& (/idim7/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_FLOAT,& (/idim7/),nvarid) #endif IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ierr=NF90_ENDDEF(nid_restart) firsttime=.false. endif ! Append "time" value ! get current length of "Time" dimension ierr=nf90_inq_dimid(nid_restart,var_name,idim1d) ierr=nf90_inquire_dimension(nid_restart,idim1d,len=timeindex) timeindex=timeindex+1 ierr=NF90_INQ_VARID(nid_restart,var_name,nvarid) ierr=NF90_PUT_VAR(nid_restart,nvarid,var,& start=(/timeindex/)) IF (ierr/=NF90_NOERR) THEN write(*,*)'put_var_rgen: problem writing Time' write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("get_var_rgen","Failed to write Time",1) ENDIF return ! nothing left to do ELSEIF (var_size==length) THEN ! We know it is a "controle" kind of 1D array idim1d=idim1 ELSEIF (var_size==nsoilmx_PEM) THEN ! We know it is an "mlayer" kind of 1D array idim1d=idim3 ELSEIF (var_size==nslope+1) THEN ! We know it is an "inter slope" kind of 1D array idim1d=idim9 ELSEIF (var_name == "nb_str_max") THEN ! We know it is a kind of stratification array idim1d = idim10 ELSE write(*,*) "put_var_rgen error : wrong dimension" write(*,*) " var_size =",var_size CALL abort_physic("get_var_rgen","Wrong variable dimension",1) ENDIF ! of IF (var_size==length) THEN ! Swich to NetCDF define mode ierr=NF90_REDEF (nid_restart) ! Define the variable #ifdef NC_DOUBLE ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_DOUBLE,(/idim1d/),nvarid) #else ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_FLOAT,(/idim1d/),nvarid) #endif ! Add a "title" attribute IF (LEN_TRIM(title)>0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title) ! Swich out of define mode ierr=NF90_ENDDEF(nid_restart) ! Write variable to file ierr=NF90_PUT_VAR(nid_restart,nvarid,var) IF (ierr/=NF90_NOERR) THEN write(*,*)'put_var_rgen: problem writing '//trim(var_name) write(*,*)trim(nf90_strerror(ierr)) CALL abort_physic("get_var_rgen","Failed writing variable",1) ENDIF ENDIF ! of IF (is_master) END SUBROUTINE put_var_rgen !======================================================================= END MODULE iostart_pem