diff --git a/GFDL_tools/fv_ada_nudge.F90 b/GFDL_tools/fv_ada_nudge.F90 index b091e556c..90cfe8510 100644 --- a/GFDL_tools/fv_ada_nudge.F90 +++ b/GFDL_tools/fv_ada_nudge.F90 @@ -38,9 +38,7 @@ module fv_ada_nudge_mod use external_sst_mod, only: i_sst, j_sst, sst_ncep, sst_anom, forecast_mode use diag_manager_mod, only: register_diag_field, send_data use constants_mod, only: pi, grav, rdgas, cp_air, kappa, cnst_radius=>radius, seconds_per_day - use fms_mod, only: write_version_number, open_namelist_file, & - check_nml_error, file_exist, close_file -!use fms_io_mod, only: field_size + use fms_mod, only: write_version_number, check_nml_error use mpp_mod, only: mpp_error, FATAL, stdlog, get_unit, mpp_pe, input_nml_file use mpp_mod, only: mpp_root_pe, stdout ! snz use mpp_mod, only: mpp_clock_id, mpp_clock_begin, mpp_clock_end @@ -63,8 +61,10 @@ module fv_ada_nudge_mod get_var3_r4, get_var2_r4, get_var1_real use fv_arrays_mod, only: fv_grid_type, fv_grid_bounds_type, fv_nest_type, R_GRID - use fms_io_mod, only: register_restart_field, restart_file_type, restore_state - use fms_io_mod, only: save_restart, get_mosaic_tile_file + use fms2_io_mod, only : register_restart_field, open_file, close_file, & + read_restart, register_field, & + register_variable_attribute, file_exists + use fv_io_mod, only : fv_io_register_axis use axis_utils_mod, only : frac_index #ifdef ENABLE_ADA @@ -235,8 +235,9 @@ module fv_ada_nudge_mod integer :: id_u_da, id_v_da, id_t_da, id_q_da, id_ps_da ! snz integer :: id_ada - type(restart_file_type) :: ada_driver_restart ! snz - character(len=*), parameter :: restart_file="ada_driver.res.nc" ! snz + type(FmsNetcdfFile_t) :: ada_driver_restart ! snz + character(len=*), parameter :: restart_file="INPUT/ada_driver.res.nc" ! snz + character(len=8), dimension(4) :: dim_names_4d #ifdef ENABLE_ADA ! snz type(model_data_type) :: Atm_var @@ -1482,7 +1483,6 @@ subroutine fv_ada_nudge_init(time, axes, npz, zvir, ak, bk, ts, phis, gridstruct integer :: id_restart !< Currently not used for anything other than a return !value. - character(len=256) :: restart_file_instance real, pointer, dimension(:,:,:) :: agrid @@ -1522,20 +1522,9 @@ subroutine fv_ada_nudge_init(time, axes, npz, zvir, ak, bk, ts, phis, gridstruct track_file_name = "No_File_specified" -#ifdef INTERNAL_FILE_NML - read(input_nml_file, nml = fv_ada_nudge_nml, iostat = io) - ierr = check_nml_error(io,'fv_ada_nudge_nml') -#else - if( file_exist( 'input.nml' ) ) then - unit = open_namelist_file () - io = 1 - do while ( io .ne. 0 ) - read( unit, nml = fv_ada_nudge_nml, iostat = io, end = 10 ) - ierr = check_nml_error(io,'fv_ada_nudge_nml') - end do -10 call close_file ( unit ) - end if -#endif + read(input_nml_file, nml = fv_ada_nudge_nml, iostat = io) + ierr = check_nml_error(io,'fv_ada_nudge_nml') + call write_version_number ( 'FV_ADA_NUDGE_MOD', version ) if ( master ) then f_unit=stdlog() @@ -1793,26 +1782,32 @@ subroutine fv_ada_nudge_init(time, axes, npz, zvir, ak, bk, ts, phis, gridstruct 'da ps', 'Pa', missing_value=missing_value ) ! snz add the following lines for recording the return values from the previous assim run - call get_mosaic_tile_file(restart_file, restart_file_instance,& - & .FALSE., domain) - - id_restart = register_restart_field(ada_driver_restart, restart_file_instance, & - & "u_adj", Atm_var%u_adj(:,:,:), domain=domain) - id_restart = register_restart_field(ada_driver_restart, restart_file_instance, & - & "v_adj", Atm_var%v_adj(:,:,:), domain=domain) - id_restart = register_restart_field(ada_driver_restart, restart_file_instance, & - & "t_adj", Atm_var%t_adj(:,:,:), domain=domain) - id_restart = register_restart_field(ada_driver_restart, restart_file_instance, & - & "q_adj", Atm_var%q_adj(:,:,:), domain=domain) - id_restart = register_restart_field(ada_driver_restart, restart_file_instance, & - & "ps_adj", Atm_var%ps_adj(:,:), domain=domain) - - if ( file_exist('INPUT/'//trim(restart_file_instance), domain=domain) ) then - if ( mpp_pe() .eq. mpp_root_pe() ) then - write (stdout_unit,*) 'Reading ada restart information from ', 'INPUT/'//trim(restart_file_instance) - end if - call restore_state(ada_driver_restart, DIRECTORY='INPUT') - end if + +! set dimensions for register restart + dim_names_4d(1) = "xaxis_1" + dim_names_4d(2) = "yaxis_1" + dim_names_4d(3) = "zaxis_1" + dim_names_4d(4) = "Time" + + if (open_file(ada_driver_restart, restart_file, "read", domain, is_restart=.true.)) + call fv_io_register_axis(ada_driver_restart, numx=1, numy=1, numz=1, zsize=(/size(Atm_var%u_adj,3)/)) + call register_restart_field(ada_driver_restart, & + & "u_adj", Atm_var%u_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "v_adj", Atm_var%v_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "t_adj", Atm_var%t_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "q_adj", Atm_var%q_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "ps_adj", Atm_var%ps_adj(:,:), dim_names_4d) + + if ( mpp_pe() .eq. mpp_root_pe() ) then + write (stdout_unit,*) 'Reading ada restart information from ', 'INPUT/'//trim(restart_file) + end if + call read_restart(ada_driver_restart) + call close_file(ada_driver_restart) + endif #endif ! snz for ENABLE_ADA @@ -1848,7 +1843,7 @@ subroutine get_ncep_analysis ( ps, u, v, t, q, zvir, ts, nfile, fname, bd, time !! model run integer :: is, ie, js, je - if( .not. file_exist(fname) ) then + if( .not. file_exists(fname) ) then call mpp_error(FATAL,'==> Error from get_ncep_analysis: file not found') else call open_ncfile( fname, ncid ) ! open the file @@ -2516,7 +2511,21 @@ subroutine fv_ada_nudge_end #ifdef ENABLE_ADA ! snz - call save_restart(ada_driver_restart) + if (open_file(ada_driver_restart, restart_file, "overwrite", domain, is_restart=.true.)) then + call fv_io_register_axis(ada_driver_restart, numx=1, numy=1, numz=1, zsize=(/size(Atm_var%u_adj,3)/)) + call register_restart_field(ada_driver_restart, & + & "u_adj", Atm_var%u_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "v_adj", Atm_var%v_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "t_adj", Atm_var%t_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "q_adj", Atm_var%q_adj(:,:,:), dim_names_4d) + call register_restart_field(ada_driver_restart, & + & "ps_adj", Atm_var%ps_adj(:,:), dim_names_4d) + call write_restart(ada_driver_restart) + call close_file(ada_driver_restart) + endif deallocate ( Atm_var%u, Atm_var%v, Atm_var%t, Atm_var%u_adj, Atm_var%v_adj, Atm_var%t_adj) deallocate ( Atm_var%q, Atm_var%ps, Atm_var%q_adj, Atm_var%ps_adj) diff --git a/GFDL_tools/fv_climate_nudge.F90 b/GFDL_tools/fv_climate_nudge.F90 index bc490228d..b392bef2e 100644 --- a/GFDL_tools/fv_climate_nudge.F90 +++ b/GFDL_tools/fv_climate_nudge.F90 @@ -20,10 +20,10 @@ !*********************************************************************** module fv_climate_nudge_mod -use fms_mod, only: open_namelist_file, check_nml_error, & - close_file, stdlog, mpp_pe, mpp_root_pe, & +use fms_mod, only: check_nml_error, & + stdlog, mpp_pe, mpp_root_pe, & write_version_number, string, error_mesg, & - FATAL, WARNING, NOTE, file_exist + FATAL, WARNING, NOTE use mpp_mod, only: input_nml_file use diag_manager_mod, only: register_diag_field, send_data, & register_static_field @@ -129,20 +129,8 @@ subroutine fv_climate_nudge_init ( Time, axes, flag ) if (module_is_initialized) return ! read namelist -#ifdef INTERNAL_FILE_NML read (input_nml_file, nml=fv_climate_nudge_nml, iostat=io) ierr = check_nml_error (io, 'fv_climate_nudge_nml') -#else - if (file_exist('input.nml') ) then - unit = open_namelist_file() - ierr=1 - do while (ierr /= 0) - read (unit, nml=fv_climate_nudge_nml, iostat=io, end=10) - ierr = check_nml_error (io, 'fv_climate_nudge_nml') - enddo -10 call close_file (unit) - endif -#endif !----- write version and namelist to log file ----- diff --git a/GFDL_tools/fv_cmip_diag.F90 b/GFDL_tools/fv_cmip_diag.F90 index 91c5d40f9..0ef4e558a 100644 --- a/GFDL_tools/fv_cmip_diag.F90 +++ b/GFDL_tools/fv_cmip_diag.F90 @@ -21,11 +21,10 @@ module fv_cmip_diag_mod use mpp_mod, only: input_nml_file -use fms_mod, only: open_namelist_file, check_nml_error, & - close_file, stdlog, mpp_pe, mpp_root_pe, & - write_version_number, file_exist, & +use fms_mod, only: check_nml_error, string, & + stdlog, mpp_pe, mpp_root_pe, & + write_version_number, & error_mesg, FATAL, WARNING, NOTE -use fms_io_mod, only: set_domain, nullify_domain, string use time_manager_mod, only: time_type use mpp_domains_mod, only: domain2d use diag_manager_mod, only: register_diag_field, diag_axis_init, & @@ -120,20 +119,8 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time ) endif !----- read namelist ----- -#ifdef INTERNAL_FILE_NML read (input_nml_file, nml=fv_cmip_diag_nml, iostat=io) ierr = check_nml_error (io, 'fv_cmip_diag_nml') -#else - if (file_exist('input.nml') ) then - iunit = open_namelist_file() - ierr=1 - do while (ierr /= 0) - read (iunit, nml=fv_cmip_diag_nml, iostat=io, end=10) - ierr = check_nml_error (io, 'fv_cmip_diag_nml') - enddo -10 call close_file (iunit) - endif -#endif !----- write version and namelist to log file ----- @@ -141,10 +128,6 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time ) call write_version_number ( 'FV_CMIP_DIAG_MOD', version ) if (mpp_pe() == mpp_root_pe()) write (iunit, nml=fv_cmip_diag_nml) - -! Set domain so that diag_manager can access tile information - call set_domain(Atm(1)%domain) - ! axis identifiers area_id = get_diag_field_id ('dynamics', 'area') if (area_id .eq. DIAG_FIELD_NOT_FOUND) call error_mesg & @@ -411,7 +394,6 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time ) call diag_field_add_attribute (id_zg1000, 'coordinates', 'p1000') !--- done --- - call nullify_domain() module_is_initialized=.true. !----------------------------------------------------------------------- @@ -455,8 +437,6 @@ subroutine fv_cmip_diag ( Atm, zvir, Time ) ngc = Atm(n)%ng npz = Atm(n)%npz - call set_domain(Atm(n)%domain) - ! set flags for computing quantities compute_rh = .false. compute_wa = .false. @@ -688,10 +668,6 @@ subroutine fv_cmip_diag ( Atm, zvir, Time ) used = send_data (id_zg1000, dat3(:,:,1), Time) endif -!---------------------------------------------------------------------- - - call nullify_domain() - !---------------------------------------------------------------------- end subroutine fv_cmip_diag diff --git a/GFDL_tools/fv_diag_column.F90 b/GFDL_tools/fv_diag_column.F90 index 7a76e27b3..25a3cf008 100644 --- a/GFDL_tools/fv_diag_column.F90 +++ b/GFDL_tools/fv_diag_column.F90 @@ -8,7 +8,6 @@ module fv_diag_column_mod use fms_mod, only: write_version_number, lowercase use mpp_mod, only: mpp_error, FATAL, stdlog, mpp_pe, mpp_root_pe, mpp_sum, & mpp_max, NOTE, input_nml_file, get_unit - use mpp_io_mod, only: mpp_flush use fv_sg_mod, only: qsmith implicit none @@ -79,20 +78,7 @@ subroutine fv_diag_column_init(Atm, yr_init, mo_init, dy_init, hr_init, do_diag_ diag_sonde_j(:) = -999 diag_sonde_tile(:) = -99 -#ifdef INTERNAL_FILE_NML read(input_nml_file, nml=fv_diag_column_nml,iostat=ios) -#else - inquire (file=trim(Atm%nml_filename), exist=exists) - if (.not. exists) then - write(errmsg,*) 'fv_diag_column_nml: namelist file ',trim(Atm%nml_filename),' does not exist' - call mpp_error(FATAL, errmsg) - else - open (unit=nlunit, file=Atm%nml_filename, READONLY, status='OLD', iostat=ios) - endif - rewind(nlunit) - read (nlunit, nml=fv_diag_column_nml, iostat=ios) - close (nlunit) -#endif if (do_diag_debug .or. do_diag_sonde) then call read_column_table @@ -394,10 +380,8 @@ subroutine debug_column(pt, delp, delz, u, v, w, q, npz, ncnst, sphum, nwat, zvi write(unit, *) '===================================================================' write(unit, *) - - call mpp_flush(unit) - - + flush(unit) + enddo end subroutine debug_column @@ -487,8 +471,6 @@ subroutine debug_column_dyn(pt, delp, delz, u, v, w, q, heat_source, cappa, akap write(unit, *) '===================================================================' write(unit, *) - call mpp_flush(unit) - enddo end subroutine debug_column_dyn @@ -582,8 +564,6 @@ subroutine sounding_column( pt, delp, delz, u, v, q, peln, pkz, thetae, phis, & enddo endif - call mpp_flush(unit) - enddo diff --git a/GFDL_tools/read_climate_nudge_data.F90 b/GFDL_tools/read_climate_nudge_data.F90 index 6122478cd..96c5baccc 100644 --- a/GFDL_tools/read_climate_nudge_data.F90 +++ b/GFDL_tools/read_climate_nudge_data.F90 @@ -1,3 +1,4 @@ + !*********************************************************************** !* GNU Lesser General Public License !* @@ -20,14 +21,15 @@ !*********************************************************************** module read_climate_nudge_data_mod -use fms_mod, only: open_namelist_file, check_nml_error, close_file, & +use fms_mod, only: check_nml_error, & stdlog, mpp_pe, mpp_root_pe, write_version_number, & - string, error_mesg, FATAL, NOTE, file_exist -use mpp_mod, only: input_nml_file -use mpp_io_mod, only: mpp_open, MPP_NETCDF, MPP_RDONLY,MPP_MULTI, MPP_SINGLE -use mpp_io_mod, only: axistype, fieldtype, mpp_get_time_axis, mpp_get_atts -use mpp_io_mod, only: mpp_get_fields, mpp_get_info, mpp_get_axes, mpp_get_times -use mpp_io_mod, only: mpp_get_axis_data, mpp_read, mpp_close, mpp_get_default_calendar + string, error_mesg, FATAL, NOTE +use fms2_io_mod, only: open_file, close_file, get_num_dimensions, & + get_dimension_names, get_dimension_size, FmsNetcdfFile_t, & + get_num_variables, get_variable_names, get_variable_size, & + get_variable_num_dimensions, get_variable_units, & + get_time_calendar, read_data, variable_att_axists +use mpp_mod, only: input_nml_file, mpp_npes, mpp_get_current_pelist use constants_mod, only: PI, GRAV, RDGAS, RVGAS implicit none @@ -50,7 +52,7 @@ module read_climate_nudge_data_mod real, parameter :: D608 = RVGAS/RDGAS - 1. integer, parameter :: NUM_REQ_AXES = 3 - integer, parameter :: INDEX_LON = 1, INDEX_LAT = 2, INDEX_LEV = 3 + integer, parameter :: INDEX_LON = 1, INDEX_LAT = 2, INDEX_LEV = 3, INDEX_TIME = 4 character(len=8), dimension(NUM_REQ_AXES) :: required_axis_names = & (/ 'lon', 'lat', 'lev' /) @@ -78,14 +80,13 @@ module read_climate_nudge_data_mod integer, allocatable :: file_index(:) type filedata_type - integer :: ncid - integer, pointer :: length_axes(:) ! length of all dimensions in file - integer :: ndim, nvar, natt, ntim, varid_time + integer, allocatable :: dim_size + integer :: ndim, nvar, ntim integer :: time_offset integer, dimension(NUM_REQ_FLDS) :: field_index ! varid for variables integer, dimension(NUM_REQ_AXES) :: axis_index ! varid for dimensions - type(axistype), dimension(NUM_REQ_FLDS) :: axes - type(fieldtype), dimension(NUM_REQ_FLDS) :: fields + character(len=32), dimension(NUM_REQ_FLDS+1) :: axes + character(len=32), dimension(NUM_REQ_FLDS) :: fields end type type(filedata_type), allocatable :: Files(:) @@ -102,10 +103,10 @@ subroutine read_climate_nudge_data_init (nlon, nlat, nlev, ntime) ! ntime number of time levels integer :: iunit, ierr, io - character(len=128) :: name integer :: istat, i, j, k, n, nd, siz(4), i1, i2 - type(axistype), allocatable :: axes(:) - type(fieldtype), allocatable :: fields(:) + character(len=32), allocatable :: fields(:) + type(FmsNetcdfFile_t) :: fileobj + integer, allocatable, dimension(:) :: pes !< Array of ther pes in the current pelist if (module_is_initialized) return ! initial file names to blanks @@ -117,20 +118,8 @@ subroutine read_climate_nudge_data_init (nlon, nlat, nlev, ntime) enddo !----- read namelist ----- -#ifdef INTERNAL_FILE_NML read (input_nml_file, nml=read_climate_nudge_data_nml, iostat=io) ierr = check_nml_error (io, 'read_climate_nudge_data_nml') -#else - if (file_exist('input.nml') ) then - iunit = open_namelist_file() - ierr=1 - do while (ierr /= 0) - read (iunit, nml=read_climate_nudge_data_nml, iostat=io, end=10) - ierr = check_nml_error (io, 'read_climate_nudge_data_nml') - enddo -10 call close_file (iunit) - endif -#endif !----- write version and namelist to log file ----- @@ -151,60 +140,63 @@ subroutine read_climate_nudge_data_init (nlon, nlat, nlev, ntime) numtime = 0 ! open input file(s) - do n = 1, nfiles - call mpp_open(Files(n)%ncid, trim(filenames(n)), form=MPP_NETCDF,action=MPP_RDONLY,threading=MPP_MULTI, & - fileset=MPP_SINGLE ) -! call error_mesg ('read_climate_nudge_data_mod', 'Buffer size for reading = '//trim(string(read_buffer_size)), NOTE) - - call mpp_get_info(Files(n)%ncid, Files(n)%ndim, Files(n)%nvar, Files(n)%natt, Files(n)%ntim) - - allocate (Files(n)%length_axes(Files(n)%ndim)) - allocate (axes(Files(n)%ndim)) - call mpp_get_axes(Files(n)%ncid,axes) - - ! inquire dimension sizes - do i = 1, Files(n)%ndim - call mpp_get_atts(axes(i), name=name, len=Files(n)%length_axes(i)) - do j = 1, NUM_REQ_AXES - if (trim(name) .eq. trim(required_axis_names(j))) then - call check_axis_size (j,Files(n)%length_axes(i)) - Files(n)%axes(j) = axes(i) - Files(n)%axis_index(j) = i - exit - endif - enddo - enddo - deallocate(axes) - ! time axis indexing - Files(n)%time_offset = numtime - numtime = numtime + Files(n)%ntim - - allocate(fields(Files(n)%nvar)) - call mpp_get_fields(Files(n)%ncid,fields) - Files(n)%field_index = 0 - do i = 1, Files(n)%nvar - call mpp_get_atts(fields(i), name=name, ndim=nd, siz=siz) - do j = 1, NUM_REQ_FLDS - if (trim(name) .eq. trim(required_field_names(j))) then - Files(n)%field_index(j) = i - Files(n)%fields(j) = fields(i) - if (j .gt. 3) then - call check_resolution (siz(1:nd)) - endif - exit - endif - enddo - - ! special case for surface geopotential (sometimes the name is PHIS) - if (trim(name) .eq. 'PHIS') then - Files(n)%field_index(INDEX_ZS) = i - Files(n)%fields(INDEX_ZS) = fields(i) - call check_resolution (siz(1:nd)) - endif - enddo - deallocate(fields) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + do n = 1, nfiles + if (open_file(fileobj, trim(filenames(n)), "read", pelist=pes)) then + Files(n)%ndim=get_num_dimensions(fileobj) + + allocate (Files(n)%dim_size(Files(n)%ndim)) + + ! inquire dimension sizes + do i = 1, Files(n)%ndim + call get_dimension_names(fileobj, Files(n)%axes) + do j = 1, NUM_REQ_AXES + if (trim(Files(n)%axes(j)) .eq. trim(required_axis_names(j))) then + call get_dimension_size(fileobj, Files(n)%axes(j), Files(n)%dim_size(j)) + call check_axis_size (j,Files(n)%dim_size(j)) + Files(n)%axis_index(j) = i + exit + endif + enddo + enddo + ! time axis indexing + call get_dimension_size(fileobj, Files(n)%axes(INDEX_TIME), Files(n)%dim_size(INDEX_TIME)) + Files(n)%ntim = Files(n)%dim_size(INDEX_TIME) + Files(n)%time_offset = numtime + numtime = numtime + Files(n)%ntim + + Files(n)%nvar = get_num_variables(fileobj) + allocate(fields(Files(n)%nvar)) + call get_variable_names(fileobj, fields) + Files(n)%field_index = 0 + do i = 1, Files(n)%nvar + nd = get_variable_num_dimensions(fileobj, fields(i)) + call get_variable_size(fileobj, fields(i), siz) + do j = 1, NUM_REQ_FLDS + if (trim(fields(i)) .eq. trim(required_field_names(j))) then + Files(n)%field_index(j) = i + Files(n)%fields(j) = fields(i) + if (j .gt. 3) then + call check_resolution (siz(1:nd)) + endif + exit + endif + enddo + + ! special case for surface geopotential (sometimes the name is PHIS) + if (trim(fields(i)) .eq. 'PHIS') then + Files(n)%field_index(INDEX_ZS) = i + Files(n)%fields(INDEX_ZS) = fields(i) + call check_resolution (siz(1:nd)) + endif + enddo + deallocate(fields) + call close_file(fileobj) + endif enddo ! "n" files loop + deallocate(pes) ! setup file indexing allocate(file_index(numtime)) @@ -234,8 +226,9 @@ subroutine read_time ( times, units, calendar ) character(len=*), intent(out) :: units, calendar integer :: istat, i1, i2, n real :: l_times(size(times)) -type(axistype), save:: time_axis -character(len=32) :: default_calendar +character(len=32), save:: time_axis +type(FmsNetcdfFile_t) :: fileobj +integer, allocatable, dimension(:) :: pes !< Array of ther pes in the current pelist if (.not.module_is_initialized) then call error_mesg ('read_climate_nudge_data_mod/read_time', & @@ -248,18 +241,27 @@ subroutine read_time ( times, units, calendar ) ! data i2 = 0 + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) do n = 1, nfiles - i1 = i2+1 - i2 = i2+Files(n)%ntim - if( n == 1) then - default_calendar = mpp_get_default_calendar() - call mpp_get_time_axis( Files(n)%ncid, time_axis) - call mpp_get_atts(time_axis, units=units, calendar=calendar) - if( trim(calendar) == trim(default_calendar)) calendar = 'gregorian ' + if (open_file(fileobj, trim(filenames(n)), "read", pelist=pes)) then + i1 = i2+1 + i2 = i2+Files(n)%ntim + if( n == 1) then + time_axis = Files(n)%axes(INDEX_TIME) + call get_variable_units(fileobj, time_axis, units) + if( variable_att_exists(fileobj, time_axis, calendar) then + call get_time_calendar(fileobj, time_axis, calendar) + else + calendar = 'gregorian ' + endif + endif + call read_data(fileobj,time_axis, l_times(i1:i2)) + times(i1:i2) = l_times(i1:i2) + call close_file(fileobj) endif - call mpp_get_times(Files(n)%ncid, l_times(i1:i2)) - times(i1:i2) = l_times(i1:i2) enddo + deallocate(pes) ! NOTE: need to do the conversion to time_type in this routine ! this will allow different units and calendars for each file @@ -273,6 +275,8 @@ subroutine read_grid ( lon, lat, ak, bk ) real :: pref integer :: istat + type(FmsNetcdfFile_t) :: fileobj + integer, allocatable, dimension(:) :: pes !< Array of ther pes in the current pelist if (.not.module_is_initialized) then call error_mesg ('read_climate_nudge_data_mod/read_grid', & @@ -281,29 +285,34 @@ subroutine read_grid ( lon, lat, ak, bk ) ! static fields from first file only - call mpp_get_axis_data(Files(1)%axes(INDEX_LON), lon) - call mpp_get_axis_data(Files(1)%axes(INDEX_LAT), lat) - - ! units are assumed to be degrees east and north - ! convert to radians - lon = lon * PI/180. - lat = lat * PI/180. - - ! vertical coodinate - if (Files(1)%field_index(INDEX_AK) .gt. 0) then - call mpp_read(Files(1)%ncid, Files(1)%fields(INDEX_AK), ak) - if (Files(1)%field_index(INDEX_P0) .gt. 0) then - call mpp_read(Files(1)%ncid, Files(1)%fields(INDEX_P0), pref) - else - pref = P0 - endif - ak = ak*pref - else - ak = 0. + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if (open_file(fileobj, trim(filenames(1)), "read", pelist=pes)) then + call read_data(fileobj, Files(1)%axes(INDEX_LON), lon) + call read_data(fileobj, Files(1)%axes(INDEX_LAT), lat) + + ! units are assumed to be degrees east and north + ! convert to radians + lon = lon * PI/180. + lat = lat * PI/180. + + ! vertical coodinate + if (Files(1)%field_index(INDEX_AK) .gt. 0) then + call read_data(fileobj, Files(1)%fields(INDEX_AK), ak) + if (Files(1)%field_index(INDEX_P0) .gt. 0) then + call read_data(fileobj, Files(1)%fields(INDEX_P0), pref) + else + pref = P0 + endif + ak = ak*pref + else + ak = 0. + endif + + call read_data(fileobj, Files(1)%fields(INDEX_BK), bk) + call close_file(fileobj) endif - - call mpp_read(Files(1)%ncid, Files(1)%fields(INDEX_BK), bk) - + deallocate(pes) end subroutine read_grid @@ -382,6 +391,8 @@ subroutine read_climate_nudge_data_2d (itime, field, dat, is, js) integer, intent(in), optional :: is, js integer :: istat, atime, n, this_index integer :: nread(4), start(4) +type(FmsNetcdfFile_t) :: fileobj +integer, allocatable, dimension(:) :: pes !< Array of ther pes in the current pelist if (.not.module_is_initialized) then call error_mesg ('read_climate_nudge_data_mod', & @@ -428,7 +439,13 @@ subroutine read_climate_nudge_data_2d (itime, field, dat, is, js) nread(1) = size(dat,1) nread(2) = size(dat,2) - call mpp_read(Files(n)%ncid, Files(n)%fields(this_index), dat, start, nread) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if (open_file(fileobj, trim(filenames(n)), "read", pelist=pes)) then + call read_data(fileobj, Files(n)%fields(this_index), dat, corner=start, edge_lengths=nread) + call close_file(fileobj) + endif + deallocate(pes) ! geopotential height (convert to m2/s2 if necessary) if (field .eq. 'phis') then @@ -449,6 +466,8 @@ subroutine read_climate_nudge_data_3d (itime, field, dat, is, js) real, intent(out), dimension(:,:,:) :: dat integer, intent(in), optional :: is, js integer :: istat, atime, n, this_index, start(4), nread(4) +type(FmsNetcdfFile_t) :: fileobj +integer, allocatable, dimension(:) :: pes !< Array of ther pes in the current pelist !logical :: convert_virt_temp = .false. if (.not.module_is_initialized) then @@ -505,7 +524,13 @@ subroutine read_climate_nudge_data_3d (itime, field, dat, is, js) nread(2) = size(dat,2) nread(3) = size(dat,3) - call mpp_read(Files(n)%ncid, Files(n)%fields(this_index), dat, start, nread) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if (open_file(fileobj, trim(filenames(n)), "read", pelist=pes)) then + call read_data(fileobj, Files(n)%fields(this_index), dat, corner=start, edge_lengths=nread) + call close_file(fileobj) + endif + deallocate(pes) ! convert virtual temp to temp ! necessary for some of the high resol AVN analyses @@ -521,9 +546,6 @@ subroutine read_climate_nudge_data_end integer :: istat, n if ( .not.module_is_initialized) return - do n = 1, nfiles - call mpp_close(Files(n)%ncid) - enddo deallocate (Files) module_is_initialized = .false. diff --git a/driver/GFDL/atmosphere.F90 b/driver/GFDL/atmosphere.F90 index 3e64add09..e32d35564 100644 --- a/driver/GFDL/atmosphere.F90 +++ b/driver/GFDL/atmosphere.F90 @@ -34,14 +34,14 @@ module atmosphere_mod use block_control_mod, only: block_control_type use constants_mod, only: cp_air, rdgas, grav, rvgas, kappa, pstd_mks use time_manager_mod, only: time_type, get_time, set_time, operator(+) -use fms_mod, only: file_exist, open_namelist_file, & - close_file, error_mesg, FATAL, & +use fms_mod, only: error_mesg, FATAL, & check_nml_error, stdlog, & write_version_number, & - mpp_pe, mpp_root_pe, set_domain, & + mpp_pe, mpp_root_pe, & mpp_clock_id, mpp_clock_begin, & mpp_clock_end, CLOCK_SUBCOMPONENT, & - clock_flag_default, nullify_domain + clock_flag_default +use fms2_io_mod, only: file_exists use mpp_mod, only: mpp_error, FATAL, NOTE, input_nml_file, & mpp_npes, mpp_get_current_pelist, & mpp_set_current_pelist, stdout, & @@ -210,7 +210,7 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Surf_diff, Grid_box) !----- initialize FV dynamical core ----- !NOTE do we still need the second file_exist call? - cold_start = (.not.file_exist('INPUT/fv_core.res.nc') .and. .not.file_exist('INPUT/fv_core.res.tile1.nc')) + cold_start = (.not.file_exists('INPUT/fv_core.res.nc') .and. .not.file_exists('INPUT/fv_core.res.tile1.nc')) call fv_control_init( Atm, dt_atmos, mygrid, grids_on_this_pe, p_split ) ! allocates Atm components; sets mygrid @@ -264,7 +264,6 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Surf_diff, Grid_box) ! Allocate grid variables to be used to calculate gradient in 2nd order flux exchange ! This data is only needed for the COARSEST grid. !call switch_current_Atm(Atm(mygrid)) - call set_domain(Atm(mygrid)%domain) allocate(Grid_box%dx ( isc:iec , jsc:jec+1)) allocate(Grid_box%dy ( isc:iec+1, jsc:jec )) @@ -318,12 +317,11 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Surf_diff, Grid_box) Atm(mygrid)%atmos_axes(4), Atm(mygrid)%coarse_graining) endif if (Atm(mygrid)%coarse_graining%write_coarse_restart_files) then - call fv_coarse_restart_init(mygrid, Atm(mygrid)%npz, Atm(mygrid)%flagstruct%nt_prog, & + call fv_coarse_restart_init(Atm(mygrid)%npz, Atm(mygrid)%flagstruct%nt_prog, & Atm(mygrid)%flagstruct%nt_phys, Atm(mygrid)%flagstruct%hydrostatic, & Atm(mygrid)%flagstruct%hybrid_z, Atm(mygrid)%flagstruct%fv_land, & Atm(mygrid)%coarse_graining%write_coarse_dgrid_vel_rst, & Atm(mygrid)%coarse_graining%write_coarse_agrid_vel_rst, & - Atm(mygrid)%coarse_graining%domain, & Atm(mygrid)%coarse_graining%restart) endif @@ -439,8 +437,6 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Surf_diff, Grid_box) call timing_off('ATMOS_INIT') - call set_domain(Atm(mygrid)%domain) - end subroutine atmosphere_init @@ -599,10 +595,6 @@ subroutine atmosphere_end (Time, Grid_box )!rab, Radiation, Physics) !rab type (radiation_type), intent(inout) :: Radiation !rab type (physics_type), intent(inout) :: Physics - ! initialize domains for writing global physics data - call set_domain ( Atm(mygrid)%domain ) - - !--- end nudging module --- #if defined (ATMOS_NUDGE) if ( Atm(mygrid)%flagstruct%nudge ) call atmos_nudge_end @@ -620,7 +612,6 @@ subroutine atmosphere_end (Time, Grid_box )!rab, Radiation, Physics) call atmos_global_diag_end call fv_cmip_diag_end - call nullify_domain ( ) call fv_end(Atm, mygrid) deallocate (Atm) @@ -916,8 +907,6 @@ subroutine atmosphere_state_update (Time, Physics_tendency, Physics, Atm_block) n = mygrid - call set_domain ( Atm(mygrid)%domain ) - !--- put u/v tendencies into haloed arrays u_dt and v_dt !$OMP parallel do default(shared) private(nb, ibs, ibe, jbs, jbe) do nb = 1,Atm_block%nblks @@ -999,7 +988,6 @@ subroutine atmosphere_state_update (Time, Physics_tendency, Physics, Atm_block) endif #endif - call nullify_domain() !---- diagnostics for FV dynamics ----- if (Atm(mygrid)%flagstruct%print_freq /= -99) then call mpp_clock_begin(id_fv_diag) diff --git a/driver/SHiELD/atmosphere.F90 b/driver/SHiELD/atmosphere.F90 index 3699357cb..869bf2cfe 100644 --- a/driver/SHiELD/atmosphere.F90 +++ b/driver/SHiELD/atmosphere.F90 @@ -34,14 +34,13 @@ module atmosphere_mod use constants_mod, only: cp_air, rdgas, grav, rvgas, kappa, pstd_mks, pi use time_manager_mod, only: time_type, get_time, set_time, operator(+), & operator(-), operator(/), time_type_to_real -use fms_mod, only: file_exist, open_namelist_file, & - close_file, error_mesg, FATAL, & +use fms_mod, only: error_mesg, FATAL, & check_nml_error, stdlog, & write_version_number, & - set_domain, & mpp_clock_id, mpp_clock_begin, & mpp_clock_end, CLOCK_SUBCOMPONENT, & - clock_flag_default, nullify_domain + clock_flag_default +use fms2_io_mod, only: file_exists use mpp_mod, only: mpp_error, stdout, FATAL, WARNING, NOTE, & input_nml_file, mpp_root_pe, & mpp_npes, mpp_pe, mpp_chksum, & @@ -184,7 +183,7 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) !----- initialize FV dynamical core ----- !NOTE do we still need the second file_exist call? - cold_start = (.not.file_exist('INPUT/fv_core.res.nc') .and. .not.file_exist('INPUT/fv_core.res.tile1.nc')) + cold_start = (.not.file_exists('INPUT/fv_core.res.nc') .and. .not.file_exists('INPUT/fv_core.res.tile1.nc')) call fv_control_init( Atm, dt_atmos, mygrid, grids_on_this_pe, p_split ) ! allocates Atm components; sets mygrid @@ -237,7 +236,6 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) ! Allocate grid variables to be used to calculate gradient in 2nd order flux exchange ! This data is only needed for the COARSEST grid. !call switch_current_Atm(Atm(mygrid)) - call set_domain(Atm(mygrid)%domain) allocate(Grid_box%dx ( isc:iec , jsc:jec+1)) allocate(Grid_box%dy ( isc:iec+1, jsc:jec )) @@ -292,12 +290,11 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) Atm(mygrid)%atmos_axes(4), Atm(mygrid)%coarse_graining) endif if (Atm(mygrid)%coarse_graining%write_coarse_restart_files) then - call fv_coarse_restart_init(mygrid, Atm(mygrid)%npz, Atm(mygrid)%flagstruct%nt_prog, & + call fv_coarse_restart_init(Atm(mygrid)%npz, Atm(mygrid)%flagstruct%nt_prog, & Atm(mygrid)%flagstruct%nt_phys, Atm(mygrid)%flagstruct%hydrostatic, & Atm(mygrid)%flagstruct%hybrid_z, Atm(mygrid)%flagstruct%fv_land, & Atm(mygrid)%coarse_graining%write_coarse_dgrid_vel_rst, & Atm(mygrid)%coarse_graining%write_coarse_agrid_vel_rst, & - Atm(mygrid)%coarse_graining%domain, & Atm(mygrid)%coarse_graining%restart) endif @@ -334,7 +331,6 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) if ( Atm(mygrid)%flagstruct%na_init>0 ) then - call nullify_domain ( ) if ( .not. Atm(mygrid)%flagstruct%hydrostatic ) then call prt_maxmin('Before adi: W', Atm(mygrid)%w, isc, iec, jsc, jec, Atm(mygrid)%ng, npz, 1.) endif @@ -350,15 +346,12 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) endif #ifdef DEBUG - call nullify_domain() call fv_diag(Atm(mygrid:mygrid), zvir, Time, -1) if (Atm(mygrid)%coarse_graining%write_coarse_diagnostics) then call fv_coarse_diag(Atm(mygrid:mygrid), fv_time) endif #endif - call set_domain(Atm(mygrid)%domain) - end subroutine atmosphere_init @@ -540,8 +533,6 @@ subroutine atmosphere_end (Time, Grid_box )!rab, Radiation, Physics) !rab type (physics_type), intent(inout) :: Physics ! initialize domains for writing global physics data - call set_domain ( Atm(mygrid)%domain ) - if ( Atm(mygrid)%flagstruct%nudge ) call fv_nwp_nudge_end @@ -549,7 +540,6 @@ subroutine atmosphere_end (Time, Grid_box )!rab, Radiation, Physics) call gfdl_mp_end ( ) endif - call nullify_domain ( ) if (first_diag) then call timing_on('FV_DIAG') call fv_diag(Atm(mygrid:mygrid), zvir, fv_time, Atm(mygrid)%flagstruct%print_freq) @@ -1091,8 +1081,6 @@ subroutine atmosphere_state_update (Time, IPD_Data, Atm_block) if( nq<3 ) call mpp_error(FATAL, 'GFS phys must have 3 interactive tracers') - call set_domain ( Atm(mygrid)%domain ) - call timing_on('GFS_TENDENCIES') call atmos_phys_qdt_diag(Atm(n)%q, Atm(n)%phys_diag, nt_dyn, dt_atmos, .true.) @@ -1270,7 +1258,6 @@ subroutine atmosphere_state_update (Time, IPD_Data, Atm_block) call twoway_nesting(Atm, ngrids, grids_on_this_pe, zvir, fv_time, mygrid) call timing_off('TWOWAY_UPDATE') endif - call nullify_domain() !---- diagnostics for FV dynamics ----- if (Atm(mygrid)%flagstruct%print_freq /= -99) then @@ -1279,7 +1266,6 @@ subroutine atmosphere_state_update (Time, IPD_Data, Atm_block) fv_time = Time_next call get_time (fv_time, seconds, days) - call nullify_domain() call timing_on('FV_DIAG') call fv_diag(Atm(mygrid:mygrid), zvir, fv_time, Atm(mygrid)%flagstruct%print_freq) if (Atm(mygrid)%coarse_graining%write_coarse_diagnostics) then diff --git a/model/fv_arrays.F90 b/model/fv_arrays.F90 index d21a1215a..370d02d49 100644 --- a/model/fv_arrays.F90 +++ b/model/fv_arrays.F90 @@ -24,7 +24,7 @@ module fv_arrays_mod #include use mpp_domains_mod, only: domain2d - use fms_io_mod, only: restart_file_type + use fms2_io_mod, only: FmsNetcdfFile_t, FmsNetcdfDomainFile_t use time_manager_mod, only: time_type use horiz_interp_type_mod, only: horiz_interp_type use mpp_mod, only: mpp_broadcast @@ -999,8 +999,9 @@ module fv_arrays_mod !These are for tracer flux BCs logical :: do_flux_BCs, do_2way_flux_BCs !pi_8, kappa, radius, grav, rdgas use field_manager_mod, only: MODEL_ATMOS - use fms_mod, only: write_version_number, open_namelist_file, & - check_nml_error, close_file, file_exist - use fms_io_mod, only: set_domain + use fms_mod, only: write_version_number, check_nml_error + use fms2_io_mod, only: file_exists use mpp_mod, only: FATAL, mpp_error, mpp_pe, stdlog, & mpp_npes, mpp_get_current_pelist, & input_nml_file, get_unit, WARNING, & @@ -382,7 +381,7 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split) else Atm(n)%nml_filename = 'input.nml' endif - if (.not. file_exist(Atm(n)%nml_filename)) then + if (.not. file_exists(Atm(n)%nml_filename)) then call mpp_error(FATAL, "Could not find nested grid namelist "//Atm(n)%nml_filename) endif enddo @@ -442,7 +441,7 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split) #ifdef INTERNAL_FILE_NML if (this_grid .gt. 1) then write(Atm(this_grid)%nml_filename,'(A4, I2.2)') 'nest', this_grid - if (.not. file_exist('input_'//trim(Atm(this_grid)%nml_filename)//'.nml')) then + if (.not. file_exists('input_'//trim(Atm(this_grid)%nml_filename)//'.nml')) then call mpp_error(FATAL, "Could not find nested grid namelist "//'input_'//trim(Atm(this_grid)%nml_filename)//'.nml') endif else @@ -538,7 +537,6 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split) Atm(this_grid)%layout,Atm(this_grid)%io_layout,Atm(this_grid)%bd,Atm(this_grid)%tile_of_mosaic, & Atm(this_grid)%gridstruct%square_domain,Atm(this_grid)%npes_per_tile,Atm(this_grid)%domain, & Atm(this_grid)%domain_for_coupler,Atm(this_grid)%num_contact,Atm(this_grid)%pelist) - call set_domain(Atm(this_grid)%domain) call broadcast_domains(Atm,Atm(this_grid)%pelist,size(Atm(this_grid)%pelist)) do n=1,ngrids tile_id = mpp_get_tile_id(Atm(n)%domain) @@ -858,16 +856,8 @@ subroutine read_namelist_nest_nml integer :: f_unit, ios, ierr, dum namelist /nest_nml/ dum ! ngrids, ntiles, nest_pes, p_split !emptied lmh 7may2019 -#ifdef INTERNAL_FILE_NML read (input_nml_file,nest_nml,iostat=ios) ierr = check_nml_error(ios,'nest_nml') -#else - f_unit=open_namelist_file() - rewind (f_unit) - read (f_unit,nest_nml,iostat=ios) - ierr = check_nml_error(ios,'nest_nml') - call close_file(f_unit) -#endif if (ierr > 0) then call mpp_error(FATAL, " &nest_nml is depreciated. Please use &fv_nest_nml instead.") endif @@ -880,16 +870,8 @@ subroutine read_namelist_fv_nest_nml namelist /fv_nest_nml/ grid_pes, num_tile_top, tile_coarse, nest_refine, & nest_ioffsets, nest_joffsets, p_split -#ifdef INTERNAL_FILE_NML read (input_nml_file,fv_nest_nml,iostat=ios) ierr = check_nml_error(ios,'fv_nest_nml') -#else - f_unit=open_namelist_file() - rewind (f_unit) - read (f_unit,fv_nest_nml,iostat=ios) - ierr = check_nml_error(ios,'fv_nest_nml') - call close_file(f_unit) -#endif end subroutine read_namelist_fv_nest_nml @@ -901,18 +883,10 @@ subroutine read_namelist_fv_grid_nml character(len=120) :: grid_file = '' namelist /fv_grid_nml/ grid_name, grid_file -#ifdef INTERNAL_FILE_NML ! Read Main namelist read (input_nml_file,fv_grid_nml,iostat=ios) ierr = check_nml_error(ios,'fv_grid_nml') -#else - f_unit=open_namelist_file() - rewind (f_unit) - ! Read Main namelist - read (f_unit,fv_grid_nml,iostat=ios) - ierr = check_nml_error(ios,'fv_grid_nml') - rewind (f_unit) -#endif + call write_version_number ( 'FV_CONTROL_MOD', version ) unit = stdlog() write(unit, nml=fv_grid_nml) @@ -968,19 +942,12 @@ subroutine read_namelist_fv_core_nml(Atm) write_coarse_agrid_vel_rst, write_coarse_dgrid_vel_rst -#ifdef INTERNAL_FILE_NML ! Read FVCORE namelist read (input_nml_file,fv_core_nml,iostat=ios) ierr = check_nml_error(ios,'fv_core_nml') ! Reset input_file_nml to default behavior (CHECK do we still need this???) !call read_input_nml -#else - f_unit = open_namelist_file(Atm%nml_filename) - ! Read FVCORE namelist - read (f_unit,fv_core_nml,iostat=ios) - ierr = check_nml_error(ios,'fv_core_nml') - call close_file(f_unit) -#endif + call write_version_number ( 'FV_CONTROL_MOD', version ) unit = stdlog() write(unit, nml=fv_core_nml) diff --git a/model/fv_regional_bc.F90 b/model/fv_regional_bc.F90 index 5d5d2c874..ea28ebe2a 100644 --- a/model/fv_regional_bc.F90 +++ b/model/fv_regional_bc.F90 @@ -34,7 +34,11 @@ module fv_regional_mod mpp_error ,mpp_pe, mpp_sync, & mpp_npes, mpp_root_pe, mpp_gather, & mpp_get_current_pelist, NULL_PE - use mpp_io_mod + use fms2_io_mod, only: open_file, close_file, register_axis, & + register_variable_attribute, & + register_global_attribute, read_data, & + register_field, FmsNetcdfFile_t, & + FmsNetcdfDomainFile_t, write_data use tracer_manager_mod,only: get_tracer_index use field_manager_mod, only: MODEL_ATMOS use time_manager_mod, only: get_time & @@ -58,7 +62,6 @@ module fv_regional_mod use fv_fill_mod, only: fillz use fv_eta_mod, only: get_eta_level use fms_mod, only: check_nml_error - use fms_io_mod, only: read_data use boundary_mod, only: fv_nest_BC_type_3D private @@ -1184,6 +1187,8 @@ subroutine start_regional_restart(Atm & !*** Local variables !--------------------- ! + type(FmsNetcdfFile_t) :: Grid_input + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist integer :: ierr, ios real, allocatable :: wk2(:,:) ! @@ -1220,7 +1225,13 @@ subroutine start_regional_restart(Atm & allocate (wk2(levp+1,2)) allocate (ak_in(levp+1)) !<-- Save the input vertical structure for allocate (bk_in(levp+1)) ! remapping BC updates during the forecast. - call read_data('INPUT/gfs_ctrl.nc','vcoord',wk2, no_domain=.TRUE.) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if (open_file(Grid_input, 'INPUT/gfs_ctrl.nc', "read", pelist=pes)) then + call read_data(Grid_input,'vcoord',wk2) + call close_file(Grid_input) + endif + deallocate(pes) ak_in(1:levp+1) = wk2(1:levp+1,1) ak_in(1) = 1.e-9 bk_in(1:levp+1) = wk2(1:levp+1,2) @@ -6001,11 +6012,8 @@ subroutine dump_field_3d (domain, name, field, isd, ied, jsd, jed, nlev, stag) integer, intent(IN) :: isd, ied, jsd, jed, nlev integer, intent(IN) :: stag - integer :: unit character(len=128) :: fname - type(axistype) :: x, y, z - type(fieldtype) :: f - type(domain1D) :: xdom, ydom + type(FmsNetcdfDomainFile_t) :: fileobj integer :: nz integer :: is, ie, js, je integer :: isg, ieg, jsg, jeg, nxg, nyg, npx, npy @@ -6015,11 +6023,15 @@ subroutine dump_field_3d (domain, name, field, isd, ied, jsd, jed, nlev, stag) integer, allocatable, dimension(:) :: pelist character(len=1) :: stagname integer :: isection_s, isection_e, jsection_s, jsection_e + character(len=8), dimension(3) :: dim_names_3d + + dim_names_3d(1) = "grid_xt" + dim_names_3d(2) = "grid_yt" + dim_names_3d(3) = "lev" write(fname,"(A,A,A,I1.1,A)") "regional_",name,".tile", 7 , ".nc" write(0,*)'dump_field_3d: file name = |', trim(fname) , '|' - call mpp_get_domain_components( domain, xdom, ydom ) call mpp_get_compute_domain( domain, is, ie, js, je ) call mpp_get_global_domain ( domain, isg, ieg, jsg, jeg, xsize=npx, ysize=npy, position=CENTER ) @@ -6064,32 +6076,50 @@ subroutine dump_field_3d (domain, name, field, isd, ied, jsd, jed, nlev, stag) call mpp_gather(isection_s,isection_e,jsection_s,jsection_e, nz, & pelist, field(isection_s:isection_e,jsection_s:jsection_e,:), glob_field, is_root_pe, halo, halo) - call mpp_open( unit, trim(fname), action=MPP_WRONLY, form=MPP_NETCDF, threading=MPP_SINGLE) - - call mpp_write_meta( unit, x, 'grid_xt', 'km', 'X distance', 'X', domain=xdom, data=(/(i*1.0,i=1,nxg)/) ) - call mpp_write_meta( unit, y, 'grid_yt', 'km', 'Y distance', 'Y', domain=ydom, data=(/(j*1.0,j=1,nyg)/) ) - call mpp_write_meta( unit, z, 'lev', 'km', 'Z distance', data=(/(i*1.0,i=1,nz)/) ) - - call mpp_write_meta( unit, f, (/x,y,z/), name, 'unit', name) - call mpp_write_meta( unit, "stretch_factor", rval=stretch_factor ) - call mpp_write_meta( unit, "target_lon", rval=target_lon ) - call mpp_write_meta( unit, "target_lat", rval=target_lat ) - call mpp_write_meta( unit, "cube_res", ival= cube_res) - call mpp_write_meta( unit, "parent_tile", ival=parent_tile ) - call mpp_write_meta( unit, "refine_ratio", ival=refine_ratio ) - call mpp_write_meta( unit, "istart_nest", ival=istart_nest ) - call mpp_write_meta( unit, "jstart_nest", ival=jstart_nest ) - call mpp_write_meta( unit, "iend_nest", ival=iend_nest ) - call mpp_write_meta( unit, "jend_nest", ival=jend_nest ) - call mpp_write_meta( unit, "ihalo_shift", ival=halo ) - call mpp_write_meta( unit, "jhalo_shift", ival=halo ) - call mpp_write_meta( unit, mpp_get_id(f), "hstagger", cval=stagname ) - call mpp_write( unit, x ) - call mpp_write( unit, y ) - call mpp_write( unit, z ) - call mpp_write( unit, f, glob_field ) - - call mpp_close( unit ) + if (open_file(fileobj, fname, "overwrite", domain)) then + call register_axis(fileobj, "grid_xt", nxg) + call register_field(fileobj, "grid_xt", "double", (/"grid_xt"/)) + call register_variable_attribute(fileobj, "grid_xt", "axis", "X", str_len=1) + call register_variable_attribute(fileobj, "grid_xt", "units", "km", str_len=len("km")) + call register_variable_attribute(fileobj, "grid_xt", "long_name", "X distance", str_len=len("X distance")) + call register_variable_attribute(fileobj, "grid_xt", "cartesian", "X", str_len=1) + call write_data(fileobj, "grid_xt", (/(i*1.0,i=1,nxg)/)) + + call register_axis(fileobj, "grid_yt", nyg) + call register_field(fileobj, "grid_yt", "double", (/"grid_yt"/)) + call register_variable_attribute(fileobj, "grid_yt", "axis", "Y", str_len=1) + call register_variable_attribute(fileobj, "grid_yt", "units", "km", str_len=len("km")) + call register_variable_attribute(fileobj, "grid_yt", "long_name", "Y distance", str_len=len("Y distance")) + call register_variable_attribute(fileobj, "grid_yt", "cartesian", "Y", str_len=1) + call write_data(fileobj, "grid_yt", (/(j*1.0,j=1,nyg)/)) + + call register_axis(fileobj, "lev", nz) + call register_field(fileobj, "lev", "double", (/"lev"/)) + call register_variable_attribute(fileobj, "lev", "axis", "Z", str_len=1) + call register_variable_attribute(fileobj, "lev", "units", "km", str_len=len("km")) + call register_variable_attribute(fileobj, "lev", "long_name", "Z distance", str_len=len("Z distance")) + call write_data(fileobj, "lev", (/(i*1.0,i=1,nz)/)) + + call register_global_attribute(fileobj, "stretch_factor", stretch_factor ) + call register_global_attribute(fileobj, "target_lon", target_lon ) + call register_global_attribute(fileobj, "target_lat", target_lat ) + call register_global_attribute(fileobj, "cube_res", cube_res) + call register_global_attribute(fileobj, "parent_tile", parent_tile ) + call register_global_attribute(fileobj, "refine_ratio", refine_ratio ) + call register_global_attribute(fileobj, "istart_nest", istart_nest ) + call register_global_attribute(fileobj, "jstart_nest", jstart_nest ) + call register_global_attribute(fileobj, "iend_nest", iend_nest ) + call register_global_attribute(fileobj, "jend_nest", jend_nest ) + call register_global_attribute(fileobj, "ihalo_shift", halo ) + call register_global_attribute(fileobj, "jhalo_shift", halo ) + call register_global_attribute(fileobj, "hstagger", stagname ) + + call register_field(fileobj, name, "double", dim_names_3d) + + call write_data(fileobj, name, glob_field) + + call close_file(fileobj) + endif end subroutine dump_field_3d @@ -6101,11 +6131,8 @@ subroutine dump_field_2d (domain, name, field, isd, ied, jsd, jed, stag) integer, intent(IN) :: isd, ied, jsd, jed integer, intent(IN) :: stag - integer :: unit character(len=128) :: fname - type(axistype) :: x, y - type(fieldtype) :: f - type(domain1D) :: xdom, ydom + type(FmsNetcdfDomainFile_t) :: fileobj integer :: is, ie, js, je integer :: isg, ieg, jsg, jeg, nxg, nyg, npx, npy integer :: i, j, halo, iext, jext @@ -6115,10 +6142,15 @@ subroutine dump_field_2d (domain, name, field, isd, ied, jsd, jed, stag) character(len=1) :: stagname integer :: isection_s, isection_e, jsection_s, jsection_e + character(len=8), dimension(3) :: dim_names_3d + + dim_names_3d(1) = "grid_xt" + dim_names_3d(2) = "grid_yt" + dim_names_3d(3) = "lev" + write(fname,"(A,A,A,I1.1,A)") "regional_",name,".tile", 7 , ".nc" ! write(0,*)'dump_field_3d: file name = |', trim(fname) , '|' - call mpp_get_domain_components( domain, xdom, ydom ) call mpp_get_compute_domain( domain, is, ie, js, je ) call mpp_get_global_domain ( domain, isg, ieg, jsg, jeg, xsize=npx, ysize=npy, position=CENTER ) @@ -6162,30 +6194,43 @@ subroutine dump_field_2d (domain, name, field, isd, ied, jsd, jed, stag) call mpp_gather(isection_s,isection_e,jsection_s,jsection_e, & pelist, field(isection_s:isection_e,jsection_s:jsection_e), glob_field, is_root_pe, halo, halo) - call mpp_open( unit, trim(fname), action=MPP_WRONLY, form=MPP_NETCDF, threading=MPP_SINGLE) - - call mpp_write_meta( unit, x, 'grid_xt', 'km', 'X distance', 'X', domain=xdom, data=(/(i*1.0,i=1,nxg)/) ) - call mpp_write_meta( unit, y, 'grid_yt', 'km', 'Y distance', 'Y', domain=ydom, data=(/(j*1.0,j=1,nyg)/) ) - - call mpp_write_meta( unit, f, (/x,y/), name, 'unit', name) - call mpp_write_meta( unit, "stretch_factor", rval=stretch_factor ) - call mpp_write_meta( unit, "target_lon", rval=target_lon ) - call mpp_write_meta( unit, "target_lat", rval=target_lat ) - call mpp_write_meta( unit, "cube_res", ival= cube_res) - call mpp_write_meta( unit, "parent_tile", ival=parent_tile ) - call mpp_write_meta( unit, "refine_ratio", ival=refine_ratio ) - call mpp_write_meta( unit, "istart_nest", ival=istart_nest ) - call mpp_write_meta( unit, "jstart_nest", ival=jstart_nest ) - call mpp_write_meta( unit, "iend_nest", ival=iend_nest ) - call mpp_write_meta( unit, "jend_nest", ival=jend_nest ) - call mpp_write_meta( unit, "ihalo_shift", ival=halo ) - call mpp_write_meta( unit, "jhalo_shift", ival=halo ) - call mpp_write_meta( unit, mpp_get_id(f), "hstagger", cval=stagname ) - call mpp_write( unit, x ) - call mpp_write( unit, y ) - call mpp_write( unit, f, glob_field ) - - call mpp_close( unit ) + if (open_file(fileobj, fname, "overwrite", domain)) then + call register_axis(fileobj, "grid_xt", nxg) + call register_field(fileobj, "grid_xt", "double", (/"grid_xt"/)) + call register_variable_attribute(fileobj, "grid_xt", "axis", "X", str_len=1) + call register_variable_attribute(fileobj, "grid_xt", "units", "km", str_len=len("km")) + call register_variable_attribute(fileobj, "grid_xt", "long_name", "X distance", str_len=len("X distance")) + call register_variable_attribute(fileobj, "grid_xt", "cartesian", "X", str_len=1) + call write_data(fileobj, "grid_xt", (/(i*1.0,i=1,nxg)/)) + + call register_axis(fileobj, "grid_yt", nyg) + call register_field(fileobj, "grid_yt", "double", (/"grid_yt"/)) + call register_variable_attribute(fileobj, "grid_yt", "axis", "Y", str_len=1) + call register_variable_attribute(fileobj, "grid_yt", "units", "km", str_len=len("km")) + call register_variable_attribute(fileobj, "grid_yt", "long_name", "Y distance", str_len=len("Y distance")) + call register_variable_attribute(fileobj, "grid_yt", "cartesian", "Y", str_len=1) + call write_data(fileobj, "grid_yt", (/(j*1.0,j=1,nyg)/)) + + call register_global_attribute(fileobj, "stretch_factor", stretch_factor ) + call register_global_attribute(fileobj, "target_lon", target_lon ) + call register_global_attribute(fileobj, "target_lat", target_lat ) + call register_global_attribute(fileobj, "cube_res", cube_res) + call register_global_attribute(fileobj, "parent_tile", parent_tile ) + call register_global_attribute(fileobj, "refine_ratio", refine_ratio ) + call register_global_attribute(fileobj, "istart_nest", istart_nest ) + call register_global_attribute(fileobj, "jstart_nest", jstart_nest ) + call register_global_attribute(fileobj, "iend_nest", iend_nest ) + call register_global_attribute(fileobj, "jend_nest", jend_nest ) + call register_global_attribute(fileobj, "ihalo_shift", halo ) + call register_global_attribute(fileobj, "jhalo_shift", halo ) + call register_global_attribute(fileobj, "hstagger", stagname ) + + call register_field(fileobj, name, "double", dim_names_3d) + + call write_data(fileobj, name, glob_field) + + call close_file(fileobj) + endif end subroutine dump_field_2d diff --git a/tools/coarse_grained_restart_files.F90 b/tools/coarse_grained_restart_files.F90 index c11bbdfe1..c94bcf6bf 100644 --- a/tools/coarse_grained_restart_files.F90 +++ b/tools/coarse_grained_restart_files.F90 @@ -7,11 +7,12 @@ module coarse_grained_restart_files_mod remap_edges_along_y, vertically_remap_field use constants_mod, only: GRAV, RDGAS, RVGAS use field_manager_mod, only: MODEL_ATMOS - use fms_io_mod, only: register_restart_field, save_restart + use fms2_io_mod, only: register_restart_field, write_restart, open_file, close_file use fv_arrays_mod, only: coarse_restart_type, fv_atmos_type - use mpp_domains_mod, only: domain2d, EAST, NORTH, mpp_update_domains + use mpp_domains_mod, only: domain2d, EAST, NORTH, CENTER, mpp_update_domains use mpp_mod, only: FATAL, mpp_error use tracer_manager_mod, only: get_tracer_names, get_tracer_index, set_tracer_profile + use fv_io_mod, only: fv_io_register_axis implicit none private @@ -25,14 +26,13 @@ module coarse_grained_restart_files_mod contains - subroutine fv_coarse_restart_init(tile_count, nz, nt_prog, & + subroutine fv_coarse_restart_init(nz, nt_prog, & nt_phys, hydrostatic, hybrid_z, fv_land, & write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst, & - coarse_domain, restart) - integer, intent(in) :: tile_count, nz, nt_prog, nt_phys + restart) + integer, intent(in) :: nz, nt_prog, nt_phys logical, intent(in) :: hydrostatic, hybrid_z, fv_land logical, intent(in) :: write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst - type(domain2d), intent(inout) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart call get_fine_array_bounds(is, ie, js, je) @@ -45,9 +45,6 @@ subroutine fv_coarse_restart_init(tile_count, nz, nt_prog, & call allocate_coarse_restart_type(hydrostatic, hybrid_z, & fv_land, write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst, & restart) - call register_coarse_restart_files(tile_count, hydrostatic, & - hybrid_z, fv_land, write_coarse_dgrid_vel_rst, & - write_coarse_agrid_vel_rst, coarse_domain, restart) end subroutine fv_coarse_restart_init subroutine fv_io_write_restart_coarse(Atm, timestamp) @@ -56,13 +53,41 @@ subroutine fv_io_write_restart_coarse(Atm, timestamp) integer :: tile_count, n_tiles + call register_coarse_restart_files(Atm%flagstruct%hydrostatic, & + Atm%flagstruct%hybrid_z, Atm%flagstruct%fv_land, & + Atm%coarse_graining%write_coarse_dgrid_vel_rst, & + Atm%coarse_graining%write_coarse_agrid_vel_rst, & + Atm%coarse_graining%domain, & + Atm%coarse_graining%restart, timestamp) + call coarse_grain_restart_data(Atm) - call save_restart(Atm%coarse_graining%restart%fv_core_coarse, timestamp) - call save_restart(Atm%coarse_graining%restart%fv_tracer_coarse, timestamp) - call save_restart(Atm%coarse_graining%restart%fv_srf_wnd_coarse, timestamp) + + if (Atm%coarse_graining%restart%fv_core_coarse_is_open) then + call write_restart(Atm%coarse_graining%restart%fv_core_coarse) + call close_file(Atm%coarse_graining%restart%fv_core_coarse) + Atm%coarse_graining%restart%fv_core_coarse_is_open=.false. + endif + if (Atm%coarse_graining%restart%fv_tracer_coarse_is_open) then + call write_restart(Atm%coarse_graining%restart%fv_tracer_coarse) + call close_file(Atm%coarse_graining%restart%fv_tracer_coarse) + Atm%coarse_graining%restart%fv_tracer_coarse_is_open=.false. + endif + if (Atm%coarse_graining%restart%fv_srf_wnd_coarse_is_open) then + call write_restart(Atm%coarse_graining%restart%fv_srf_wnd_coarse) + call close_file(Atm%coarse_graining%restart%fv_srf_wnd_coarse) + Atm%coarse_graining%restart%fv_srf_wnd_coarse_is_open=.false. + endif if (Atm%flagstruct%fv_land) then - call save_restart(Atm%coarse_graining%restart%mg_drag_coarse, timestamp) - call save_restart(Atm%coarse_graining%restart%fv_land_coarse, timestamp) + if (Atm%coarse_graining%restart%mg_drag_coarse_is_open) then + call write_restart(Atm%coarse_graining%restart%mg_drag_coarse) + call close_file(Atm%coarse_graining%restart%mg_drag_coarse) + Atm%coarse_graining%restart%mg_drag_coarse_is_open=.false. + endif + if (Atm%coarse_graining%restart%fv_land_coarse_is_open) then + call write_restart(Atm%coarse_graining%restart%fv_land_coarse) + call close_file(Atm%coarse_graining%restart%fv_land_coarse) + Atm%coarse_graining%restart%fv_land_coarse_is_open=.false. + endif endif end subroutine fv_io_write_restart_coarse @@ -123,148 +148,228 @@ subroutine deallocate_coarse_restart_type(restart) if (allocated(restart%ze0)) deallocate(restart%ze0) end subroutine deallocate_coarse_restart_type - subroutine register_coarse_restart_files(tile_count, hydrostatic, & + subroutine register_coarse_restart_files(hydrostatic, & hybrid_z, fv_land, write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst, & - coarse_domain, restart) - integer, intent(in) :: tile_count + coarse_domain, restart, timestamp) logical, intent(in) :: hydrostatic, hybrid_z, fv_land logical, intent(in) :: write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart + character(len=*), optional, intent(in) :: timestamp - call register_fv_core_coarse(tile_count, hydrostatic, hybrid_z, & + call register_fv_core_coarse(hydrostatic, hybrid_z, & write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst, & - coarse_domain, restart) - call register_fv_tracer_coarse(tile_count, coarse_domain, restart) - call register_fv_srf_wnd_coarse(tile_count, coarse_domain, restart) + coarse_domain, restart, timestamp) + call register_fv_tracer_coarse(coarse_domain, restart, timestamp) + call register_fv_srf_wnd_coarse(coarse_domain, restart, timestamp) if (fv_land) then - call register_mg_drag_coarse(tile_count, coarse_domain, restart) - call register_fv_land_coarse(tile_count, coarse_domain, restart) + call register_mg_drag_coarse(coarse_domain, restart, timestamp) + call register_fv_land_coarse(coarse_domain, restart, timestamp) endif end subroutine register_coarse_restart_files - subroutine register_fv_core_coarse(tile_count, hydrostatic, hybrid_z, & + subroutine register_fv_core_coarse(hydrostatic, hybrid_z, & write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst, coarse_domain, & - restart) - integer, intent(in) :: tile_count + restart, timestamp) logical, intent(in) :: hydrostatic, hybrid_z logical, intent(in) :: write_coarse_dgrid_vel_rst, write_coarse_agrid_vel_rst type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart + character(len=*), optional, intent(in) :: timestamp + character(len=8), dimension(4) :: dim_names_4d, dim_names_4d2, dim_names_4d3 + character(len=8), dimension(3) :: dim_names_3d character(len=64) :: filename - integer :: id_restart - - filename = 'fv_core_coarse.res.nc' - - if (write_coarse_dgrid_vel_rst) then - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'u', restart%u, domain=coarse_domain, position=NORTH, & - tile_count=tile_count) - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'v', restart%v, domain=coarse_domain, position=EAST, & - tile_count=tile_count) - endif - - if (write_coarse_agrid_vel_rst) then - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'ua', restart%ua, domain=coarse_domain, tile_count=tile_count) - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'va', restart%va, domain=coarse_domain, tile_count=tile_count) + integer, dimension(1) :: zsize + + dim_names_4d(1) = "xaxis_1" + dim_names_4d(2) = "yaxis_1" + dim_names_4d(3) = "zaxis_1" + dim_names_4d(4) = "Time" + dim_names_4d2 = dim_names_4d + dim_names_4d2(1) = "xaxis_2" + dim_names_4d2(2) = "yaxis_2" + dim_names_4d3 = dim_names_4d + dim_names_4d3(2) = "yaxis_2" + dim_names_3d(1) = "xaxis_1" + dim_names_3d(2) = "yaxis_2" + dim_names_3d(3) = "Time" + + if (present(timestamp)) then + filename = 'RESTART/'//trim(timestamp)//'.fv_core_coarse.res.nc' + else + filename = 'RESTART/fv_core_coarse.res.nc' endif - if (.not. hydrostatic) then - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'W', restart%w, domain=coarse_domain, mandatory=.false., tile_count=tile_count) - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'DZ', restart%delz, domain=coarse_domain, mandatory=.false., tile_count=tile_count) - if (hybrid_z) then - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'ZE0', restart%ze0, domain=coarse_domain, mandatory=.false., tile_count=tile_count) - endif + restart%fv_core_coarse_is_open = open_file(restart%fv_core_coarse, filename, & + "overwrite", coarse_domain, is_restart=.true.) + if (restart%fv_core_coarse_is_open) then + zsize = (/size(restart%u,3)/) + call fv_io_register_axis(restart%fv_core_coarse, numx=2 ,numy=2, xpos=(/CENTER, EAST/) ,ypos=(/NORTH, CENTER/) ,numz=1, zsize=zsize) + + if (write_coarse_dgrid_vel_rst) then + call register_restart_field(restart%fv_core_coarse, & + 'u', restart%u, dim_names_4d) + call register_restart_field(restart%fv_core_coarse, & + 'v', restart%v, dim_names_4d2) + endif + + if (write_coarse_agrid_vel_rst) then + call register_restart_field(restart%fv_core_coarse, & + 'ua', restart%ua, dim_names_4d3) + call register_restart_field(restart%fv_core_coarse, & + 'va', restart%va, dim_names_4d3) + endif + + if (.not. hydrostatic) then + call register_restart_field(restart%fv_core_coarse, & + 'W', restart%w, dim_names_4d3, is_optional=.true.) + call register_restart_field(restart%fv_core_coarse, & + 'DZ', restart%delz, dim_names_4d3, is_optional=.true.) + if (hybrid_z) then + call register_restart_field(restart%fv_core_coarse, & + 'ZE0', restart%ze0, dim_names_4d3, is_optional=.false.) + endif + endif + + call register_restart_field(restart%fv_core_coarse, & + 'T', restart%pt, dim_names_4d3) + call register_restart_field(restart%fv_core_coarse, & + 'delp', restart%delp, dim_names_4d3) + call register_restart_field(restart%fv_core_coarse, & + 'phis', restart%phis, dim_names_3d) endif - - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'T', restart%pt, domain=coarse_domain, tile_count=tile_count) - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'delp', restart%delp, domain=coarse_domain, tile_count=tile_count) - id_restart = register_restart_field(restart%fv_core_coarse, & - filename, 'phis', restart%phis, domain=coarse_domain, tile_count=tile_count) end subroutine register_fv_core_coarse - subroutine register_fv_tracer_coarse(tile_count, coarse_domain, restart) - integer, intent(in) :: tile_count + subroutine register_fv_tracer_coarse(coarse_domain, restart, timestamp) type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart + character(len=*), optional, intent(in) :: timestamp + character(len=8), dimension(4) :: dim_names_4d character(len=64) :: filename, tracer_name - integer :: id_restart, n_tracer + integer :: n_tracer + integer, dimension(1) :: zsize - filename = 'fv_tracer_coarse.res.nc' + dim_names_4d(1) = "xaxis_1" + dim_names_4d(2) = "yaxis_1" + dim_names_4d(3) = "zaxis_1" + dim_names_4d(4) = "Time" - do n_tracer = 1, n_prognostic_tracers - call get_tracer_names(MODEL_ATMOS, n_tracer, tracer_name) - call set_tracer_profile(MODEL_ATMOS, n_tracer, restart%q(:,:,:,n_tracer)) - id_restart = register_restart_field(restart%fv_tracer_coarse, & - filename, tracer_name, restart%q(:,:,:,n_tracer), domain=coarse_domain, & - mandatory=.false., tile_count=tile_count) - enddo + if (present(timestamp)) then + filename = 'RESTART/'//trim(timestamp)//'.fv_tracer_coarse.res.nc' + else + filename = 'RESTART/fv_tracer_coarse.res.nc' + endif - do n_tracer = n_prognostic_tracers + 1, n_tracers - call get_tracer_names(MODEL_ATMOS, n_tracer, tracer_name) - call set_tracer_profile(MODEL_ATMOS, n_tracer, restart%qdiag(:,:,:,n_tracer)) - id_restart = register_restart_field(restart%fv_tracer_coarse, & - filename, tracer_name, restart%qdiag(:,:,:,n_tracer), domain=coarse_domain, & - mandatory=.false., tile_count=tile_count) - enddo + restart%fv_tracer_coarse_is_open = open_file(restart%fv_tracer_coarse, filename, & + "overwrite", coarse_domain, is_restart=.true.) + if (restart%fv_tracer_coarse_is_open) then + zsize=(/size(restart%q,3)/) + call fv_io_register_axis(restart%fv_tracer_coarse, numx=1 ,numy=1, xpos=(/CENTER/) ,ypos=(/CENTER/) ,numz=1, zsize=zsize) + + do n_tracer = 1, n_prognostic_tracers + call get_tracer_names(MODEL_ATMOS, n_tracer, tracer_name) + call set_tracer_profile(MODEL_ATMOS, n_tracer, restart%q(:,:,:,n_tracer)) + call register_restart_field(restart%fv_tracer_coarse, & + tracer_name, restart%q(:,:,:,n_tracer), dim_names_4d, & + is_optional=.true.) + enddo + + do n_tracer = n_prognostic_tracers + 1, n_tracers + call get_tracer_names(MODEL_ATMOS, n_tracer, tracer_name) + call set_tracer_profile(MODEL_ATMOS, n_tracer, restart%qdiag(:,:,:,n_tracer)) + call register_restart_field(restart%fv_tracer_coarse, & + tracer_name, restart%qdiag(:,:,:,n_tracer), dim_names_4d, & + is_optional=.true.) + enddo + endif end subroutine register_fv_tracer_coarse - subroutine register_fv_srf_wnd_coarse(tile_count, coarse_domain, restart) - integer, intent(in) :: tile_count + subroutine register_fv_srf_wnd_coarse(coarse_domain, restart, timestamp) type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart + character(len=*), optional, intent(in) :: timestamp + character(len=8), dimension(3) :: dim_names_3d character(len=64) :: filename - integer :: id_restart - filename = 'fv_srf_wnd_coarse.res.nc' + dim_names_3d(1) = "xaxis_1" + dim_names_3d(2) = "yaxis_1" + dim_names_3d(3) = "Time" + + if (present(timestamp)) then + filename = 'RESTART/'//trim(timestamp)//'.fv_srf_wnd_coarse.res.nc' + else + filename = 'RESTART/fv_srf_wnd_coarse.res.nc' + endif + + restart%fv_srf_wnd_coarse_is_open = open_file(restart%fv_srf_wnd_coarse, filename, & + "overwrite", coarse_domain, is_restart=.true.) + if (restart%fv_srf_wnd_coarse_is_open) then + call fv_io_register_axis(restart%fv_srf_wnd_coarse, numx=1 ,numy=1, xpos=(/CENTER/) ,ypos=(/CENTER/)) - id_restart = register_restart_field(restart%fv_srf_wnd_coarse, & - filename, 'u_srf', restart%u_srf, domain=coarse_domain, & - tile_count=tile_count) - id_restart = register_restart_field(restart%fv_srf_wnd_coarse, & - filename, 'v_srf', restart%v_srf, domain=coarse_domain, & - tile_count=tile_count) + call register_restart_field(restart%fv_srf_wnd_coarse, & + 'u_srf', restart%u_srf, dim_names_3d) + call register_restart_field(restart%fv_srf_wnd_coarse, & + 'v_srf', restart%v_srf, dim_names_3d) + endif end subroutine register_fv_srf_wnd_coarse - subroutine register_mg_drag_coarse(tile_count, coarse_domain, restart) - integer, intent(in) :: tile_count + subroutine register_mg_drag_coarse(coarse_domain, restart, timestamp) type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(out) :: restart + character(len=*), optional, intent(in) :: timestamp + character(len=8), dimension(3) :: dim_names_3d character(len=64) :: filename - integer :: id_restart - filename = 'mg_drag_coarse.res.nc' + dim_names_3d(1) = "xaxis_1" + dim_names_3d(2) = "yaxis_1" + dim_names_3d(3) = "Time" + + if (present(timestamp)) then + filename = 'RESTART/'//trim(timestamp)//'.mg_drag_coarse.res.nc' + else + filename = 'RESTART/mg_drag_coarse.res.nc' + endif + + restart%mg_drag_coarse_is_open = open_file(restart%mg_drag_coarse, filename, & + "overwrite", coarse_domain, is_restart=.true.) + if (restart%mg_drag_coarse_is_open) then + call fv_io_register_axis(restart%mg_drag_coarse, numx=1, numy=1, xpos=(/CENTER/), ypos=(/CENTER/)) - id_restart = register_restart_field(restart%mg_drag_coarse, & - filename, 'ghprime', restart%sgh, domain=coarse_domain, & - tile_count=tile_count) + call register_restart_field(restart%mg_drag_coarse, & + 'ghprime', restart%sgh, dim_names_3d) + endif end subroutine register_mg_drag_coarse - subroutine register_fv_land_coarse(tile_count, coarse_domain, restart) - integer, intent(in) :: tile_count + subroutine register_fv_land_coarse(coarse_domain, restart, timestamp) type(domain2d), intent(in) :: coarse_domain type(coarse_restart_type), intent(inout) :: restart + character(len=*), optional, intent(in) :: timestamp + character(len=8), dimension(3) :: dim_names_3d character(len=64) :: filename - integer :: id_restart - filename = 'fv_land_coarse.res.nc' + dim_names_3d(1) = "xaxis_1" + dim_names_3d(2) = "yaxis_1" + dim_names_3d(3) = "Time" + + if (present(timestamp)) then + filename = 'RESTART/'//trim(timestamp)//'.fv_land_coarse.res.nc' + else + filename = 'RESTART/fv_land_coarse.res.nc' + endif + + restart%fv_land_coarse_is_open = open_file(restart%fv_land_coarse, filename, & + "overwrite", coarse_domain, is_restart=.true.) + if (restart%fv_land_coarse_is_open) then + call fv_io_register_axis(restart%fv_land_coarse, numx=1, numy=1, xpos=(/CENTER/), ypos=(/CENTER/)) - id_restart = register_restart_field(restart%fv_land_coarse, & - filename, 'oro', restart%oro, domain=coarse_domain, & - tile_count=tile_count) + call register_restart_field(restart%fv_land_coarse, & + 'oro', restart%oro, dim_names_3d) + endif end subroutine register_fv_land_coarse subroutine coarse_grain_restart_data(Atm) diff --git a/tools/coarse_graining.F90 b/tools/coarse_graining.F90 index 347ef29a7..4584d4f0e 100644 --- a/tools/coarse_graining.F90 +++ b/tools/coarse_graining.F90 @@ -1,6 +1,6 @@ module coarse_graining_mod - use fms_mod, only: check_nml_error, close_file, open_namelist_file + use fms_mod, only: check_nml_error use mpp_domains_mod, only: domain2d, mpp_define_io_domain, mpp_define_mosaic, mpp_get_compute_domain use fv_mapz_mod, only: mappm use mpp_mod, only: FATAL, input_nml_file, mpp_error, mpp_npes diff --git a/tools/external_ic.F90 b/tools/external_ic.F90 index 70176a521..7db34ee42 100644 --- a/tools/external_ic.F90 +++ b/tools/external_ic.F90 @@ -27,14 +27,13 @@ module external_ic_mod use external_sst_mod, only: i_sst, j_sst, sst_ncep - use fms_mod, only: file_exist, read_data, field_exist, write_version_number - use fms_mod, only: open_namelist_file, check_nml_error, close_file - use fms_mod, only: get_mosaic_tile_file, read_data, error_mesg - use fms_io_mod, only: get_tile_string, field_size, free_restart_type - use fms_io_mod, only: restart_file_type, register_restart_field - use fms_io_mod, only: save_restart, restore_state, set_filename_appendix, get_global_att_value + use fms_mod, only: write_version_number, check_nml_error + use fms2_io_mod, only: file_exists, open_file, close_file, read_data, variable_exists, & + get_variable_size, get_global_attribute, global_att_exists, & + FmsNetcdfFile_t, FmsNetcdfDomainFile_t, read_restart, & + register_restart_field, register_axis use mpp_mod, only: mpp_error, FATAL, NOTE, mpp_pe, mpp_root_pe - use mpp_mod, only: stdlog, input_nml_file + use mpp_mod, only: stdlog, input_nml_file, mpp_npes, mpp_get_current_pelist use mpp_parameter_mod, only: AGRID_PARAM=>AGRID use mpp_domains_mod, only: mpp_get_tile_id, domain2d, mpp_update_domains, NORTH, EAST use tracer_manager_mod, only: get_tracer_names, get_number_tracers, get_tracer_index @@ -221,10 +220,9 @@ end subroutine get_external_ic subroutine get_cubed_sphere_terrain( Atm, fv_domain ) type(fv_atmos_type), intent(inout), target :: Atm type(domain2d), intent(inout) :: fv_domain + type(FmsNetcdfDomainFile_t) :: Fv_core integer :: tile_id(1) character(len=64) :: fname - character(len=7) :: gn - integer :: n=1 integer :: jbeg, jend real ftop real, allocatable :: g_dat2(:,:,:) @@ -244,22 +242,15 @@ subroutine get_cubed_sphere_terrain( Atm, fv_domain ) jed = Atm%bd%jed ng = Atm%bd%ng - if (Atm%grid_number > 1) then - !write(gn,'(A2, I1)') ".g", Atm%grid_number - write(gn,'(A5, I2.2)') ".nest", Atm%grid_number - else - gn = '' - end if - tile_id = mpp_get_tile_id( fv_domain ) - call get_tile_string(fname, 'INPUT/fv_core.res'//trim(gn)//'.tile', tile_id(n), '.nc' ) + fname = 'INPUT/fv_core.res.nc' call mpp_error(NOTE, 'external_ic: looking for '//fname) - if( file_exist(fname) ) then - call read_data(fname, 'phis', Atm%phis(is:ie,js:je), & - domain=fv_domain, tile_count=n) + if( open_file(Fv_core, fname, "read", fv_domain, is_restart=.true.) ) then + call read_data(Fv_core, 'phis', Atm%phis(is:ie,js:je)) + call close_file(Fv_core) else call mpp_error(NOTE, fname//' not found; generating terrain from USGS data') call surfdrv( Atm%npx, Atm%npy, Atm%gridstruct%grid_64, Atm%gridstruct%agrid_64, & @@ -326,19 +317,25 @@ subroutine get_nggps_ic (Atm, fv_domain) integer :: is, ie, js, je integer :: isd, ied, jsd, jed integer :: ios, ierr, unit, id_res - type (restart_file_type) :: ORO_restart, SFC_restart, GFS_restart - character(len=6) :: gn, stile_name + type(FmsNetcdfDomainFile_t) :: ORO_restart, SFC_restart, GFS_restart + type(FmsNetcdfFile_t) :: Gfs_ctl + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist + character(len=8), dimension(2) :: dim_names_2d + character(len=8), dimension(3) :: dim_names_3d, dim_names_3d2, dim_names_3d3, dim_names_3d4 + character(len=6) :: stile_name character(len=64) :: tracer_name - character(len=64) :: fn_gfs_ctl = 'gfs_ctrl.nc' - character(len=64) :: fn_gfs_ics = 'gfs_data.nc' - character(len=64) :: fn_sfc_ics = 'sfc_data.nc' - character(len=64) :: fn_oro_ics = 'oro_data.nc' + character(len=64) :: fn_gfs_ctl = 'INPUT/gfs_ctrl.nc' + character(len=64) :: fn_gfs_ics = 'INPUT/gfs_data.nc' + character(len=64) :: fn_sfc_ics = 'INPUT/sfc_data.nc' + character(len=64) :: fn_oro_ics = 'INPUT/oro_data.nc' logical :: remap logical :: filtered_terrain = .true. logical :: gfs_dwinds = .true. integer :: levp = 64 logical :: checker_tr = .false. integer :: nt_checker = 0 + character(len=20) :: suffix + character(len=1) :: tile_num real(kind=R_GRID), dimension(2):: p1, p2, p3 real(kind=R_GRID), dimension(3):: e1, e2, ex, ey integer:: i,j,k,nts, ks @@ -351,15 +348,8 @@ subroutine get_nggps_ic (Atm, fv_domain) call mpp_error(NOTE,'Using external_IC::get_nggps_ic which is valid only for data which has been & &horizontally interpolated to the current cubed-sphere grid') -#ifdef INTERNAL_FILE_NML read (input_nml_file,external_ic_nml,iostat=ios) ierr = check_nml_error(ios,'external_ic_nml') -#else - unit=open_namelist_file() - read (unit,external_ic_nml,iostat=ios) - ierr = check_nml_error(ios,'external_ic_nml') - call close_file(unit) -#endif unit = stdlog() call write_version_number ( 'EXTERNAL_IC_MOD::get_nggps_ic', version ) @@ -396,22 +386,18 @@ subroutine get_nggps_ic (Atm, fv_domain) call get_number_tracers(MODEL_ATMOS, num_tracers=ntracers, num_prog=ntprog) ntdiag = ntracers-ntprog -!--- set the 'nestXX' appendix for all files using fms_io - if (Atm%grid_number > 1) then - write(gn,'(A4, I2.2)') "nest", Atm%grid_number + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Gfs_ctl, fn_gfs_ctl, "read", pelist=pes) ) then +!--- read in the number of tracers in the NCEP NGGPS ICs + call read_data (Gfs_ctl, 'ntrac', ntrac) + call close_file(Gfs_ctl) else - gn = '' - end if - call set_filename_appendix('') - -!--- test for existence of the GFS control file - if (.not. file_exist('INPUT/'//trim(fn_gfs_ctl), no_domain=.TRUE.)) then call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: file '//trim(fn_gfs_ctl)//' for NGGPS IC does not exist') endif + deallocate(pes) call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using control file '//trim(fn_gfs_ctl)//' for NGGPS IC') -!--- read in the number of tracers in the NCEP NGGPS ICs - call read_data ('INPUT/'//trim(fn_gfs_ctl), 'ntrac', ntrac, no_domain=.TRUE.) if (ntrac > ntracers) call mpp_error(FATAL,'==> External_ic::get_nggps_ic: more NGGPS tracers & &than defined in field_table '//trim(fn_gfs_ctl)//' for NGGPS IC') @@ -421,11 +407,10 @@ subroutine get_nggps_ic (Atm, fv_domain) !--- read in the number of levp - call open_ncfile( 'INPUT/'//trim(fn_gfs_ctl), ncid ) ! open the file + call open_ncfile(fn_gfs_ctl, ncid ) ! open the file call get_ncdim1( ncid, 'levsp', levsp ) call close_ncfile( ncid ) - ! read in gfs_data. If levp = 66, read only the lowest 65 level if (levsp .eq. 66) then call mpp_error(NOTE,'==> External_ic::get_nggps_ic: Correcting BAD IC') @@ -445,55 +430,56 @@ subroutine get_nggps_ic (Atm, fv_domain) enddo endif -!--- test for existence of the GFS orography and surface files - if (.not. file_exist('INPUT/'//trim(fn_oro_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_oro_ics)//' for NGGPS IC does not exist') - endif - call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_oro_ics)//' for NGGPS IC') - - if (.not. file_exist('INPUT/'//trim(fn_sfc_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_sfc_ics)//' for NGGPS IC does not exist') - endif - call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_sfc_ics)//' for NGGPS IC') - + ! set dimensions for register restart + dim_names_2d(1) = "lat" + dim_names_2d(2) = "lon" !--- read in surface temperature (k) and land-frac ! surface skin temperature - id_res = register_restart_field (SFC_restart, fn_sfc_ics, 'tsea', Atm%ts, domain=Atm%domain) + if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(SFC_restart, "lat", "y") + call register_axis(SFC_restart, "lon", "x") + call register_restart_field(SFC_restart, 'tsea', Atm%ts, dim_names_2d) + call read_restart(SFC_restart) + call close_file(SFC_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_sfc_ics)//' for NGGPS IC does not exist') + endif + call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_sfc_ics)//' for NGGPS IC') ! terrain surface height -- (needs to be transformed into phis = zs*grav) - if (filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_filt', Atm%phis, domain=Atm%domain) - elseif (.not. filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_raw', Atm%phis, domain=Atm%domain) - endif + if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(ORO_restart, "lat", "y") + call register_axis(ORO_restart, "lon", "x") + if (filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_filt', Atm%phis, dim_names_2d) + elseif (.not. filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_raw', Atm%phis, dim_names_2d) + endif - if ( Atm%flagstruct%full_zs_filter) then - allocate (oro_g(isd:ied,jsd:jed)) - oro_g = 0. - ! land-frac - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'land_frac', oro_g, domain=Atm%domain) - call mpp_update_domains(oro_g, Atm%domain) - if (Atm%neststruct%nested) then - call extrapolation_BC(oro_g, 0, 0, Atm%npx, Atm%npy, Atm%bd, .true.) + if ( Atm%flagstruct%full_zs_filter) then + allocate (oro_g(isd:ied,jsd:jed)) + oro_g = 0. + ! land-frac + call register_restart_field(ORO_restart, 'land_frac', oro_g, dim_names_2d) + call mpp_update_domains(oro_g, Atm%domain) + if (Atm%neststruct%nested) then + call extrapolation_BC(oro_g, 0, 0, Atm%npx, Atm%npy, Atm%bd, .true.) + endif endif - endif - if ( Atm%flagstruct%fv_land ) then - ! stddev - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'stddev', Atm%sgh, domain=Atm%domain) - ! land-frac - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'land_frac', Atm%oro, domain=Atm%domain) + if ( Atm%flagstruct%fv_land ) then + ! stddev + call register_restart_field(ORO_restart, 'stddev', Atm%sgh, dim_names_2d) + ! land-frac + call register_restart_field(ORO_restart, 'land_frac', Atm%oro, dim_names_2d) + endif + call read_restart(ORO_restart) + call close_file(ORO_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_oro_ics)//' for NGGPS IC does not exist') endif - - - ! read in the restart - call restore_state (ORO_restart) - call restore_state (SFC_restart) - ! free the restart type to be re-used by the nest - call free_restart_type(ORO_restart) - call free_restart_type(SFC_restart) - + call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_oro_ics)//' for NGGPS IC') ! initialize all tracers to default values prior to being input do nt = 1, ntprog @@ -747,16 +733,16 @@ subroutine read_gfs_data_original() allocate (ak(levp+1)) allocate (bk(levp+1)) - call read_data('INPUT/'//trim(fn_gfs_ctl),'vcoord',wk2, no_domain=.TRUE.) - ak(1:levp+1) = wk2(1:levp+1,1) - bk(1:levp+1) = wk2(1:levp+1,2) - deallocate (wk2) - - - if (.not. file_exist('INPUT/'//trim(fn_gfs_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_gfs_ics)//' for NGGPS IC does not exist') + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Gfs_ctl, fn_gfs_ctl, "read", pelist=pes) ) then + call read_data(Gfs_ctl,'vcoord',wk2) + ak(1:levp+1) = wk2(1:levp+1,1) + bk(1:levp+1) = wk2(1:levp+1,2) + deallocate (wk2) + call close_file(Gfs_ctl) endif - call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_gfs_ics)//' for NGGPS IC') + deallocate(pes) allocate (zh(is:ie,js:je,levp+1)) ! SJL allocate (ps(is:ie,js:je)) @@ -768,39 +754,58 @@ subroutine read_gfs_data_original() allocate ( v_s(is:ie, js:je+1, 1:levp) ) if (trim(source) == source_fv3gfs) allocate (temp(is:ie,js:je,1:levp)) + ! initialize dim_names for register restart + dim_names_3d(1) = "lev" + dim_names_3d(2) = "lat" + dim_names_3d(3) = "lonp" + dim_names_3d2 = dim_names_3d + dim_names_3d2(2) = "latp" + dim_names_3d2(3) = "lon" + dim_names_3d3 = dim_names_3d + dim_names_3d3(3) = "lon" + dim_names_3d4 = dim_names_3d3 + dim_names_3d4(1) = "levp" ! surface pressure (Pa) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ps', ps, domain=Atm%domain) - - ! D-grid west face tangential wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'u_w', u_w, domain=Atm%domain,position=EAST) - ! D-grid west face normal wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'v_w', v_w, domain=Atm%domain,position=EAST) - ! D-grid south face tangential wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'u_s', u_s, domain=Atm%domain,position=NORTH) - ! D-grid south face normal wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'v_s', v_s, domain=Atm%domain,position=NORTH) - - ! vertical velocity 'omega' (Pa/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'w', omga, domain=Atm%domain) - ! GFS grid height at edges (including surface height) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ZH', zh, domain=Atm%domain) - - ! real temperature (K) - if (trim(source) == source_fv3gfs) id_res = register_restart_field (GFS_restart, fn_gfs_ics, 't', temp, mandatory=.false., & - domain=Atm%domain) - ! prognostic tracers - do nt = 1, ntracers - q(:,:,:,nt) = -999.99 - call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, trim(tracer_name), q(:,:,:,nt), & - mandatory=.false.,domain=Atm%domain) - enddo - - ! read in the gfs_data and free the restart type to be re-used by the nest - call restore_state(GFS_restart) - call free_restart_type(GFS_restart) + if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(GFS_restart, "lat", "y") + call register_axis(GFS_restart, "lon", "x") + call register_axis(GFS_restart, "lonp", "x", domain_position=east) + call register_axis(GFS_restart, "latp", "y", domain_position=north) + call register_axis(GFS_restart, "lev", size(u_w,3)) + call register_axis(GFS_restart, 'levp', size(zh,3)) + + call register_restart_field(GFS_restart, 'ps', ps, dim_names_2d) + ! D-grid west face tangential wind component (m/s) + call register_restart_field(GFS_restart, 'u_w', u_w, dim_names_3d) + ! D-grid west face normal wind component (m/s) + call register_restart_field(GFS_restart, 'v_w', v_w, dim_names_3d) + ! D-grid south face tangential wind component (m/s) + call register_restart_field(GFS_restart, 'u_s', u_s, dim_names_3d2) + ! D-grid south face normal wind component (m/s) + call register_restart_field(GFS_restart, 'v_s', v_s, dim_names_3d2) + ! vertical velocity 'omega' (Pa/s) + call register_restart_field(GFS_restart, 'w', omga, dim_names_3d3) + ! GFS grid height at edges (including surface height) + call register_restart_field(GFS_restart, 'zh', zh, dim_names_3d4) + + ! real temperature (K) + if (trim(source) == source_fv3gfs) call register_restart_field(GFS_restart, 't', temp, dim_names_3d3, is_optional=.true.) + + ! prognostic tracers + do nt = 1, ntracers + q(:,:,:,nt) = -999.99 + call get_tracer_names(MODEL_ATMOS, nt, tracer_name) + call register_restart_field(GFS_restart, trim(tracer_name), q(:,:,:,nt), dim_names_3d3, is_optional=.true.) + enddo + ! read in the gfs_data and free the restart type to be re-used by the nest + call read_restart(GFS_restart) + call close_file(GFS_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_gfs_ics)//' for NGGPS IC does not exist') + endif + call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_gfs_ics)//' for NGGPS IC') endsubroutine read_gfs_data_original @@ -840,50 +845,60 @@ subroutine read_gfs_data_bad() !--- read in ak and bk from the gfs control file using fms_io read_data --- ! ! put the lowest 64 levels into ak and bk - call read_data('INPUT/'//trim(fn_gfs_ctl),'vcoord',wk2_tmp, no_domain=.TRUE.) - ak(1:levp+1) = wk2_tmp(2:levsp,1) - bk(1:levp+1) = wk2_tmp(2:levsp,2) - - deallocate (wk2_tmp) - - - if (.not. file_exist('INPUT/'//trim(fn_gfs_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_gfs_ics)//' for NGGPS IC does not exist') + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Gfs_ctl, fn_gfs_ctl, "read", pelist=pes) ) then + call read_data(Gfs_ctl,'vcoord',wk2_tmp) + ak(1:levp+1) = wk2_tmp(2:levsp,1) + bk(1:levp+1) = wk2_tmp(2:levsp,2) + + deallocate (wk2_tmp) + call close_file(Gfs_ctl) endif - call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_gfs_ics)//' for NGGPS IC') + deallocate(pes) ! surface pressure (Pa) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ps', ps, domain=Atm%domain) - - - ! D-grid west face tangential wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'u_w', u_w_tmp, domain=Atm%domain,position=EAST) - ! D-grid west face normal wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'v_w', v_w_tmp, domain=Atm%domain,position=EAST) - ! D-grid south face tangential wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'u_s', u_s_tmp, domain=Atm%domain,position=NORTH) - ! D-grid south face normal wind component (m/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'v_s', v_s_tmp, domain=Atm%domain,position=NORTH) - ! vertical velocity 'omega' (Pa/s) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'w', omga_tmp, domain=Atm%domain) - ! GFS grid height at edges (including surface height) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ZH', zh_tmp, domain=Atm%domain) - ! real temperature (K) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 't', temp_tmp, mandatory=.false., & - domain=Atm%domain) - - ! prognostic tracers - do nt = 1, ntracers - call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, trim(tracer_name), q_tmp(:,:,:,nt), & - mandatory=.false.,domain=Atm%domain) - enddo + if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(GFS_restart, "lat", "y") + call register_axis(GFS_restart, "lon", "x") + call register_axis(GFS_restart, "lonp", "x", domain_position=east) + call register_axis(GFS_restart, "latp", "y", domain_position=north) + call register_axis(GFS_restart, "lev", size(u_w_tmp,3)) + call register_axis(GFS_restart, "levp", size(zh_tmp,3)) + + + + call register_restart_field(GFS_restart, 'ps', ps, dim_names_2d) + ! D-grid west face tangential wind component (m/s) + call register_restart_field(GFS_restart, 'u_w', u_w_tmp, dim_names_3d) + ! D-grid west face normal wind component (m/s) + call register_restart_field(GFS_restart, 'v_w', v_w_tmp, dim_names_3d) + ! D-grid south face tangential wind component (m/s) + call register_restart_field(GFS_restart, 'u_s', u_s_tmp, dim_names_3d2) + ! D-grid south face normal wind component (m/s) + call register_restart_field(GFS_restart, 'v_s', v_s_tmp, dim_names_3d2) + ! vertical velocity 'omega' (Pa/s) + call register_restart_field(GFS_restart, 'w', omga_tmp, dim_names_3d3) + ! GFS grid height at edges (including surface height) + call register_restart_field(GFS_restart, 'zh', zh_tmp, dim_names_3d4) + ! real temperature (K) + call register_restart_field(GFS_restart, 't', temp_tmp, dim_names_3d3, is_optional=.true.) + + ! Prognostic tracers + do nt = 1, ntracers + call get_tracer_names(MODEL_ATMOS, nt, tracer_name) + call register_restart_field(GFS_restart, trim(tracer_name), q_tmp(:,:,:,nt), dim_names_3d3, is_optional=.true.) + enddo - ! read in the gfs_data and free the restart type to be re-used by the nest - call restore_state(GFS_restart) - call free_restart_type(GFS_restart) + ! read in the gfs_data and free the restart type to be re-used by the nest + call read_restart(GFS_restart) + call close_file(GFS_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_gfs_ics)//' for NGGPS IC does not exist') + endif + call mpp_error(NOTE,'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_gfs_ics)//' for NGGPS IC') ! extract and return the lowest 64 levels of data do nt = 1, ntracers @@ -947,19 +962,25 @@ subroutine get_hrrr_ic (Atm, fv_domain) integer :: is, ie, js, je integer :: isd, ied, jsd, jed integer :: ios, ierr, unit, id_res - type (restart_file_type) :: ORO_restart, SFC_restart, HRRR_restart - character(len=6) :: gn, stile_name + type (FmsNetcdfDomainFile_t) :: ORO_restart, SFC_restart, HRRR_restart + type(FmsNetcdfFile_t) :: Hrr_ctl + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist + character(len=8), dimension(2) :: dim_names_2d + character(len=8), dimension(3) :: dim_names_3d, dim_names_3d2, dim_names_3d3, dim_names_3d4 + character(len=6) :: stile_name character(len=64) :: tracer_name - character(len=64) :: fn_hrr_ctl = 'hrrr_ctrl.nc' - character(len=64) :: fn_hrr_ics = 'hrrr_data.nc' - character(len=64) :: fn_sfc_ics = 'sfc_data.nc' - character(len=64) :: fn_oro_ics = 'oro_data.nc' + character(len=64) :: fn_hrr_ctl = 'INPUT/hrrr_ctrl.nc' + character(len=64) :: fn_hrr_ics = 'INPUT/hrrr_data.nc' + character(len=64) :: fn_sfc_ics = 'INPUT/sfc_data.nc' + character(len=64) :: fn_oro_ics = 'INPUT/oro_data.nc' logical :: remap logical :: filtered_terrain = .true. logical :: gfs_dwinds = .true. integer :: levp = 50 logical :: checker_tr = .false. integer :: nt_checker = 0 + character(len=20) :: suffix + character(len=1) :: tile_num real(kind=R_GRID), dimension(2):: p1, p2, p3 real(kind=R_GRID), dimension(3):: e1, e2, ex, ey integer:: i,j,k,nts, ks @@ -969,15 +990,8 @@ subroutine get_hrrr_ic (Atm, fv_domain) call mpp_error(NOTE,'Using external_IC::get_hrrr_ic which is valid only for data which has been & &horizontally interpolated to the current lambert grid') -#ifdef INTERNAL_FILE_NML read (input_nml_file,external_ic_nml,iostat=ios) ierr = check_nml_error(ios,'external_ic_nml') -#else - unit=open_namelist_file() - read (unit,external_ic_nml,iostat=ios) - ierr = check_nml_error(ios,'external_ic_nml') - call close_file(unit) -#endif unit = stdlog() call write_version_number ( 'EXTERNAL_IC_MOD::get_hrrr_ic', version ) @@ -1006,48 +1020,28 @@ subroutine get_hrrr_ic (Atm, fv_domain) call get_number_tracers(MODEL_ATMOS, num_tracers=ntracers, num_prog=ntprog) ntdiag = ntracers-ntprog -!--- set the 'nestXX' appendix for all files using fms_io - if (Atm%grid_number > 1) then - write(gn,'(A4, I2.2)') "nest", Atm%grid_number - else - gn = '' - end if - call set_filename_appendix('') - -!--- test for existence of the HRRR control file - if (.not. file_exist('INPUT/'//trim(fn_hrr_ctl), no_domain=.TRUE.)) then - call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: file '//trim(fn_hrr_ctl)//' for HRRR IC does not exist') - endif - call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using control file '//trim(fn_hrr_ctl)//' for HRRR IC') - + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Hrr_ctl, fn_hrr_ctl, "read", pelist=pes) ) then !--- read in the number of tracers in the HRRR ICs - call read_data ('INPUT/'//trim(fn_hrr_ctl), 'ntrac', ntrac, no_domain=.TRUE.) - if (ntrac > ntracers) call mpp_error(FATAL,'==> External_ic::get_hrrr_ic: more HRRR tracers & - &than defined in field_table '//trim(fn_hrr_ctl)//' for HRRR IC') + call read_data (Hrr_ctl, 'ntrac', ntrac) + if (ntrac > ntracers) call mpp_error(FATAL,'==> External_ic::get_hrrr_ic: more HRRR tracers & + &than defined in field_table '//trim(fn_hrr_ctl)//' for HRRR IC') !--- read in ak and bk from the HRRR control file using fms_io read_data --- - allocate (wk2(levp+1,2)) - allocate (ak(levp+1)) - allocate (bk(levp+1)) - call read_data('INPUT/'//trim(fn_hrr_ctl),'vcoord',wk2, no_domain=.TRUE.) - ak(1:levp+1) = wk2(1:levp+1,1) - bk(1:levp+1) = wk2(1:levp+1,2) - deallocate (wk2) - - if (.not. file_exist('INPUT/'//trim(fn_oro_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_oro_ics)//' for HRRR IC does not exist') - endif - call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_oro_ics)//' for HRRR IC') - - if (.not. file_exist('INPUT/'//trim(fn_sfc_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_sfc_ics)//' for HRRR IC does not exist') - endif - call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_sfc_ics)//' for HRRR IC') - - if (.not. file_exist('INPUT/'//trim(fn_hrr_ics), domain=Atm%domain)) then - call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_hrr_ics)//' for HRRR IC does not exist') + allocate (wk2(levp+1,2)) + allocate (ak(levp+1)) + allocate (bk(levp+1)) + call read_data(Hrr_ctl,'vcoord',wk2) + ak(1:levp+1) = wk2(1:levp+1,1) + bk(1:levp+1) = wk2(1:levp+1,2) + deallocate (wk2) + call close_file(Hrr_ctl) + else + call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: file '//trim(fn_hrr_ctl)//' for HRRR IC does not exist') endif - call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_hrr_ics)//' for HRRR IC') + deallocate(pes) + call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using control file '//trim(fn_hrr_ctl)//' for HRRR IC') allocate (zh(is:ie,js:je,levp+1)) allocate (ps(is:ie,js:je)) @@ -1069,63 +1063,107 @@ subroutine get_hrrr_ic (Atm, fv_domain) enddo enddo endif + ! set dimensions for register restart + dim_names_2d(1) = "lat" + dim_names_2d(2) = "lon" !--- read in surface temperature (k) and land-frac ! surface skin temperature - id_res = register_restart_field (SFC_restart, fn_sfc_ics, 'tsea', Atm%ts, domain=Atm%domain) + if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(SFC_restart, "lat", "y") + call register_axis(SFC_restart, "lon", "x") + call register_restart_field(SFC_restart, 'tsea', Atm%ts, dim_names_2d) + call read_restart(SFC_restart) + call close_file(SFC_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_sfc_ics)//' for HRRR IC does not exist') + endif + call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_sfc_ics)//' for HRRR IC') ! terrain surface height -- (needs to be transformed into phis = zs*grav) - if (filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_filt', Atm%phis, domain=Atm%domain) - elseif (.not. filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_raw', Atm%phis, domain=Atm%domain) - endif + if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(ORO_restart, "lat", "y") + call register_axis(ORO_restart, "lon", "x") + if (filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_filt', Atm%phis, dim_names_2d) + elseif (.not. filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_raw', Atm%phis, dim_names_2d) + endif - if ( Atm%flagstruct%full_zs_filter) then - allocate (oro_g(isd:ied,jsd:jed)) - oro_g = 0. - ! land-frac - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'land_frac', oro_g, domain=Atm%domain) - call mpp_update_domains(oro_g, Atm%domain) - if (Atm%neststruct%nested) then - call extrapolation_BC(oro_g, 0, 0, Atm%npx, Atm%npy, Atm%bd, .true.) + if ( Atm%flagstruct%full_zs_filter) then + allocate (oro_g(isd:ied,jsd:jed)) + oro_g = 0. + ! land-frac + call register_restart_field(ORO_restart, 'land_frac', oro_g, dim_names_2d) + call mpp_update_domains(oro_g, Atm%domain) + if (Atm%neststruct%nested) then + call extrapolation_BC(oro_g, 0, 0, Atm%npx, Atm%npy, Atm%bd, .true.) + endif endif - endif - if ( Atm%flagstruct%fv_land ) then - ! stddev - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'stddev', Atm%sgh, domain=Atm%domain) - ! land-frac - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'land_frac', Atm%oro, domain=Atm%domain) + if ( Atm%flagstruct%fv_land ) then + ! stddev + call register_restart_field(ORO_restart, 'stddev', Atm%sgh, dim_names_2d) + ! land-frac + call register_restart_field(ORO_restart, 'land_frac', Atm%oro, (/"lat", "lon"/)) + endif + call read_restart(ORO_restart) + call close_file(ORO_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_oro_ics)//' for HRRR IC does not exist') endif + call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_oro_ics)//' for HRRR IC') + + ! initialize dim_names for register restart + dim_names_3d(1) = "lev" + dim_names_3d(2) = "lat" + dim_names_3d(3) = "lonp" + dim_names_3d2 = dim_names_3d + dim_names_3d2(2) = "latp" + dim_names_3d2(3) = "lon" + dim_names_3d3 = dim_names_3d + dim_names_3d3(3) = "lon" + dim_names_3d4 = dim_names_3d3 + dim_names_3d4(1) = "levp" ! edge pressure (Pa) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'ps', ps, domain=Atm%domain) - - ! physical temperature (K) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'pt', t, domain=Atm%domain) - - ! D-grid west face tangential wind component (m/s) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'u_w', u_w, domain=Atm%domain,position=EAST) - ! D-grid west face normal wind component (m/s) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'v_w', v_w, domain=Atm%domain,position=EAST) - ! D-grid south face tangential wind component (m/s) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'u_s', u_s, domain=Atm%domain,position=NORTH) - ! D-grid south face normal wind component (m/s) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'v_s', v_s, domain=Atm%domain,position=NORTH) - + if( open_file(HRRR_restart, fn_hrr_ics, "read", Atm%domain,is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(HRRR_restart, "lat", "y") + call register_axis(HRRR_restart, "lon", "x") + call register_axis(HRRR_restart, "lonp", "x", domain_position=east) + call register_axis(HRRR_restart, "latp", "y", domain_position=north) + call register_axis(HRRR_restart, "lev", size(u_w,3)) + call register_axis(HRRR_restart, "levp", size(zh,3)) + + call register_restart_field(HRRR_restart, 'ps', ps, dim_names_2d) + ! physical temperature (K) + call register_restart_field(HRRR_restart, 'pt', t, dim_names_2d) + ! D-grid west face tangential wind component (m/s) + call register_restart_field(HRRR_restart, 'u_w', u_w, dim_names_3d) + ! D-grid west face normal wind component (m/s) + call register_restart_field(HRRR_restart, 'v_w', v_w, dim_names_3d) + ! D-grid south face tangential wind component (m/s) + call register_restart_field(HRRR_restart, 'u_s', u_s, dim_names_3d2) + ! D-grid south face normal wind component (m/s) + call register_restart_field(HRRR_restart, 'v_s', v_s, dim_names_3d2) + ! vertical velocity 'omega' (m/s) + call register_restart_field(HRRR_restart, 'w', w, dim_names_3d3) + ! GFS grid height at edges (including surface height) + call register_restart_field(HRRR_restart, 'zh', zh, dim_names_3d4) + + ! prognostic tracers + do nt = 1, ntracers + call get_tracer_names(MODEL_ATMOS, nt, tracer_name) + call register_restart_field(HRRR_restart, trim(tracer_name), q(:,:,:,nt), dim_names_3d3, is_optional=.true.) + enddo + call read_restart(HRRR_restart) + call close_file(HRRR_restart) + else + call mpp_error(FATAL,'==> Error in External_ic::get_hrrr_ic: tiled file '//trim(fn_hrr_ics)//' for HRRR IC does not exist') + endif + call mpp_error(NOTE,'==> External_ic::get_hrrr_ic: using tiled data file '//trim(fn_hrr_ics)//' for HRRR IC') - ! vertical velocity (m/s) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'w', w, domain=Atm%domain) - ! GFS grid height at edges (including surface height) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, 'ZH', zh, domain=Atm%domain) - ! prognostic tracers - do nt = 1, ntracers - call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - id_res = register_restart_field (HRRR_restart, fn_hrr_ics, trim(tracer_name), q(:,:,:,nt), & - mandatory=.false.,domain=Atm%domain) - enddo ! initialize all tracers to default values prior to being input do nt = 1, ntprog @@ -1140,16 +1178,6 @@ subroutine get_hrrr_ic (Atm, fv_domain) enddo - ! read in the restart - call restore_state (ORO_restart) - call restore_state (SFC_restart) - call restore_state (HRRR_restart) - ! free the restart type to be re-used by the nest - call free_restart_type(ORO_restart) - call free_restart_type(SFC_restart) - call free_restart_type(HRRR_restart) - - ! multiply static terrain 'phis' by gravity to be true geopotential Atm%phis = Atm%phis*grav @@ -1396,7 +1424,7 @@ subroutine get_ncep_ic( Atm, fv_domain, nq ) fname = Atm%flagstruct%res_latlon_dynamics - if( file_exist(fname) ) then + if( file_exists(fname) ) then call open_ncfile( fname, ncid ) ! open the file call get_ncdim1( ncid, 'lon', tsize(1) ) call get_ncdim1( ncid, 'lat', tsize(2) ) @@ -1872,6 +1900,8 @@ subroutine get_ecmwf_ic( Atm, fv_domain ) 0.97771, 0.98608, 0.99347, 1./ character(len=128) :: fname + character(len=8), dimension(2) :: dim_names_2d + character(len=8), dimension(3) :: dim_names_3d3, dim_names_3d4 real, allocatable:: wk2(:,:) real(kind=4), allocatable:: wk2_r4(:,:) real, dimension(:,:,:), allocatable:: ud, vd @@ -1910,10 +1940,14 @@ subroutine get_ecmwf_ic( Atm, fv_domain ) integer :: id_res, ntprog, ntracers, ks, iq, nt character(len=64) :: tracer_name integer :: levp_gfs = 64 - type (restart_file_type) :: ORO_restart, GFS_restart - character(len=64) :: fn_oro_ics = 'oro_data.nc' - character(len=64) :: fn_gfs_ics = 'gfs_data.nc' - character(len=64) :: fn_gfs_ctl = 'gfs_ctrl.nc' + type(FmsNetcdfDomainFile_t) :: ORO_restart, GFS_restart + type(FmsNetcdfFile_t) :: Gfs_ctl + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist + character(len=64) :: fn_oro_ics = 'INPUT/oro_data.nc' + character(len=64) :: fn_gfs_ics = 'INPUT/gfs_data.nc' + character(len=64) :: fn_gfs_ctl = 'INPUT/gfs_ctrl.nc' + character(len=20) :: suffix + character(len=1) :: tile_num logical :: filtered_terrain = .true. namelist /external_ic_nml/ filtered_terrain @@ -1966,14 +2000,35 @@ subroutine get_ecmwf_ic( Atm, fv_domain ) !!$ call set_eta(npz, ks, Atm%ptop, Atm%ak, Atm%bk, Atm%flagstruct%npz_type) !!$ endif + !!! If a nested grid, add "nestXX.tileX" to the filename + if (Atm%neststruct%nested) then + write(tile_num,'(I1)') Atm%neststruct%nestupdate + suffix = ".tile"//trim(tile_num)//"" + fn_oro_ics = fn_oro_ics(1:len_trim(fn_oro_ics)-2)//trim(suffix)//".nc" + fn_gfs_ics = fn_gfs_ics(1:len_trim(fn_gfs_ics)-2)//trim(suffix)//".nc" + endif + +! initialize dim_names for register_restart + dim_names_2d(1) = "lat" + dim_names_2d(2) = "lon" + dim_names_3d3(1) = "lev" + dim_names_3d3(2) = "lat" + dim_names_3d3(3) = "lon" + dim_names_3d4 = dim_names_3d3 + dim_names_3d4(1) = "levp" + !! Read in model terrain from oro_data.tile?.nc - if (filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_filt', Atm%phis, domain=Atm%domain) - elseif (.not. filtered_terrain) then - id_res = register_restart_field (ORO_restart, fn_oro_ics, 'orog_raw', Atm%phis, domain=Atm%domain) + if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(ORO_restart, "lat", "y") + call register_axis(ORO_restart, "lon", "x") + if (filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_filt', Atm%phis, dim_names_2d) + elseif (.not. filtered_terrain) then + call register_restart_field(ORO_restart, 'orog_raw', Atm%phis, dim_names_2d) + endif + call read_restart(ORO_restart) + call close_file(ORO_restart) endif - call restore_state (ORO_restart) - call free_restart_type(ORO_restart) Atm%phis = Atm%phis*grav if(is_master()) write(*,*) 'done reading model terrain from oro_data.nc' call mpp_update_domains( Atm%phis, Atm%domain ) @@ -1983,19 +2038,29 @@ subroutine get_ecmwf_ic( Atm, fv_domain ) allocate (ps_gfs(is:ie,js:je)) allocate (zh_gfs(is:ie,js:je,levp_gfs+1)) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'o3mr', o3mr_gfs, & - mandatory=.false.,domain=Atm%domain) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ps', ps_gfs, domain=Atm%domain) - id_res = register_restart_field (GFS_restart, fn_gfs_ics, 'ZH', zh_gfs, domain=Atm%domain) - call restore_state (GFS_restart) - call free_restart_type(GFS_restart) - + if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then + call register_axis(GFS_restart, "lat", "y") + call register_axis(GFS_restart, "lon", "x") + call register_axis(GFS_restart, "lev", size(o3mr_gfs,3)) + call register_axis(GFS_restart, "levp", size(zh_gfs,3)) + call register_restart_field(GFS_restart, 'o3mr', o3mr_gfs, dim_names_3d3, is_optional=.true.) + call register_restart_field(GFS_restart, 'ps', ps_gfs, dim_names_2d) + call register_restart_field(GFS_restart, 'ZH', zh_gfs, dim_names_3d4) + call read_restart(GFS_restart) + call close_file(GFS_restart) + endif ! Get GFS ak, bk for o3mr vertical interpolation allocate (wk2(levp_gfs+1,2)) allocate (ak_gfs(levp_gfs+1)) allocate (bk_gfs(levp_gfs+1)) - call read_data('INPUT/'//trim(fn_gfs_ctl),'vcoord',wk2, no_domain=.TRUE.) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Gfs_ctl, fn_gfs_ctl, "read", pelist=pes) ) then + call read_data(Gfs_ctl,'vcoord',wk2) + call close_file(Gfs_ctl) + endif + deallocate(pes) ak_gfs(1:levp_gfs+1) = wk2(1:levp_gfs+1,1) bk_gfs(1:levp_gfs+1) = wk2(1:levp_gfs+1,2) deallocate (wk2) @@ -2014,7 +2079,7 @@ subroutine get_ecmwf_ic( Atm, fv_domain ) !! Start to read EC data fname = Atm%flagstruct%res_latlon_dynamics - if( file_exist(fname) ) then + if( file_exists(fname) ) then call open_ncfile( fname, ncid ) ! open the file call get_ncdim1( ncid, 'longitude', tsize(1) ) @@ -2457,6 +2522,8 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) type(domain2d), intent(inout) :: fv_domain integer, intent(in):: nq + type(FmsNetcdfFile_t) :: Latlon_dyn, Latlon_tra + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist character(len=128) :: fname, tracer_name real, allocatable:: ps0(:,:), gz0(:,:), u0(:,:,:), v0(:,:,:), t0(:,:,:), dp0(:,:,:), q0(:,:,:,:) real, allocatable:: ua(:,:,:), va(:,:,:) @@ -2473,9 +2540,10 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) ! Read in lat-lon FV core restart file fname = Atm%flagstruct%res_latlon_dynamics - - if( file_exist(fname) ) then - call field_size(fname, 'T', tsize, field_found=found) + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if( open_file(Latlon_dyn, fname, "read", is_restart=.true., pelist=pes) ) then + call get_variable_size(Latlon_dyn, 'T', tsize) if(is_master()) write(*,*) 'Using lat-lon FV restart:', fname if ( found ) then @@ -2506,13 +2574,13 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) allocate ( t0(1:im,1:jm,1:km) ) allocate ( dp0(1:im,1:jm,1:km) ) - call read_data (fname, 'ak', ak0) - call read_data (fname, 'bk', bk0) - call read_data (fname, 'Surface_geopotential', gz0) - call read_data (fname, 'U', u0) - call read_data (fname, 'V', v0) - call read_data (fname, 'T', t0) - call read_data (fname, 'DELP', dp0) + call read_data (Latlon_dyn, 'ak', ak0) + call read_data (Latlon_dyn, 'bk', bk0) + call read_data (Latlon_dyn, 'Surface_geopotential', gz0) + call read_data (Latlon_dyn, 'U', u0) + call read_data (Latlon_dyn, 'V', v0) + call read_data (Latlon_dyn, 'T', t0) + call read_data (Latlon_dyn, 'DELP', dp0) ! Share the load if(is_master()) call pmaxmin( 'ZS_data', gz0, im, jm, 1./grav) @@ -2520,7 +2588,7 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) if(mpp_pe()==1) call pmaxmin( 'V_data', v0, im*jm, km, 1.) if(mpp_pe()==2) call pmaxmin( 'T_data', t0, im*jm, km, 1.) if(mpp_pe()==3) call pmaxmin( 'DEL-P', dp0, im*jm, km, 0.01) - + call close_file(Latlon_dyn) else call mpp_error(FATAL,'==> Error in get_external_ic: Expected file '//trim(fname)//' for dynamics does not exist') @@ -2528,8 +2596,7 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) ! Read in tracers: only AM2 "physics tracers" at this point fname = Atm%flagstruct%res_latlon_tracers - - if( file_exist(fname) ) then + if( open_file(Latlon_tra, fname, "read", is_restart=.true., pelist=pes) ) then if(is_master()) write(*,*) 'Using lat-lon tracer restart:', fname allocate ( q0(im,jm,km,Atm%ncnst) ) @@ -2537,15 +2604,17 @@ subroutine get_fv_ic( Atm, fv_domain, nq ) do tr_ind = 1, nq call get_tracer_names(MODEL_ATMOS, tr_ind, tracer_name) - if (field_exist(fname,tracer_name)) then - call read_data(fname, tracer_name, q0(1:im,1:jm,1:km,tr_ind)) + if (variable_exists(Latlon_tra,tracer_name)) then + call read_data(Latlon_tra, tracer_name, q0(1:im,1:jm,1:km,tr_ind)) call mpp_error(NOTE,'==> Have read tracer '//trim(tracer_name)//' from '//trim(fname)) cycle endif enddo + call close_file(Latlon_tra) else call mpp_error(FATAL,'==> Error in get_external_ic: Expected file '//trim(fname)//' for tracers does not exist') endif + deallocate(pes) ! D to A transform on lat-lon grid: allocate ( ua(im,jm,km) ) @@ -4472,14 +4541,27 @@ subroutine get_data_source(source,regional) character (len = 80) :: source integer :: ncids,sourceLength logical :: lstatus,regional + type(FmsNetcdfFile_t) :: Gfs_data + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist ! ! Use the fms call here so we can actually get the return code value. ! + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) if (regional) then - lstatus = get_global_att_value('INPUT/gfs_data.nc',"source", source) + if (open_file(Gfs_data , 'INPUT/gfs_data.nc', "read", pelist=pes)) then + lstatus = global_att_exists(Gfs_data, "source") + if(lstatus) call get_global_attribute(Gfs_data, "source", source) + call close_file(Gfs_data) + endif else - lstatus = get_global_att_value('INPUT/gfs_data.tile1.nc',"source", source) + if (open_file(Gfs_data , 'INPUT/gfs_data.tile1.nc', "read", pelist=pes)) then + lstatus = global_att_exists(Gfs_data, "source") + if(lstatus) call get_global_attribute(Gfs_data, "source", source) + call close_file(Gfs_data) + endif endif + deallocate(pes) if (.not. lstatus) then if (mpp_pe() == 0) write(0,*) 'INPUT source not found ',lstatus,' set source=No Source Attribute' source='No Source Attribute' diff --git a/tools/fv_diagnostics.F90 b/tools/fv_diagnostics.F90 index 88bba0f29..bcf064352 100644 --- a/tools/fv_diagnostics.F90 +++ b/tools/fv_diagnostics.F90 @@ -27,7 +27,6 @@ module fv_diagnostics_mod use constants_mod, only: grav, rdgas, rvgas, pi=>pi_8, radius, kappa, WTMAIR, WTMCO2, & omega, hlv, cp_air, cp_vapor, TFREEZE use fms_mod, only: write_version_number - use fms_io_mod, only: set_domain, nullify_domain use time_manager_mod, only: time_type, get_date, get_time use mpp_domains_mod, only: domain2d, mpp_update_domains, DGRID_NE, NORTH, EAST use diag_manager_mod, only: diag_axis_init, register_diag_field, & @@ -45,7 +44,6 @@ module fv_diagnostics_mod use tracer_manager_mod, only: get_tracer_names, get_number_tracers, get_tracer_index use field_manager_mod, only: MODEL_ATMOS use mpp_mod, only: mpp_error, FATAL, stdlog, mpp_pe, mpp_root_pe, mpp_sum, mpp_max, NOTE, input_nml_file - use mpp_io_mod, only: mpp_flush use sat_vapor_pres_mod, only: compute_qs, lookup_es use fv_arrays_mod, only: max_step @@ -177,8 +175,6 @@ subroutine fv_diag_init(Atm, axes, Time, npx, npy, npz, p_ref) ncnst = Atm(1)%ncnst m_calendar = Atm(1)%flagstruct%moist_phys - call set_domain(Atm(1)%domain) ! Set domain so that diag_manager can access tile information - sphum = get_tracer_index (MODEL_ATMOS, 'sphum') liq_wat = get_tracer_index (MODEL_ATMOS, 'liq_wat') ice_wat = get_tracer_index (MODEL_ATMOS, 'ice_wat') @@ -1311,8 +1307,6 @@ subroutine fv_diag_init(Atm, axes, Time, npx, npy, npz, p_ref) yr_init = 0 ; mo_init = 0 ; hr_init = 0 ; mn_init = 0 endif - call nullify_domain() ! Nullify set_domain info - module_is_initialized=.true. istep = 0 #ifndef GFS_PHYS @@ -1468,7 +1462,6 @@ subroutine fv_diag(Atm, zvir, Time, print_freq) endif fv_time = Time - call set_domain(Atm(1)%domain) if ( m_calendar ) then call get_date(fv_time, yr, mon, dd, hr, mn, seconds) @@ -3780,8 +3773,6 @@ subroutine fv_diag(Atm, zvir, Time, print_freq) if (allocated(dmmr)) deallocate(dmmr) if (allocated(dvmr)) deallocate(dvmr) - call nullify_domain() - end subroutine fv_diag subroutine wind_max(isc, iec, jsc, jec ,isd, ied, jsd, jed, us, vs, ws_max, domain) diff --git a/tools/fv_grid_tools.F90 b/tools/fv_grid_tools.F90 index 36f5ef7ae..d1087b918 100644 --- a/tools/fv_grid_tools.F90 +++ b/tools/fv_grid_tools.F90 @@ -40,7 +40,6 @@ module fv_grid_tools_mod mpp_get_data_domain, mpp_get_compute_domain, & mpp_get_global_domain, mpp_global_sum, mpp_global_max, mpp_global_min use mpp_domains_mod, only: domain2d - use mpp_io_mod, only: mpp_get_att_value use mpp_parameter_mod, only: AGRID_PARAM=>AGRID, & DGRID_NE_PARAM=>DGRID_NE, & @@ -50,9 +49,9 @@ module fv_grid_tools_mod BGRID_SW_PARAM=>BGRID_SW, & SCALAR_PAIR, & CORNER, CENTER, XUPDATE - use fms_mod, only: get_mosaic_tile_grid - use fms_io_mod, only: file_exist, field_exist, read_data, & - get_global_att_value, get_var_att_value + use fms2_io_mod, only: file_exists, variable_exists, open_file, read_data, & + get_global_attribute, get_variable_attribute, & + close_file, get_mosaic_tile_grid, FmsNetcdfFile_t use mosaic_mod, only : get_mosaic_ntiles use mpp_mod, only: mpp_transmit, mpp_recv @@ -86,6 +85,7 @@ subroutine read_grid(Atm, grid_file, ndims, nregions, ng) integer, intent(IN) :: nregions integer, intent(IN) :: ng + type(FmsNetcdfFile_t) :: Grid_input real, allocatable, dimension(:,:) :: tmpx, tmpy real(kind=R_GRID), pointer, dimension(:,:,:) :: grid character(len=128) :: units = "" @@ -108,31 +108,34 @@ subroutine read_grid(Atm, grid_file, ndims, nregions, ng) jed = Atm%bd%jed grid => Atm%gridstruct%grid_64 - if(.not. file_exist(grid_file)) call mpp_error(FATAL, 'fv_grid_tools(read_grid): file '// & + if(.not. file_exists(grid_file)) call mpp_error(FATAL, 'fv_grid_tools(read_grid): file '// & trim(grid_file)//' does not exist') !--- make sure the grid file is mosaic file. - if( field_exist(grid_file, 'atm_mosaic_file') .OR. field_exist(grid_file, 'gridfiles') ) then - stdunit = stdout() - write(stdunit,*) '==>Note from fv_grid_tools_mod(read_grid): read atmosphere grid from mosaic version grid' - else - call mpp_error(FATAL, 'fv_grid_tools(read_grid): neither atm_mosaic_file nor gridfiles exists in file ' & - //trim(grid_file)) - endif + if( open_file(Grid_input, grid_file, "read") ) then + if( variable_exists(Grid_input, 'atm_mosaic_file') .OR. variable_exists(Grid_input, 'gridfiles') ) then + stdunit = stdout() + write(stdunit,*) '==>Note from fv_grid_tools_mod(read_grid): read atmosphere grid from mosaic version grid' + else + call mpp_error(FATAL, 'fv_grid_tools(read_grid): neither atm_mosaic_file nor gridfiles exists in file ' & + //trim(grid_file)) + endif - if(field_exist(grid_file, 'atm_mosaic_file')) then - call read_data(grid_file, "atm_mosaic_file", atm_mosaic) - atm_mosaic = "INPUT/"//trim(atm_mosaic) - else - atm_mosaic = trim(grid_file) + if(variable_exists(Grid_input, 'atm_mosaic_file') ) then + call read_data(Grid_input, "atm_mosaic_file", atm_mosaic) + atm_mosaic = "INPUT/"//trim(atm_mosaic) + else + atm_mosaic = trim(grid_file) + endif + call close_file(Grid_input) endif call get_mosaic_tile_grid(atm_hgrid, atm_mosaic, Atm%domain) grid_form = "none" - if( get_global_att_value(atm_hgrid, "history", attvalue) ) then + if (open_file(Grid_input, atm_hgrid, "read")) then + call get_global_attribute(Grid_input, "history", attvalue) if( index(attvalue, "gnomonic_ed") > 0) grid_form = "gnomonic_ed" - endif if(grid_form .NE. "gnomonic_ed") call mpp_error(FATAL, & "fv_grid_tools(read_grid): the grid should be 'gnomonic_ed' when reading from grid file, contact developer") @@ -144,22 +147,24 @@ subroutine read_grid(Atm, grid_file, ndims, nregions, ng) 'fv_grid_tools(read_grid): nregions should be 6 when reading from mosaic file '//trim(grid_file) ) endif - call get_var_att_value(atm_hgrid, 'x', 'units', units) + call get_variable_attribute(Grid_input, 'x', 'units', units) - !--- get the geographical coordinates of super-grid. - isc2 = 2*is-1; iec2 = 2*ie+1 - jsc2 = 2*js-1; jec2 = 2*je+1 - if( Atm%gridstruct%bounded_domain ) then - isc2 = 2*(isd+halo)-1; iec2 = 2*(ied+1+halo)-1 ! For the regional domain the cell corner locations must be transferred - jsc2 = 2*(jsd+halo)-1; jec2 = 2*(jed+1+halo)-1 ! from the entire supergrid to the compute grid, including the halo region. + !--- get the geographical coordinates of super-grid. + isc2 = 2*is-1; iec2 = 2*ie+1 + jsc2 = 2*js-1; jec2 = 2*je+1 + if( Atm%gridstruct%bounded_domain ) then + isc2 = 2*(isd+halo)-1; iec2 = 2*(ied+1+halo)-1 ! For the regional domain the cell corner locations must be transferred + jsc2 = 2*(jsd+halo)-1; jec2 = 2*(jed+1+halo)-1 ! from the entire supergrid to the compute grid, including the halo region. + endif + allocate(tmpx(isc2:iec2, jsc2:jec2) ) + allocate(tmpy(isc2:iec2, jsc2:jec2) ) + start = 1; nread = 1 + start(1) = isc2; nread(1) = iec2 - isc2 + 1 + start(2) = jsc2; nread(2) = jec2 - jsc2 + 1 + call read_data(Grid_input, 'x', tmpx, corner=start, edge_lengths=nread) !<-- tmpx (lon, deg east) is on the supergrid + call read_data(Grid_input, 'y', tmpy, corner=start, edge_lengths=nread) !<-- tmpy (lat, deg) is on the supergrid + call close_file(Grid_input) endif - allocate(tmpx(isc2:iec2, jsc2:jec2) ) - allocate(tmpy(isc2:iec2, jsc2:jec2) ) - start = 1; nread = 1 - start(1) = isc2; nread(1) = iec2 - isc2 + 1 - start(2) = jsc2; nread(2) = jec2 - jsc2 + 1 - call read_data(atm_hgrid, 'x', tmpx, start, nread, no_domain=.TRUE.) !<-- tmpx (lon, deg east) is on the supergrid - call read_data(atm_hgrid, 'y', tmpy, start, nread, no_domain=.TRUE.) !<-- tmpy (lat, deg) is on the supergrid !--- geographic grid at cell corner grid(isd: is-1, jsd:js-1,1:ndims)=0. @@ -1237,6 +1242,7 @@ subroutine setup_orthogonal_grid(npx, npy, bd, grid_file) ! real, pointer, dimension(:,:,:) :: e1, e2 + type(FmsNetcdfFile_t) :: Grid_input character(len=256) :: atm_mosaic, atm_hgrid real, allocatable, dimension(:,:) :: tmpx, tmpy, tmpu, tmpv, tmpa @@ -1262,23 +1268,26 @@ subroutine setup_orthogonal_grid(npx, npy, bd, grid_file) jed = bd%jed - if(.not. file_exist(grid_file)) call mpp_error(FATAL, 'fv_grid_tools(read_grid): file '// & + if(.not. file_exists(grid_file)) call mpp_error(FATAL, 'fv_grid_tools(read_grid): file '// & trim(grid_file)//' does not exist') !--- make sure the grid file is mosaic file. - if( field_exist(grid_file, 'atm_mosaic_file') .OR. field_exist(grid_file, 'gridfiles') ) then - stdunit = stdout() - write(stdunit,*) '==>Note from fv_grid_tools_mod(read_grid): read atmosphere grid from mosaic version grid' - else - call mpp_error(FATAL, 'fv_grid_tools(read_grid): neither atm_mosaic_file nor gridfiles exists in file ' & - //trim(grid_file)) - endif + if( open_file(Grid_input, grid_file, "read") ) then + if( variable_exists(Grid_input, 'atm_mosaic_file') .OR. variable_exists(Grid_input, 'gridfiles') ) then + stdunit = stdout() + write(stdunit,*) '==>Note from fv_grid_tools_mod(read_grid): read atmosphere grid from mosaic version grid' + else + call mpp_error(FATAL, 'fv_grid_tools(read_grid): neither atm_mosaic_file nor gridfiles exists in file ' & + //trim(grid_file)) + endif - if(field_exist(grid_file, 'atm_mosaic_file')) then - call read_data(grid_file, "atm_mosaic_file", atm_mosaic) - atm_mosaic = "INPUT/"//trim(atm_mosaic) - else - atm_mosaic = trim(grid_file) + if(variable_exists(Grid_input, 'atm_mosaic_file') ) then + call read_data(Grid_input, "atm_mosaic_file", atm_mosaic) + atm_mosaic = "INPUT/"//trim(atm_mosaic) + else + atm_mosaic = trim(grid_file) + endif + call close_file(Grid_input) endif call get_mosaic_tile_grid(atm_hgrid, atm_mosaic, Atm%domain) @@ -1295,69 +1304,73 @@ subroutine setup_orthogonal_grid(npx, npy, bd, grid_file) start = 1; nread = 1 start(1) = isc2; nread(1) = iec2 - isc2 + 1 start(2) = jsc2; nread(2) = jec2 - jsc2 + 1 - call read_data(atm_hgrid, 'x', tmpx, start, nread, no_domain=.TRUE.) !<-- tmpx (lon, deg east) is on the supergrid - call read_data(atm_hgrid, 'y', tmpy, start, nread, no_domain=.TRUE.) !<-- tmpy (lat, deg) is on the supergrid + if (open_file(Grid_input, atm_hgrid, "read")) then + call read_data(Grid_input, 'x', tmpx, corner=start, edge_lengths=nread) !<-- tmpx (lon, deg east) is on the supergrid + call read_data(Grid_input, 'y', tmpy, corner=start, edge_lengths=nread) !<-- tmpy (lat, deg) is on the supergrid - !--- geographic grid at cell corner - grid(isd: is-1, jsd:js-1,1:ndims)=0. - grid(isd: is-1, je+2:jed+1,1:ndims)=0. - grid(ie+2:ied+1,jsd:js-1,1:ndims)=0. - grid(ie+2:ied+1,je+2:jed+1,1:ndims)=0. + !--- geographic grid at cell corner + grid(isd: is-1, jsd:js-1,1:ndims)=0. + grid(isd: is-1, je+2:jed+1,1:ndims)=0. + grid(ie+2:ied+1,jsd:js-1,1:ndims)=0. + grid(ie+2:ied+1,je+2:jed+1,1:ndims)=0. - do j = jsd, jed+1 - do i = isd, ied+1 - grid(i,j,1) = tmpx(2*i+halo+2,2*j+halo+2)*pi/180. - grid(i,j,2) = tmpy(2*i+halo+2,2*j+halo+2)*pi/180. + do j = jsd, jed+1 + do i = isd, ied+1 + grid(i,j,1) = tmpx(2*i+halo+2,2*j+halo+2)*pi/180. + grid(i,j,2) = tmpy(2*i+halo+2,2*j+halo+2)*pi/180. + enddo enddo - enddo - call mpp_update_domains( grid, Atm%domain, position=CORNER) + call mpp_update_domains( grid, Atm%domain, position=CORNER) - iec2 = 2*(ied+1+halo)-2 ! For the regional domain the cell corner locations must be transferred - jec2 = 2*(jed+1+halo)-1 ! from the entire supergrid to the compute grid, including the halo region. + iec2 = 2*(ied+1+halo)-2 ! For the regional domain the cell corner locations must be transferred + jec2 = 2*(jed+1+halo)-1 ! from the entire supergrid to the compute grid, including the halo region. - allocate(tmpu(isc2:iec2, jsc2:jec2) ) + allocate(tmpu(isc2:iec2, jsc2:jec2) ) - nread(1) = iec2 - isc2 + 1 - nread(2) = jec2 - jsc2 + 1 - call read_data(atm_hgrid, 'dx', tmpu, start, nread, no_domain=.TRUE.) + nread(1) = iec2 - isc2 + 1 + nread(2) = jec2 - jsc2 + 1 + call read_data(Grid_input, 'dx', tmpu, corner=start, edge_lengths=nread) - do j = jsd, jed+1 - do i = isd, ied - dx(i,j) = tmpu(2*i+halo+2,2*j+halo+2) + tmpu(2*i+halo+3,2*j+halo+2) + do j = jsd, jed+1 + do i = isd, ied + dx(i,j) = tmpu(2*i+halo+2,2*j+halo+2) + tmpu(2*i+halo+3,2*j+halo+2) + enddo enddo - enddo - iec2 = 2*(ied+1+halo)-1 ! For the regional domain the cell corner locations must be transferred - jec2 = 2*(jed+1+halo)-2 ! from the entire supergrid to the compute grid, including the halo region. + iec2 = 2*(ied+1+halo)-1 ! For the regional domain the cell corner locations must be transferred + jec2 = 2*(jed+1+halo)-2 ! from the entire supergrid to the compute grid, including the halo region. - allocate(tmpv(isc2:iec2, jsc2:jec2) ) + allocate(tmpv(isc2:iec2, jsc2:jec2) ) - nread(1) = iec2 - isc2 + 1 - nread(2) = jec2 - jsc2 + 1 - call read_data(atm_hgrid, 'dy', tmpv, start, nread, no_domain=.TRUE.) + nread(1) = iec2 - isc2 + 1 + nread(2) = jec2 - jsc2 + 1 + call read_data(Grid_input, 'dy', tmpv, corner=start, edge_lengths=nread) - do j = jsd, jed - do i = isd, ied+1 - dy(i,j) = tmpv(2*i+halo+2,2*j+halo+2) + tmpv(2*i+halo+2,2*j+halo+3) + do j = jsd, jed + do i = isd, ied+1 + dy(i,j) = tmpv(2*i+halo+2,2*j+halo+2) + tmpv(2*i+halo+2,2*j+halo+3) + enddo enddo - enddo - call mpp_update_domains( dy, dx, Atm%domain, flags=SCALAR_PAIR, & - gridtype=CGRID_NE_PARAM, complete=.true.) + call mpp_update_domains( dy, dx, Atm%domain, flags=SCALAR_PAIR, & + gridtype=CGRID_NE_PARAM, complete=.true.) + + iec2 = 2*(ied+1+halo)-2 ! For the regional domain the cell corner locations must be transferred + jec2 = 2*(jed+1+halo)-2 ! from the entire supergrid to the compute grid, including the halo region. - iec2 = 2*(ied+1+halo)-2 ! For the regional domain the cell corner locations must be transferred - jec2 = 2*(jed+1+halo)-2 ! from the entire supergrid to the compute grid, including the halo region. + allocate(tmpa(isc2:iec2, jsc2:jec2) ) - allocate(tmpa(isc2:iec2, jsc2:jec2) ) + nread(1) = iec2 - isc2 + 1 + nread(2) = jec2 - jsc2 + 1 + call read_data(Grid_input, 'area', tmpa, corner=start, edge_lengths=nread) !<-- tmpx (lon, deg east) is on the supergrid + call close_file(Grid_input) + endif - nread(1) = iec2 - isc2 + 1 - nread(2) = jec2 - jsc2 + 1 - call read_data(atm_hgrid, 'area', tmpa, start, nread, no_domain=.TRUE.) !<-- tmpx (lon, deg east) is on the supergrid !agrid(:,:,:) = -1.e25 area_c(:,:) = -missing ! To prevent divide by zero error diff --git a/tools/fv_io.F90 b/tools/fv_io.F90 index 0d72dff7a..525abac57 100644 --- a/tools/fv_io.F90 +++ b/tools/fv_io.F90 @@ -33,16 +33,15 @@ module fv_io_mod ! for the model. ! - use fms_mod, only: file_exist - use fms_io_mod, only: fms_io_exit, get_tile_string, & - restart_file_type, register_restart_field, & - save_restart, restore_state, & - set_domain, nullify_domain, set_filename_appendix, & - get_mosaic_tile_file, get_instance_filename, & - save_restart_border, restore_state_border, free_restart_type, & - field_exist + use fms2_io_mod, only: FmsNetcdfFile_t, FmsNetcdfDomainFile_t, & + register_restart_field, register_axis, unlimited, & + open_file, read_restart, read_restart_bc, write_restart, & + write_restart_bc, close_file, register_field, write_data, & + get_global_io_domain_indices, register_variable_attribute, & + variable_exists, read_data, set_filename_appendix use mpp_mod, only: mpp_error, FATAL, NOTE, WARNING, mpp_root_pe, & - mpp_sync, mpp_pe, mpp_declare_pelist + mpp_sync, mpp_pe, mpp_declare_pelist, mpp_get_current_pelist, & + mpp_npes use mpp_domains_mod, only: domain2d, EAST, WEST, NORTH, CENTER, SOUTH, CORNER, & mpp_get_compute_domain, mpp_get_data_domain, & mpp_get_layout, mpp_get_ntile_count, & @@ -57,7 +56,6 @@ module fv_io_mod use fv_eta_mod, only: set_external_eta use fv_mp_mod, only: mp_gather, is_master - use fms_io_mod, only: set_domain use fv_treat_da_inc_mod, only: read_da_inc implicit none @@ -67,6 +65,7 @@ module fv_io_mod public :: fv_io_read_tracers, fv_io_register_restart, fv_io_register_nudge_restart public :: fv_io_register_restart_BCs public :: fv_io_write_BCs, fv_io_read_BCs + public :: fv_io_register_axis logical :: module_is_initialized = .FALSE. @@ -102,6 +101,226 @@ end subroutine fv_io_exit ! NAME="fv_io_exit" + !##################################################################### + ! + ! + ! + ! Register the fv axis for new fms2 io + ! + ! + subroutine fv_io_register_axis(file_obj, numx, xpos, numy, ypos, numz, zsize) + type(FmsNetcdfDomainFile_t), intent(inout) :: file_obj + integer, intent(in), optional :: numx, numy, numz + integer, dimension(:), intent(in), optional :: xpos, ypos, zsize + + integer :: i, ie, is, j + integer, dimension(:), allocatable :: buffer + character(len=1) :: suffix + character(len=7) :: axisname + + if (present(numx)) then + do i=1, numx + write(suffix,'(I1)') i + axisname = 'xaxis_'//suffix + call register_axis(file_obj, axisname, 'X', domain_position=xpos(i)) + if (.not. file_obj%is_readonly) then !if writing file + call register_field(file_obj, axisname, "double", (/axisname/)) + call register_variable_attribute(file_obj,axisname, "axis", "X", str_len=1) + call get_global_io_domain_indices(file_obj, axisname, is, ie, buffer) + call write_data(file_obj, axisname, buffer) + endif + end do + endif + + if (present(numy)) then + do i=1, numy + write(suffix,'(I1)') i + axisname = 'yaxis_'//suffix + call register_axis(file_obj, axisname, 'Y', domain_position=ypos(i)) + if (.not. file_obj%is_readonly) then !if writing file + call register_field(file_obj, axisname, "double", (/axisname/)) + call register_variable_attribute(file_obj,axisname, "axis", "Y", str_len=1) + call get_global_io_domain_indices(file_obj, axisname, is, ie, buffer) + call write_data(file_obj, axisname, buffer) + endif + end do + endif + + if (present(numz)) then + do i= 1, numz + write(suffix,'(I1)') i + axisname = 'zaxis_'//suffix + call register_axis(file_obj, axisname, zsize(i)) + if (.not. file_obj%is_readonly) then !if writing file + call register_field(file_obj, axisname, "double", (/axisname/)) + call register_variable_attribute(file_obj,axisname, "axis", "Z", str_len=1) + if (allocated(buffer)) deallocate(buffer) + allocate(buffer(zsize(i))) + do j = 1, zsize(i) + buffer(j) = j + end do + call write_data(file_obj, axisname, buffer) + deallocate(buffer) + endif + end do + endif + + call register_axis(file_obj, "Time", unlimited) + if (.not. file_obj%is_readonly) then !if writing file + call register_field(file_obj, "Time", "double", (/"Time"/)) + call register_variable_attribute(file_obj, "Time", "cartesian_axis", "T", str_len=1) + call register_variable_attribute(file_obj, "Time", "units", "time level", & + str_len=len("time level")) + call register_variable_attribute(file_obj, "Time", "long_name", "Time", & + str_len=len("Time")) + call write_data(file_obj, "Time", 1) + endif + + end subroutine fv_io_register_axis + ! NAME="fv_io_register_axis" + +!##################################################################### + ! + ! + ! + ! register restart field to be written out to restart file. + ! + subroutine fv_io_register_restart(Atm) + + type(fv_atmos_type), intent(inout) :: Atm + character(len=64) :: tracer_name + character(len=8), dimension(1) :: dim_names + character(len=8), dimension(2) :: dim_names_2d + character(len=8), dimension(4) :: dim_names_4d, dim_names_4d2, dim_names_4d3 + character(len=8), dimension(3) :: dim_names_3d, dim_names_3d2 + integer :: i, j + integer :: nt, ntracers, ntprog, ntdiag + integer, dimension(:), allocatable :: buffer + integer, parameter :: numx=1, numx_2d=2, numy=1, numy_2d=2, numz=1 + integer, dimension(1) :: xpos + integer, dimension(2) :: xpos_2d + integer, dimension(1) :: ypos + integer, dimension(2) :: ypos_2d + integer, dimension(numz) :: zsize + + dim_names_2d(1) = "xaxis_1" + dim_names_2d(2) = "Time" + dim_names_3d(1) = "xaxis_1" + dim_names_3d(2) = "yaxis_2" + dim_names_3d(3) = "Time" + dim_names_3d2 = dim_names_3d + dim_names_3d2(2) = "yaxis_1" + dim_names_4d(1) = "xaxis_1" + dim_names_4d(2) = "yaxis_1" + dim_names_4d(3) = "zaxis_1" + dim_names_4d(4) = "Time" + dim_names_4d2 = dim_names_4d + dim_names_4d2(1) = "xaxis_2" + dim_names_4d2(2) = "yaxis_2" + dim_names_4d3 = dim_names_4d + dim_names_4d3(2) = "yaxis_2" + + ntprog = size(Atm%q,4) + ntdiag = size(Atm%qdiag,4) + ntracers = ntprog+ntdiag + + xpos = (/CENTER/) + xpos_2d = (/CENTER, EAST/) + ypos = (/CENTER/) + ypos_2d = (/NORTH, CENTER/) + + ! fname = 'fv_core.res.nc' + if (Atm%Fv_restart_is_open) then + call register_axis(Atm%Fv_restart, "xaxis_1", size(Atm%ak(:), 1)) + call register_axis(Atm%Fv_restart, "Time", unlimited) + if (.not. Atm%Fv_restart%is_readonly) then !if writing file + call register_field(Atm%Fv_restart, "xaxis_1", "double", (/"xaxis_1"/)) + call register_variable_attribute(Atm%Fv_restart,"xaxis_1", "axis", "X", str_len=1) + if (allocated(buffer)) deallocate(buffer) + allocate(buffer(size(Atm%ak(:), 1))) + do j = 1, size(Atm%ak(:), 1) + buffer(j) = j + end do + call write_data(Atm%Fv_restart, "xaxis_1", buffer) + deallocate(buffer) + call register_field(Atm%Fv_restart, "Time", "double", (/"Time"/)) + call register_variable_attribute(Atm%Fv_restart, dim_names_2d(2), "cartesian_axis", "T", str_len=1) + call register_variable_attribute(Atm%Fv_restart, dim_names_2d(2), "units", "time level", str_len=len("time level")) + call register_variable_attribute(Atm%Fv_restart, dim_names_2d(2), "long_name", dim_names_2d(2), str_len=len(dim_names_2d(2))) + call write_data(Atm%Fv_restart, "Time", 1) + endif + call register_restart_field (Atm%Fv_restart, 'ak', Atm%ak(:), dim_names_2d) + call register_restart_field (Atm%Fv_restart, 'bk', Atm%bk(:), dim_names_2d) + + ! fname= 'fv_core.res'//trim(stile_name)//'.nc' + elseif (Atm%Fv_restart_tile_is_open) then + zsize = (/size(Atm%u,3)/) + call fv_io_register_axis(Atm%Fv_restart_tile, numx=numx_2d, numy=numy_2d, xpos=xpos_2d, ypos=ypos_2d, numz=numz, zsize=zsize) + call register_restart_field(Atm%Fv_restart_tile, 'u', Atm%u, dim_names_4d) + call register_restart_field(Atm%Fv_restart_tile, 'v', Atm%v, dim_names_4d2) + + if (.not.Atm%flagstruct%hydrostatic) then + call register_restart_field(Atm%Fv_restart_tile, 'W', Atm%w, dim_names_4d3) + call register_restart_field(Atm%Fv_restart_tile, 'DZ', Atm%delz, dim_names_4d3) + if ( Atm%flagstruct%hybrid_z ) then + call register_restart_field(Atm%Fv_restart_tile, 'ZE0', Atm%ze0, dim_names_4d3) + endif + endif + call register_restart_field(Atm%Fv_restart_tile, 'T', Atm%pt, dim_names_4d3) + call register_restart_field(Atm%Fv_restart_tile, 'delp', Atm%delp, dim_names_4d3) + call register_restart_field(Atm%Fv_restart_tile, 'phis', Atm%phis, dim_names_3d) + + !--- include agrid winds in restarts for use in data assimilation + if (Atm%flagstruct%agrid_vel_rst) then + call register_restart_field(Atm%Fv_restart_tile, 'ua', Atm%ua, dim_names_4d3) + call register_restart_field(Atm%Fv_restart_tile, 'va', Atm%va, dim_names_4d3) + endif + + ! fname = 'fv_srf_wnd.res'//trim(stile_name)//'.nc + elseif (Atm%Rsf_restart_is_open) then + call fv_io_register_axis(Atm%Rsf_restart, numx=numx, numy=numy, xpos=xpos, ypos=ypos) + call register_restart_field(Atm%Rsf_restart, 'u_srf', Atm%u_srf, dim_names_3d2) + call register_restart_field(Atm%Rsf_restart, 'v_srf', Atm%v_srf, dim_names_3d2) +#ifdef SIM_PHYS + call register_restart_field(Atm%Rsf_restart, 'ts', Atm%ts, dim_names_3d2) +#endif + + ! fname = 'mg_drag.res'//trim(stile_name)//'.nc' + elseif (Atm%Mg_restart_is_open) then + call fv_io_register_axis(Atm%Mg_restart, numx=numx, numy=numy, xpos=xpos, ypos=ypos) + call register_restart_field (Atm%Mg_restart, 'ghprime', Atm%sgh, dim_names_3d2) + + ! fname = 'fv_land.res'//trim(stile_name)//'.nc' + elseif (Atm%Lnd_restart_is_open) then + call fv_io_register_axis(Atm%Lnd_restart, numx=numx, numy=numy, xpos=xpos, ypos=ypos) + call register_restart_field (Atm%Lnd_restart, 'oro', Atm%oro, dim_names_3d2) + + ! fname = 'fv_tracer.res'//trim(stile_name)//'.nc' + elseif (Atm%Tra_restart_is_open) then + zsize = (/size(Atm%q,3)/) + call fv_io_register_axis(Atm%Tra_restart, numx=numx, numy=numy, xpos=xpos, ypos=ypos, numz=numz, zsize=zsize) + do nt = 1, ntprog + call get_tracer_names(MODEL_ATMOS, nt, tracer_name) + if(Atm%Tra_restart%is_readonly) then !if reading file (don't do this if writing) + ! set all tracers to an initial profile value + call set_tracer_profile (MODEL_ATMOS, nt, Atm%q(:,:,:,nt) ) + endif + call register_restart_field(Atm%Tra_restart, tracer_name, Atm%q(:,:,:,nt), & + dim_names_4d, is_optional=.true.) + enddo + do nt = ntprog+1, ntracers + call get_tracer_names(MODEL_ATMOS, nt, tracer_name) + if(Atm%Tra_restart%is_readonly) then !if reading file (don't do this if writing) + ! set all tracers to an initial profile value + call set_tracer_profile (MODEL_ATMOS, nt, Atm%qdiag(:,:,:,nt) ) + endif + call register_restart_field(Atm%Tra_restart, tracer_name, Atm%qdiag(:,:,:,nt), & + dim_names_4d, is_optional=.true.) + enddo + endif + end subroutine fv_io_register_restart + ! NAME="fv_io_register_restart" + !##################################################################### ! @@ -113,18 +332,32 @@ subroutine fv_io_read_restart(fv_domain,Atm) type(domain2d), intent(inout) :: fv_domain type(fv_atmos_type), intent(inout) :: Atm(:) - character(len=64) :: fname, tracer_name - character(len=6) :: stile_name + character(len=64) :: tracer_name integer :: isc, iec, jsc, jec, n, nt, nk, ntracers integer :: ntileMe integer :: ks, ntiles real :: ptop - character(len=128) :: tracer_longname, tracer_units - - ntileMe = size(Atm(:)) ! This will need mods for more than 1 tile per pe + character(len=128) :: tracer_longname, tracer_units + character(len=120) :: fname + character(len=20) :: suffix + character(len=1) :: tile_num + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist + + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + + suffix = '' + fname = 'INPUT/fv_core.res.nc' + Atm(1)%Fv_restart_is_open = open_file(Atm(1)%Fv_restart,fname,"read", is_restart=.true., pelist=pes) + if (Atm(1)%Fv_restart_is_open) then + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Fv_restart) + call close_file(Atm(1)%Fv_restart) + Atm(1)%Fv_restart_is_open = .false. + endif + deallocate(pes) - call restore_state(Atm(1)%Fv_restart) if (Atm(1)%flagstruct%external_eta) then call set_external_eta(Atm(1)%ak, Atm(1)%bk, Atm(1)%ptop, Atm(1)%ks) endif @@ -134,53 +367,71 @@ subroutine fv_io_read_restart(fv_domain,Atm) !call restore_state(Atm(1)%SST_restart) endif -! fix for single tile runs where you need fv_core.res.nc and fv_core.res.tile1.nc ntiles = mpp_get_ntile_count(fv_domain) - if(ntiles == 1 .and. .not. Atm(1)%neststruct%nested) then - stile_name = '.tile1' - else - stile_name = '' + !If the number of tiles is equal to 1, and it is not a nested case add the ".tile1" suffix to the filename + if (ntiles == 1 .and. .not. Atm(1)%neststruct%nested) then + suffix = ''//trim(suffix)//'.tile1' endif - do n = 1, ntileMe - call restore_state(Atm(n)%Fv_tile_restart) + fname = 'INPUT/fv_core.res'//trim(suffix)//'.nc' + Atm(1)%Fv_restart_tile_is_open = open_file(Atm(1)%Fv_restart_tile, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Fv_restart_tile_is_open) then + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Fv_restart_tile) + call close_file(Atm(1)%Fv_restart_tile) + Atm(1)%Fv_restart_tile_is_open = .false. + endif !--- restore data for fv_tracer - if it exists - fname = 'INPUT/fv_tracer.res'//trim(stile_name)//'.nc' - if (file_exist(fname)) then - call restore_state(Atm(n)%Tra_restart) - else - call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') - endif + fname = 'INPUT/fv_tracer.res'//trim(suffix)//'.nc' + Atm(1)%Tra_restart_is_open = open_file(Atm(1)%Tra_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Tra_restart_is_open) then + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Tra_restart) + call close_file(Atm(1)%Tra_restart) + Atm(1)%Tra_restart_is_open = .false. + else + call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') + endif !--- restore data for surface winds - if it exists - fname = 'INPUT/fv_srf_wnd.res'//trim(stile_name)//'.nc' - if (file_exist(fname)) then - call restore_state(Atm(n)%Rsf_restart) - Atm(n)%flagstruct%srf_init = .true. - else - call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') - Atm(n)%flagstruct%srf_init = .false. - endif + fname = 'INPUT/fv_srf_wnd.res'//trim(suffix)//'.nc' + Atm(1)%Rsf_restart_is_open = open_file(Atm(1)%Rsf_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Rsf_restart_is_open) then + Atm(1)%flagstruct%srf_init = .true. + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Rsf_restart) + call close_file(Atm(1)%Rsf_restart) + Atm(1)%Rsf_restart_is_open = .false. + else + call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') + Atm(1)%flagstruct%srf_init = .false. + endif - if ( Atm(n)%flagstruct%fv_land ) then + if ( Atm(1)%flagstruct%fv_land ) then !--- restore data for mg_drag - if it exists - fname = 'INPUT/mg_drag.res'//trim(stile_name)//'.nc' - if (file_exist(fname)) then - call restore_state(Atm(n)%Mg_restart) + fname = 'INPUT/mg_drag.res'//trim(suffix)//'.nc' + Atm(1)%Mg_restart_is_open = open_file(Atm(1)%Mg_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Mg_restart_is_open) then + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Mg_restart) + call close_file(Atm(1)%Mg_restart) + Atm(1)%Mg_restart_is_open = .false. else call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') endif !--- restore data for fv_land - if it exists - fname = 'INPUT/fv_land.res'//trim(stile_name)//'.nc' - if (file_exist(fname)) then - call restore_state(Atm(n)%Lnd_restart) + fname = 'INPUT/fv_land.res'//trim(suffix)//'.nc' + Atm(1)%Lnd_restart_is_open = open_file(Atm(1)%Lnd_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Lnd_restart_is_open) then + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Lnd_restart) + call close_file(Atm(1)%Lnd_restart) + Atm(1)%Lnd_restart_is_open = .false. else call mpp_error(NOTE,'==> Warning from fv_read_restart: Expected file '//trim(fname)//' does not exist') endif - endif - - end do + endif return @@ -192,17 +443,16 @@ end subroutine fv_io_read_restart subroutine fv_io_read_tracers(fv_domain,Atm) type(domain2d), intent(inout) :: fv_domain type(fv_atmos_type), intent(inout) :: Atm(:) - integer :: n, ntracers, ntprog, nt, isc, iec, jsc, jec, id_restart + integer :: ntracers, ntprog, nt, isc, iec, jsc, jec character(len=6) :: stile_name character(len=64):: fname, tracer_name - type(restart_file_type) :: Tra_restart_r + type(FmsNetcdfDomainFile_t) :: Tra_restart_r integer :: ntiles - n = 1 - isc = Atm(n)%bd%isc - iec = Atm(n)%bd%iec - jsc = Atm(n)%bd%jsc - jec = Atm(n)%bd%jec + isc = Atm(1)%bd%isc + iec = Atm(1)%bd%iec + jsc = Atm(1)%bd%jsc + jec = Atm(1)%bd%jec call get_number_tracers(MODEL_ATMOS, num_tracers=ntracers, num_prog=ntprog) ! fix for single tile runs where you need fv_core.res.nc and fv_core.res.tile1.nc @@ -213,26 +463,27 @@ subroutine fv_io_read_tracers(fv_domain,Atm) stile_name = '' endif - fname = 'fv_tracer.res'//trim(stile_name)//'.nc' + fname = 'INPUT/fv_tracer.res'//trim(stile_name)//'.nc' + + if (open_file(Tra_restart_r,fname,"read",fv_domain, is_restart=.true.)) then do nt = 2, ntprog call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - call set_tracer_profile (MODEL_ATMOS, nt, Atm(n)%q(isc:iec,jsc:jec,:,nt) ) - id_restart = register_restart_field(Tra_restart_r, fname, tracer_name, Atm(n)%q(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) + call set_tracer_profile (MODEL_ATMOS, nt, Atm(1)%q(isc:iec,jsc:jec,:,nt) ) + if (variable_exists(Tra_restart_r, tracer_name)) then + call read_data(Tra_restart_r, tracer_name, Atm(1)%q(:,:,:,nt)) + endif enddo do nt = ntprog+1, ntracers call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - call set_tracer_profile (MODEL_ATMOS, nt, Atm(n)%qdiag(isc:iec,jsc:jec,:,nt) ) - id_restart = register_restart_field(Tra_restart_r, fname, tracer_name, Atm(n)%qdiag(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) + call set_tracer_profile (MODEL_ATMOS, nt, Atm(1)%qdiag(isc:iec,jsc:jec,:,nt) ) + if (variable_exists(Tra_restart_r, tracer_name)) then + call read_data (Tra_restart_r, tracer_name, Atm(1)%qdiag(:,:,:,nt)) + endif enddo - if (file_exist('INPUT'//trim(fname))) then - call restore_state(Tra_restart_r) - call free_restart_type(Tra_restart_r) + call close_file(Tra_restart_r) else call mpp_error(NOTE,'==> Warning from fv_io_read_tracers: Expected file '//trim(fname)//' does not exist') endif - return end subroutine fv_io_read_tracers @@ -249,8 +500,10 @@ subroutine remap_restart(fv_domain,Atm) integer :: isc, iec, jsc, jec, n, nt, nk, ntracers, ntprog, ntdiag integer :: isd, ied, jsd, jed integer :: ntiles - type(restart_file_type) :: FV_restart_r, FV_tile_restart_r, Tra_restart_r - integer :: id_restart + + type(FmsNetcdfDomainFile_t) :: FV_tile_restart_r, Tra_restart_r + type(FmsNetcdfFile_t) :: Fv_restart_r + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist ! !------------------------------------------------------------------------- @@ -298,11 +551,15 @@ subroutine remap_restart(fv_domain,Atm) allocate ( ze0_r(isc:iec, jsc:jec, npz_rst+1) ) endif - fname = 'fv_core.res.nc' - id_restart = register_restart_field(Fv_restart_r, fname, 'ak', ak_r(:), no_domain=.true.) - id_restart = register_restart_field(Fv_restart_r, fname, 'bk', bk_r(:), no_domain=.true.) - call restore_state(Fv_restart_r) - call free_restart_type(Fv_restart_r) + fname = 'INPUT/fv_core.res.nc' + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + if (open_file(Fv_restart_r,fname,"read", is_restart=.true., pelist=pes)) then + call read_data(Fv_restart_r, 'ak', ak_r(:)) + call read_data(Fv_restart_r, 'bk', bk_r(:)) + call close_file(Fv_restart_r) + endif + deallocate(pes) ! fix for single tile runs where you need fv_core.res.nc and fv_core.res.tile1.nc ntiles = mpp_get_ntile_count(fv_domain) @@ -316,85 +573,85 @@ subroutine remap_restart(fv_domain,Atm) !!! file_exist() needs the full relative path, including INPUT/ !!! But register_restart_field ONLY looks in INPUT/ and so JUST needs the file name!! -! do n = 1, ntileMe - n = 1 - fname = 'fv_core.res'//trim(stile_name)//'.nc' - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'u', u_r, & - domain=fv_domain, position=NORTH,tile_count=n) - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'v', v_r, & - domain=fv_domain, position=EAST,tile_count=n) - if (.not.Atm(n)%flagstruct%hydrostatic) then - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'W', w_r, & - domain=fv_domain, mandatory=.false., tile_count=n) - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'DZ', delz_r, & - domain=fv_domain, mandatory=.false., tile_count=n) - if ( Atm(n)%flagstruct%hybrid_z ) then - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'ZE0', ze0_r, & - domain=fv_domain, mandatory=.false., tile_count=n) + fname = 'INPUT/fv_core.res'//trim(stile_name)//'.nc' + if (open_file(Fv_tile_restart_r, fname, "read", fv_domain, is_restart=.true.)) then + call read_data(Fv_tile_restart_r, 'u', u_r) + call read_data(Fv_tile_restart_r, 'v', v_r) + if (.not.Atm(1)%flagstruct%hydrostatic) then + call read_data(Fv_tile_restart_r, 'W', w_r) + call read_data(Fv_tile_restart_r, 'DZ', delz_r) + if ( Atm(1)%flagstruct%hybrid_z ) then + call read_data(Fv_tile_restart_r, 'ZE0', ze0_r) + endif endif + call read_data(Fv_tile_restart_r, 'T', pt_r) + call read_data(Fv_tile_restart_r, 'delp', delp_r) + call read_data(Fv_tile_restart_r, 'phis', Atm(1)%phis) + call close_file(FV_tile_restart_r) endif - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'T', pt_r, & - domain=fv_domain, tile_count=n) - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'delp', delp_r, & - domain=fv_domain, tile_count=n) - id_restart = register_restart_field(Fv_tile_restart_r, fname, 'phis', Atm(n)%phis, & - domain=fv_domain, tile_count=n) - call restore_state(FV_tile_restart_r) - call free_restart_type(FV_tile_restart_r) - fname = 'fv_srf_wnd.res'//trim(stile_name)//'.nc' - if (file_exist('INPUT/'//fname)) then - call restore_state(Atm(n)%Rsf_restart) - Atm(n)%flagstruct%srf_init = .true. - else + + fname = 'INPUT/fv_srf_wnd.res'//trim(stile_name)//'.nc' + Atm(1)%Rsf_restart_is_open = open_file(Atm(1)%Rsf_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Rsf_restart_is_open) then + Atm%flagstruct%srf_init = .true. + call fv_io_register_restart(Atm(1)) + call read_restart(Atm(1)%Rsf_restart) + call close_file(Atm(1)%Rsf_restart) + Atm(1)%Rsf_restart_is_open = .false. call mpp_error(NOTE,'==> Warning from remap_restart: Expected file '//trim(fname)//' does not exist') - Atm(n)%flagstruct%srf_init = .false. + Atm%flagstruct%srf_init = .false. endif - if ( Atm(n)%flagstruct%fv_land ) then + if ( Atm(1)%flagstruct%fv_land ) then !--- restore data for mg_drag - if it exists - fname = 'mg_drag.res'//trim(stile_name)//'.nc' - if (file_exist('INPUT/'//fname)) then - call restore_state(Atm(n)%Mg_restart) + fname = 'INPUT/mg_drag.res'//trim(stile_name)//'.nc' + Atm(1)%Mg_restart_is_open = open_file(Atm(1)%Mg_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Mg_restart_is_open) then + call read_data(Atm(1)%Mg_restart, 'ghprime', Atm(1)%sgh) + call close_file(Atm(1)%Mg_restart) else call mpp_error(NOTE,'==> Warning from remap_restart: Expected file '//trim(fname)//' does not exist') endif !--- restore data for fv_land - if it exists - fname = 'fv_land.res'//trim(stile_name)//'.nc' - if (file_exist('INPUT/'//fname)) then - call restore_state(Atm(n)%Lnd_restart) + fname = 'INPUT/fv_land.res'//trim(stile_name)//'.nc' + Atm(1)%Lnd_restart_is_open = open_file(Atm(1)%Lnd_restart, fname, "read", fv_domain, is_restart=.true.) + if (Atm(1)%Lnd_restart_is_open) then + call read_data(Atm(1)%Lnd_restart, 'oro', Atm(1)%oro) + call close_file(Atm(1)%Lnd_restart) else call mpp_error(NOTE,'==> Warning from remap_restart: Expected file '//trim(fname)//' does not exist') endif endif - fname = 'fv_tracer.res'//trim(stile_name)//'.nc' - if (file_exist('INPUT/'//fname)) then + fname = 'INPUT/fv_tracer.res'//trim(stile_name)//'.nc' + if (open_file(Tra_restart_r, fname, "read", fv_domain, is_restart=.true.)) then do nt = 1, ntprog call get_tracer_names(MODEL_ATMOS, nt, tracer_name) call set_tracer_profile (MODEL_ATMOS, nt, q_r(isc:iec,jsc:jec,:,nt) ) - id_restart = register_restart_field(Tra_restart_r, fname, tracer_name, q_r(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) + if (variable_exists(Tra_restart_r, tracer_name)) then + call read_data(Tra_restart_r, tracer_name, q_r(:,:,:,nt)) + endif enddo do nt = ntprog+1, ntracers call get_tracer_names(MODEL_ATMOS, nt, tracer_name) call set_tracer_profile (MODEL_ATMOS, nt, qdiag_r(isc:iec,jsc:jec,:,nt) ) - id_restart = register_restart_field(Tra_restart_r, fname, tracer_name, qdiag_r(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) + if (variable_exists(Tra_restart_r, tracer_name)) then + call read_data (Tra_restart_r, tracer_name, qdiag_r(:,:,:,nt)) + endif enddo - call restore_state(Tra_restart_r) - call free_restart_type(Tra_restart_r) + call close_file(Tra_restart_r) else call mpp_error(NOTE,'==> Warning from remap_restart: Expected file '//trim(fname)//' does not exist') endif ! ====== PJP added DA functionality ====== - if (Atm(n)%flagstruct%read_increment) then + if (Atm(1)%flagstruct%read_increment) then ! print point in middle of domain for a sanity check i = (isc + iec)/2 j = (jsc + jec)/2 k = npz_rst/2 if( is_master() ) write(*,*) 'Calling read_da_inc',pt_r(i,j,k) - call read_da_inc(Atm(n), Atm(n)%domain, Atm(n)%bd, npz_rst, ntprog, & + call read_da_inc(Atm(1), Atm(1)%domain, Atm(1)%bd, npz_rst, ntprog, & u_r, v_r, q_r, delp_r, pt_r, isc, jsc, iec, jec ) if( is_master() ) write(*,*) 'Back from read_da_inc',pt_r(i,j,k) endif @@ -439,144 +696,18 @@ subroutine fv_io_register_nudge_restart(Atm) ! use_ncep_sst may not be initialized at this point? call mpp_error(NOTE, 'READING FROM SST_restart DISABLED') - if ( use_ncep_sst .or. Atm(1)%flagstruct%nudge .or. Atm(1)%flagstruct%ncep_ic ) then -! if ( Atm(1)%nudge .or. Atm(1)%ncep_ic ) then - fname = 'sst_ncep.res.nc' - id_restart = register_restart_field(Atm(1)%SST_restart, fname, 'sst_ncep', sst_ncep) - id_restart = register_restart_field(Atm(1)%SST_restart, fname, 'sst_anom', sst_anom) - endif - - end subroutine fv_io_register_nudge_restart - ! NAME="fv_io_register_nudge_restart" - - - !##################################################################### - ! - ! - ! - ! register restart field to be written out to restart file. - ! - subroutine fv_io_register_restart(fv_domain,Atm) - type(domain2d), intent(inout) :: fv_domain - type(fv_atmos_type), intent(inout) :: Atm(:) - - character(len=64) :: fname, tracer_name - character(len=6) :: gn, stile_name - integer :: id_restart - integer :: n, nt, ntracers, ntprog, ntdiag, ntileMe, ntiles - - ntileMe = size(Atm(:)) - ntprog = size(Atm(1)%q,4) - ntdiag = size(Atm(1)%qdiag,4) - ntracers = ntprog+ntdiag - -!--- set the 'nestXX' appendix for all files using fms_io - if (Atm(1)%grid_number > 1) then - write(gn,'(A4, I2.2)') "nest", Atm(1)%grid_number - else - gn = '' - end if - call set_filename_appendix(gn) - -!--- fix for single tile runs where you need fv_core.res.nc and fv_core.res.tile1.nc - ntiles = mpp_get_ntile_count(fv_domain) - if(ntiles == 1 .and. .not. Atm(1)%neststruct%nested) then - stile_name = '.tile1' - else - stile_name = '' - endif - -! use_ncep_sst may not be initialized at this point? -#ifndef DYCORE_SOLO -! call mpp_error(NOTE, 'READING FROM SST_RESTART DISABLED') -!!$ if ( use_ncep_sst .or. Atm(1)%flagstruct%nudge .or. Atm(1)%flagstruct%ncep_ic ) then -!!$ fname = 'sst_ncep'//trim(gn)//'.res.nc' +!!$ if ( use_ncep_sst .or. Atm(1)%flagstruct%nudge .or. Atm(1)%flagstruct%ncep_ic ) then +!!$ if ( Atm(1)%nudge .or. Atm(1)%ncep_ic ) then +!!$ fname = 'sst_ncep.res.nc' !!$ id_restart = register_restart_field(Atm(1)%SST_restart, fname, 'sst_ncep', sst_ncep) !!$ id_restart = register_restart_field(Atm(1)%SST_restart, fname, 'sst_anom', sst_anom) -!!$ endif -#endif - - fname = 'fv_core.res.nc' - id_restart = register_restart_field(Atm(1)%Fv_restart, fname, 'ak', Atm(1)%ak(:), no_domain=.true.) - id_restart = register_restart_field(Atm(1)%Fv_restart, fname, 'bk', Atm(1)%bk(:), no_domain=.true.) - - do n = 1, ntileMe - fname = 'fv_core.res'//trim(stile_name)//'.nc' - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'u', Atm(n)%u, & - domain=fv_domain, position=NORTH,tile_count=n) - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'v', Atm(n)%v, & - domain=fv_domain, position=EAST,tile_count=n) - if (.not.Atm(n)%flagstruct%hydrostatic) then - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'W', Atm(n)%w, & - domain=fv_domain, mandatory=.false., tile_count=n) - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'DZ', Atm(n)%delz, & - domain=fv_domain, mandatory=.false., tile_count=n) - if ( Atm(n)%flagstruct%hybrid_z ) then - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'ZE0', Atm(n)%ze0, & - domain=fv_domain, mandatory=.false., tile_count=n) - endif - endif - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'T', Atm(n)%pt, & - domain=fv_domain, tile_count=n) - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'delp', Atm(n)%delp, & - domain=fv_domain, tile_count=n) - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'phis', Atm(n)%phis, & - domain=fv_domain, tile_count=n) - - !--- include agrid winds in restarts for use in data assimilation - if (Atm(n)%flagstruct%agrid_vel_rst) then - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'ua', Atm(n)%ua, & - domain=fv_domain, tile_count=n, mandatory=.false.) - id_restart = register_restart_field(Atm(n)%Fv_tile_restart, fname, 'va', Atm(n)%va, & - domain=fv_domain, tile_count=n, mandatory=.false.) - endif - - fname = 'fv_srf_wnd.res'//trim(stile_name)//'.nc' - id_restart = register_restart_field(Atm(n)%Rsf_restart, fname, 'u_srf', Atm(n)%u_srf, & - domain=fv_domain, tile_count=n) - id_restart = register_restart_field(Atm(n)%Rsf_restart, fname, 'v_srf', Atm(n)%v_srf, & - domain=fv_domain, tile_count=n) -#ifdef SIM_PHYS - id_restart = register_restart_field(Rsf_restart(n), fname, 'ts', Atm(n)%ts, & - domain=fv_domain, tile_count=n) -#endif - - if ( Atm(n)%flagstruct%fv_land ) then - !------------------------------------------------------------------------------------------------- - ! Optional terrain deviation (sgh) and land fraction (oro) - fname = 'mg_drag.res'//trim(stile_name)//'.nc' - id_restart = register_restart_field(Atm(n)%Mg_restart, fname, 'ghprime', Atm(n)%sgh, & - domain=fv_domain, tile_count=n) - - fname = 'fv_land.res'//trim(stile_name)//'.nc' - id_restart = register_restart_field(Atm(n)%Lnd_restart, fname, 'oro', Atm(n)%oro, & - domain=fv_domain, tile_count=n) - endif +!!$ endif - fname = 'fv_tracer.res'//trim(stile_name)//'.nc' - do nt = 1, ntprog - call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - ! set all tracers to an initial profile value - call set_tracer_profile (MODEL_ATMOS, nt, Atm(n)%q(:,:,:,nt) ) - id_restart = register_restart_field(Atm(n)%Tra_restart, fname, tracer_name, Atm(n)%q(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) - enddo - do nt = ntprog+1, ntracers - call get_tracer_names(MODEL_ATMOS, nt, tracer_name) - ! set all tracers to an initial profile value - call set_tracer_profile (MODEL_ATMOS, nt, Atm(n)%qdiag(:,:,:,nt) ) - id_restart = register_restart_field(Atm(n)%Tra_restart, fname, tracer_name, Atm(n)%qdiag(:,:,:,nt), & - domain=fv_domain, mandatory=.false., tile_count=n) - enddo + end subroutine fv_io_register_nudge_restart + ! NAME="fv_io_register_nudge_restart" - if ( Atm(n)%neststruct%nested ) then - call fv_io_register_restart_BCs(Atm(n)) !TODO put into fv_io_register_restart - endif - enddo - end subroutine fv_io_register_restart - ! NAME="fv_io_register_restart" @@ -590,6 +721,14 @@ subroutine fv_io_write_restart(Atm, timestamp) type(fv_atmos_type), intent(inout) :: Atm character(len=*), optional, intent(in) :: timestamp + integer :: ntiles + logical :: tile_file_exists + type(domain2d) :: fv_domain + character(len=120) :: fname + character(len=20) :: suffix + character(len=1) :: tile_num + integer, allocatable, dimension(:) :: pes !< Array of the pes in the current pelist + fv_domain = Atm%domain !!$ if ( use_ncep_sst .or. Atm%flagstruct%nudge .or. Atm%flagstruct%ncep_ic ) then !!$ call mpp_error(NOTE, 'READING FROM SST_RESTART DISABLED') @@ -597,27 +736,106 @@ subroutine fv_io_write_restart(Atm, timestamp) !!$ endif if ( (use_ncep_sst .or. Atm%flagstruct%nudge) .and. .not. Atm%gridstruct%nested ) then - call save_restart(Atm%SST_restart, timestamp) + !call save_restart(Atm%SST_restart, timestamp) endif - call save_restart(Atm%Fv_restart, timestamp) - call save_restart(Atm%Fv_tile_restart, timestamp) - call save_restart(Atm%Rsf_restart, timestamp) + suffix = '' + if (present(timestamp)) then + fname = 'RESTART/'//trim(timestamp)//'.fv_core.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/fv_core.res'//trim(suffix)//'.nc' + endif + allocate(pes(mpp_npes())) + call mpp_get_current_pelist(pes) + Atm%Fv_restart_is_open = open_file(Atm%Fv_restart, fname, "overwrite", is_restart=.true., pelist=pes) + if (Atm%Fv_restart_is_open) then + call fv_io_register_restart(Atm) + call write_restart(Atm%Fv_restart) + call close_file(Atm%Fv_restart) + Atm%Fv_restart_is_open = .false. + endif + deallocate(pes) - if ( Atm%flagstruct%fv_land ) then - call save_restart(Atm%Mg_restart, timestamp) - call save_restart(Atm%Lnd_restart, timestamp) + ntiles = mpp_get_ntile_count(fv_domain) + !If the number of tiles is equal to 1, and it is not a nested case add the ".tile1" suffix to the filename + if (ntiles == 1 .and. .not. Atm%neststruct%nested) then + suffix = ''//trim(suffix)//'.tile1' + endif + + if (present(timestamp)) then + fname = 'RESTART/'//trim(timestamp)//'.fv_core.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/fv_core.res'//trim(suffix)//'.nc' + endif + + Atm%Fv_restart_tile_is_open = open_file(Atm%Fv_restart_tile, fname, "overwrite", fv_domain, is_restart=.true.) + if (Atm%Fv_restart_tile_is_open) then + call fv_io_register_restart(Atm) + call write_restart (Atm%Fv_restart_tile) + call close_file (Atm%Fv_restart_tile) + Atm%Fv_restart_tile_is_open = .false. endif - call save_restart(Atm%Tra_restart, timestamp) + if (present(timestamp)) then + fname = 'RESTART/'//trim(timestamp)//'.fv_srf_wnd.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/fv_srf_wnd.res'//trim(suffix)//'.nc' + endif + Atm%Rsf_restart_is_open = open_file(Atm%Rsf_restart, fname, "overwrite", fv_domain, is_restart=.true.) + if (Atm%Rsf_restart_is_open) then + call fv_io_register_restart(Atm) + call write_restart (Atm%Rsf_restart) + call close_file (Atm%Rsf_restart) + Atm%Rsf_restart_is_open = .false. + endif + + if ( Atm%flagstruct%fv_land ) then + if (present(timestamp)) then + fname = 'RESTART/'//trim(timestamp)//'.mg_drag.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/mg_drag.res'//trim(suffix)//'.nc' + endif + Atm%Mg_restart_is_open = open_file(Atm%Mg_restart, fname, "overwrite", fv_domain, is_restart=.true.) + if (Atm%Mg_restart_is_open) then + call fv_io_register_restart(Atm) + call write_restart(Atm%Mg_restart) + call close_file(Atm%Mg_restart) + Atm%Mg_restart_is_open = .false. + endif + + if (present(timestamp)) then + fname = 'RESTART'//trim(timestamp)//'./fv_land.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/fv_land.res'//trim(suffix)//'.nc' + endif + Atm%Lnd_restart_is_open = open_file(Atm%Lnd_restart, fname, "overwrite", fv_domain, is_restart=.true.) + if (Atm%Lnd_restart_is_open) then + call fv_io_register_restart(Atm) + call write_restart(Atm%Lnd_restart) + call close_file(Atm%Lnd_restart) + Atm%Lnd_restart_is_open = .false. + endif + endif + if (present(timestamp)) then + fname = 'RESTART/'//trim(timestamp)//'.fv_tracer.res'//trim(suffix)//'.nc' + else + fname = 'RESTART/fv_tracer.res'//trim(suffix)//'.nc' + endif + Atm%Tra_restart_is_open = open_file(Atm%Tra_restart, fname, "overwrite", fv_domain, is_restart=.true.) + if (Atm%Tra_restart_is_open) then + call fv_io_register_restart(Atm) + call write_restart(Atm%Tra_restart) + call close_file(Atm%Tra_restart) + Atm%Tra_restart_is_open = .false. + endif end subroutine fv_io_write_restart subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & var_name, var, var_bc, istag, jstag) type(fv_atmos_type), intent(in) :: Atm - type(restart_file_type), intent(inout) :: BCfile_ne, BCfile_sw + type(FmsNetcdfFile_t), intent(inout) :: BCfile_ne, BCfile_sw character(len=120), intent(in) :: fname_ne, fname_sw character(len=*), intent(in) :: var_name real, dimension(:,:), intent(in), optional :: var @@ -675,13 +893,13 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (is.eq.1 .and. js.eq.1) is_root_pe = .TRUE. !register west halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var_bc) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_west_t1', & var_bc%west_t1, & indices, global_size, y2_pelist, & is_root_pe, jshift=y_halo) !register west prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_west', & var, indices, global_size, & y2_pelist, is_root_pe, jshift=y_halo) @@ -690,7 +908,7 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (ie.eq.npx-1 .and. je.eq.npy-1) is_root_pe = .TRUE. !register east halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var_bc) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_east_t1', & var_bc%east_t1, & indices, global_size, y1_pelist, & @@ -700,7 +918,7 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & indices(1) = ied-x_halo+1+i_stag indices(2) = ied+i_stag !register east prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_east', & var, indices, global_size, & y1_pelist, is_root_pe, jshift=y_halo, & @@ -724,13 +942,13 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (is.eq.1 .and. js.eq.1) is_root_pe = .TRUE. !register south halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var_bc) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_south_t1', & var_bc%south_t1, & indices, global_size, x2_pelist, & is_root_pe, x_halo=x_halo_ns) !register south prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_south', & var, indices, global_size, & x2_pelist, is_root_pe, x_halo=x_halo_ns) @@ -739,7 +957,7 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (ie.eq.npx-1 .and. je.eq.npy-1) is_root_pe = .TRUE. !register north halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var_bc) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_north_t1', & var_bc%north_t1, & indices, global_size, x1_pelist, & @@ -749,19 +967,24 @@ subroutine register_bcs_2d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & indices(3) = jed-y_halo+1+j_stag indices(4) = jed+j_stag !register north prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_north', & var, indices, global_size, & x1_pelist, is_root_pe, x_halo=x_halo_ns, & y_halo=(size(var,2)-y_halo), jshift=-(je+j_stag)) + deallocate (x1_pelist) + deallocate (y1_pelist) + deallocate (x2_pelist) + deallocate (y2_pelist) + end subroutine register_bcs_2d subroutine register_bcs_3d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & var_name, var, var_bc, istag, jstag, mandatory) type(fv_atmos_type), intent(in) :: Atm - type(restart_file_type), intent(inout) :: BCfile_ne, BCfile_sw + type(FmsNetcdfFile_t), intent(inout) :: BCfile_ne, BCfile_sw character(len=120), intent(in) :: fname_ne, fname_sw character(len=*), intent(in) :: var_name real, dimension(:,:,:), intent(in), optional :: var @@ -776,6 +999,10 @@ subroutine register_bcs_3d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & integer, allocatable, dimension(:) :: x1_pelist, y1_pelist integer, allocatable, dimension(:) :: x2_pelist, y2_pelist logical :: is_root_pe + logical :: mandatory_flag + + mandatory_flag = .true. + if (present(mandatory)) mandatory_flag = mandatory i_stag = 0 j_stag = 0 @@ -821,36 +1048,36 @@ subroutine register_bcs_3d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (is.eq.1 .and. js.eq.1) is_root_pe = .TRUE. !register west halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var_bc) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_west_t1', & var_bc%west_t1, & indices, global_size, y2_pelist, & - is_root_pe, jshift=y_halo, mandatory=mandatory) + is_root_pe, jshift=y_halo, is_optional=.not.mandatory_flag) !register west prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_west', & var, indices, global_size, & - y2_pelist, is_root_pe, jshift=y_halo, mandatory=mandatory) + y2_pelist, is_root_pe, jshift=y_halo, is_optional=.not.mandatory_flag) !define east root_pe is_root_pe = .FALSE. if (ie.eq.npx-1 .and. je.eq.npy-1) is_root_pe = .TRUE. !register east halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var_bc) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_east_t1', & var_bc%east_t1, & indices, global_size, y1_pelist, & - is_root_pe, jshift=y_halo, mandatory=mandatory) + is_root_pe, jshift=y_halo, is_optional=.not.mandatory_flag) !reset indices for prognostic variables in the east halo indices(1) = ied-x_halo+1+i_stag indices(2) = ied+i_stag !register east prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_east', & var, indices, global_size, & y1_pelist, is_root_pe, jshift=y_halo, & - x_halo=(size(var,1)-x_halo), ishift=-(ie+i_stag), mandatory=mandatory) + x_halo=(size(var,1)-x_halo), ishift=-(ie+i_stag), is_optional=.not.mandatory_flag) !NORTH & SOUTH !set defaults for north/south halo regions @@ -871,41 +1098,45 @@ subroutine register_bcs_3d(Atm, BCfile_ne, BCfile_sw, fname_ne, fname_sw, & is_root_pe = .FALSE. if (is.eq.1 .and. js.eq.1) is_root_pe = .TRUE. !register south halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var_bc) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_south_t1', & var_bc%south_t1, & indices, global_size, x2_pelist, & - is_root_pe, x_halo=x_halo_ns, mandatory=mandatory) + is_root_pe, x_halo=x_halo_ns, is_optional=.not.mandatory_flag) !register south prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_sw, trim(fname_sw), & + if (present(var) .and. Atm%neststruct%BCfile_sw_is_open) call register_restart_field(BCfile_sw, & trim(var_name)//'_south', & var, indices, global_size, & - x2_pelist, is_root_pe, x_halo=x_halo_ns, mandatory=mandatory) + x2_pelist, is_root_pe, x_halo=x_halo_ns, is_optional=.not.mandatory_flag) !define north root_pe is_root_pe = .FALSE. if (ie.eq.npx-1 .and. je.eq.npy-1) is_root_pe = .TRUE. !register north halo data in t1 - if (present(var_bc)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var_bc) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_north_t1', & var_bc%north_t1, & indices, global_size, x1_pelist, & - is_root_pe, x_halo=x_halo_ns, mandatory=mandatory) + is_root_pe, x_halo=x_halo_ns, is_optional=.not.mandatory_flag) !reset indices for prognostic variables in the north halo indices(3) = jed-y_halo+1+j_stag indices(4) = jed+j_stag !register north prognostic halo data - if (present(var)) id_restart = register_restart_field(BCfile_ne, trim(fname_ne), & + if (present(var) .and. Atm%neststruct%BCfile_ne_is_open) call register_restart_field(BCfile_ne, & trim(var_name)//'_north', & var, indices, global_size, & x1_pelist, is_root_pe, x_halo=x_halo_ns, & - y_halo=(size(var,2)-y_halo), jshift=-(je+j_stag), mandatory=mandatory) + y_halo=(size(var,2)-y_halo), jshift=-(je+j_stag), is_optional=.not.mandatory_flag) + deallocate (x1_pelist) + deallocate (y1_pelist) + deallocate (x2_pelist) + deallocate (y2_pelist) end subroutine register_bcs_3d - ! NAME="fv_io_regsiter_restart_BCs" + ! NAME="fv_io_register_restart_BCs" !##################################################################### subroutine fv_io_register_restart_BCs(Atm) @@ -921,8 +1152,6 @@ subroutine fv_io_register_restart_BCs(Atm) ntdiag=size(Atm%qdiag,4) ntracers=ntprog+ntdiag - call set_domain(Atm%domain) - call register_bcs_2d(Atm, Atm%neststruct%BCfile_ne, Atm%neststruct%BCfile_sw, & fname_ne, fname_sw, 'phis', var=Atm%phis) call register_bcs_3d(Atm, Atm%neststruct%BCfile_ne, Atm%neststruct%BCfile_sw, & @@ -973,11 +1202,40 @@ end subroutine fv_io_register_restart_BCs subroutine fv_io_write_BCs(Atm, timestamp) - type(fv_atmos_type), intent(inout) :: Atm + type(fv_atmos_type), intent(inout) :: Atm character(len=*), intent(in), optional :: timestamp + integer, allocatable, dimension(:) :: all_pelist + integer, dimension(2) :: layout + integer :: n + character(len=1) :: tile_num + character(len=120) :: fname_ne, fname_sw + + if (present(timestamp)) then + fname_ne = 'RESTART/'//trim(timestamp)//'fv_BC_ne.res.nc' + fname_sw = 'RESTART/'//trim(timestamp)//'fv_BC_sw.res.nc' + else + fname_ne = 'RESTART/fv_BC_ne.res.nc' + fname_sw = 'RESTART/fv_BC_sw.res.nc' + endif + + allocate(all_pelist(mpp_npes())) + call mpp_get_current_pelist(all_pelist) + + Atm%neststruct%BCfile_sw_is_open = open_file(Atm%neststruct%BCfile_sw, fname_sw, "overwrite", is_restart=.true., pelist=all_pelist) + Atm%neststruct%BCfile_ne_is_open = open_file(Atm%neststruct%BCfile_ne, fname_ne, "overwrite", is_restart=.true., pelist=all_pelist) + call fv_io_register_restart_BCs(Atm) - call save_restart_border(Atm%neststruct%BCfile_ne, timestamp) - call save_restart_border(Atm%neststruct%BCfile_sw, timestamp) + if (Atm%neststruct%BCfile_sw_is_open) then + call write_restart_bc(Atm%neststruct%BCfile_sw) + call close_file(Atm%neststruct%BCfile_sw) + endif + + if (Atm%neststruct%BCfile_ne_is_open) then + call write_restart_bc(Atm%neststruct%BCfile_ne) + call close_file(Atm%neststruct%BCfile_ne) + endif + + deallocate(all_pelist) return end subroutine fv_io_write_BCs @@ -985,9 +1243,32 @@ end subroutine fv_io_write_BCs subroutine fv_io_read_BCs(Atm) type(fv_atmos_type), intent(inout) :: Atm + integer, allocatable, dimension(:) :: all_pelist + integer, dimension(2) :: layout + integer :: n + character(len=1) :: tile_num + character(len=120) :: fname_ne, fname_sw + + fname_ne = 'RESTART/fv_BC_ne.res.nc' + fname_sw = 'RESTART/fv_BC_sw.res.nc' + + allocate(all_pelist(mpp_npes())) + call mpp_get_current_pelist(all_pelist) + + Atm%neststruct%BCfile_sw_is_open = open_file(Atm%neststruct%BCfile_sw, fname_sw, "read", is_restart=.true., pelist=all_pelist) + Atm%neststruct%BCfile_ne_is_open = open_file(Atm%neststruct%BCfile_ne, fname_ne, "read", is_restart=.true., pelist=all_pelist) + call fv_io_register_restart_BCs(Atm) + + if (Atm%neststruct%BCfile_sw_is_open) then + call read_restart_bc(Atm%neststruct%BCfile_sw) + call close_file(Atm%neststruct%BCfile_sw) + endif + + if (Atm%neststruct%BCfile_ne_is_open) then + call read_restart_bc(Atm%neststruct%BCfile_ne) + call close_file(Atm%neststruct%BCfile_ne) + endif - call restore_state_border(Atm%neststruct%BCfile_ne) - call restore_state_border(Atm%neststruct%BCfile_sw) !These do not work yet !need to modify register_bcs_?d to get ids for registered variables, and then use query_initialized_id @@ -996,6 +1277,8 @@ subroutine fv_io_read_BCs(Atm) !Atm%neststruct%delz_BC%initialized = field_exist(fname_ne, 'delz_north_t1', Atm%domain) !if (is_master()) print*, ' BCs: ', Atm%neststruct%divg_BC%initialized, Atm%neststruct%w_BC%initialized, Atm%neststruct%delz_BC%initialized + deallocate(all_pelist) + return end subroutine fv_io_read_BCs diff --git a/tools/fv_mp_mod.F90 b/tools/fv_mp_mod.F90 index 7827caa71..fb7dcef0c 100644 --- a/tools/fv_mp_mod.F90 +++ b/tools/fv_mp_mod.F90 @@ -26,7 +26,7 @@ module fv_mp_mod #if defined(SPMD) ! !USES: - use fms_mod, only : fms_init, fms_end + use fms_mod, only : fms_end use mpp_mod, only : FATAL, MPP_DEBUG, NOTE, MPP_CLOCK_SYNC,MPP_CLOCK_DETAILED, WARNING use mpp_mod, only : mpp_pe, mpp_npes, mpp_root_pe, mpp_error, mpp_set_warn_level use mpp_mod, only : mpp_declare_pelist, mpp_set_current_pelist, mpp_sync @@ -50,7 +50,6 @@ module fv_mp_mod use mpp_domains_mod, only: nest_domain_type use mpp_parameter_mod, only : WUPDATE, EUPDATE, SUPDATE, NUPDATE, XUPDATE, YUPDATE use fv_arrays_mod, only: fv_atmos_type, fv_grid_bounds_type - use fms_io_mod, only: set_domain use mpp_mod, only : mpp_get_current_pelist, mpp_set_current_pelist use mpp_domains_mod, only : mpp_get_domain_shift use ensemble_manager_mod, only : get_ensemble_id @@ -911,9 +910,6 @@ subroutine switch_current_domain(new_domain,new_domain_for_coupler) ! if (debug .AND. (gid==masterproc)) write(*,200) tile, is, ie, js, je !200 format('New domain: ', i4.4, ' ', i4.4, ' ', i4.4, ' ', i4.4, ' ', i4.4, ' ') - call set_domain(new_domain) - - end subroutine switch_current_domain !depreciated diff --git a/tools/fv_nggps_diag.F90 b/tools/fv_nggps_diag.F90 index 27bd72dfa..8e68b2e44 100644 --- a/tools/fv_nggps_diag.F90 +++ b/tools/fv_nggps_diag.F90 @@ -45,10 +45,6 @@ module fv_nggps_diags_mod ! MODEL_ATMOS ! ! -! fms_io_mod -! set_domain, nullify_domain -! -! ! fv_arrays_mod ! fv_atmos_type ! @@ -68,7 +64,6 @@ module fv_nggps_diags_mod use mpp_mod, only: mpp_pe, mpp_root_pe,FATAL,mpp_error use constants_mod, only: grav, rdgas - use fms_io_mod, only: set_domain, nullify_domain use time_manager_mod, only: time_type, get_time use diag_manager_mod, only: register_diag_field, send_data use diag_axis_mod, only: get_axis_global_length, get_diag_axis, get_diag_axis_name @@ -165,8 +160,6 @@ subroutine fv_nggps_diag_init(Atm, axes, Time) jsdo = Atm(n)%bd%jsd; jedo = Atm(n)%bd%jed hydrostatico = Atm(n)%flagstruct%hydrostatic - call set_domain(Atm(1)%domain) ! Set domain so that diag_manager can access tile information - sphum = get_tracer_index (MODEL_ATMOS, 'sphum') liq_wat = get_tracer_index (MODEL_ATMOS, 'liq_wat') ice_wat = get_tracer_index (MODEL_ATMOS, 'ice_wat') @@ -729,8 +722,6 @@ subroutine fv_nggps_diag(Atm, zvir, Time) call store_data(id_uhmin25, uhmin25, Time, kstt_uhmin25, kend_uhmin25) endif - call nullify_domain() - end subroutine fv_nggps_diag subroutine fv_nggps_tavg(Atm, Time_step_atmos,avg_max_length,zvir) diff --git a/tools/fv_nudge.F90 b/tools/fv_nudge.F90 index 75a92bf5a..dd13dd974 100644 --- a/tools/fv_nudge.F90 +++ b/tools/fv_nudge.F90 @@ -29,10 +29,9 @@ module fv_nwp_nudge_mod use external_sst_mod, only: i_sst, j_sst, sst_ncep, sst_anom, forecast_mode use diag_manager_mod, only: register_diag_field, send_data use constants_mod, only: pi=>pi_8, grav, rdgas, cp_air, kappa, cnst_radius =>radius - use fms_mod, only: write_version_number, open_namelist_file, & - check_nml_error, file_exist, close_file -!use fms_io_mod, only: field_size - use mpp_mod, only: mpp_error, FATAL, stdlog, get_unit, mpp_pe + use fms_mod, only: write_version_number, check_nml_error + use fms2_io_mod, only: file_exists + use mpp_mod, only: mpp_error, FATAL, stdlog, get_unit, mpp_pe, input_nml_file use mpp_domains_mod, only: mpp_update_domains, domain2d use time_manager_mod, only: time_type, get_time, get_date @@ -1248,14 +1247,9 @@ subroutine fv_nwp_nudge_init(time, axes, npz, zvir, ak, bk, ts, phis, gridstruct track_file_name = "No_File_specified" - if( file_exist( 'input.nml' ) ) then - unit = open_namelist_file () - io = 1 - do while ( io .ne. 0 ) - read( unit, nml = fv_nwp_nudge_nml, iostat = io, end = 10 ) - ierr = check_nml_error(io,'fv_nwp_nudge_nml') - end do -10 call close_file ( unit ) + if( file_exists( 'input.nml' ) ) then + read( input_nml_file, nml = fv_nwp_nudge_nml, iostat = io ) + ierr = check_nml_error(io,'fv_nwp_nudge_nml') end if call write_version_number ( 'FV_NUDGE_MOD', version ) if ( master ) then @@ -1443,7 +1437,7 @@ subroutine get_ncep_analysis ( ps, u, v, t, q, zvir, ts, nfile, fname, bd ) #include - if( .not. file_exist(fname) ) then + if( .not. file_exists(fname) ) then call mpp_error(FATAL,'==> Error from get_ncep_analysis: file not found: '//fname) else call open_ncfile( fname, ncid ) ! open the file diff --git a/tools/fv_restart.F90 b/tools/fv_restart.F90 index 551194842..115a3762c 100644 --- a/tools/fv_restart.F90 +++ b/tools/fv_restart.F90 @@ -32,7 +32,7 @@ module fv_restart_mod use constants_mod, only: kappa, pi=>pi_8, omega, rdgas, grav, rvgas, cp_air, radius use fv_arrays_mod, only: fv_atmos_type, fv_nest_type, fv_grid_bounds_type, R_GRID use fv_io_mod, only: fv_io_init, fv_io_read_restart, fv_io_write_restart, & - remap_restart, fv_io_register_restart, fv_io_register_nudge_restart, & + remap_restart, fv_io_register_nudge_restart, & fv_io_register_restart_BCs, fv_io_write_BCs, fv_io_read_BCs use fv_grid_utils_mod, only: ptop_min, fill_ghost, g_sum, & make_eta_level, cubed_to_latlon, great_circle_dist @@ -58,8 +58,9 @@ module fv_restart_mod use mpp_mod, only: mpp_send, mpp_recv, mpp_sync_self, mpp_set_current_pelist, mpp_get_current_pelist, mpp_npes, mpp_pe, mpp_sync use mpp_domains_mod, only: CENTER, CORNER, NORTH, EAST, mpp_get_C2F_index, WEST, SOUTH use mpp_domains_mod, only: mpp_global_field - use fms_mod, only: file_exist use fv_treat_da_inc_mod, only: read_da_inc + use fms2_io_mod, only: file_exists, set_filename_appendix, FmsNetcdfFile_t, open_file, close_file + use fms_io_mod, only: fmsset_filename_appendix=> set_filename_appendix use coarse_grained_restart_files_mod, only: fv_io_write_restart_coarse implicit none @@ -118,12 +119,14 @@ subroutine fv_restart(fv_domain, Atm, dt_atmos, seconds, days, cold_start, grid_ character(len=128):: tname, errstring, fname, tracer_name character(len=120):: fname_ne, fname_sw character(len=3) :: gn + character(len=6) :: gnn integer :: npts, sphum integer, allocatable :: pelist(:), smoothed_topo(:) real :: sumpertn real :: zvir + type(FmsNetcdfFile_t) :: fileobj logical :: do_read_restart = .false. logical :: do_read_restart_bc = .false. integer, allocatable :: ideal_test_case(:), new_nest_topo(:) @@ -167,8 +170,8 @@ subroutine fv_restart(fv_domain, Atm, dt_atmos, seconds, days, cold_start, grid_ write(fname_ne,'(A, I2.2, A)') 'INPUT/fv_BC_ne.res.nest', Atm(n)%grid_number, '.nc' write(fname_sw,'(A, I2.2, A)') 'INPUT/fv_BC_sw.res.nest', Atm(n)%grid_number, '.nc' if (is_master()) print*, 'Searching for nested grid BC files ', trim(fname_ne), ' ', trim (fname_sw) - do_read_restart = file_exist(fname, Atm(n)%domain) - do_read_restart_bc = file_exist(fname_ne, Atm(n)%domain) .and. file_exist(fname_sw, Atm(n)%domain) + do_read_restart = file_exists(fname) + do_read_restart_bc = file_exists(fname_ne) .and. file_exists(fname_sw) if (is_master()) then print*, 'FV_RESTART: ', n, do_read_restart, do_read_restart_bc if (.not. do_read_restart_bc) write(*,*) 'BC files not found, re-generating nested grid boundary conditions' @@ -176,16 +179,20 @@ subroutine fv_restart(fv_domain, Atm, dt_atmos, seconds, days, cold_start, grid_ Atm(N)%neststruct%first_step = .not. do_read_restart_bc else fname='INPUT/fv_core.res.nc' - do_read_restart = file_exist('INPUT/fv_core.res.nc') .or. file_exist('INPUT/fv_core.res.tile1.nc') + do_read_restart = open_file(fileobj, fname, "read", is_restart=.true.) + if (do_read_restart) call close_file(fileobj) if (is_master()) print*, 'FV_RESTART: ', n, do_read_restart, do_read_restart_bc endif !2. Register restarts - !--- call fv_io_register_restart to register restart field to be written out in fv_io_write_restart - if ( n==this_grid ) call fv_io_register_restart(Atm(n)%domain,Atm(n:n)) - - !if (Atm(n)%neststruct%nested) call fv_io_register_restart_BCs(Atm(n)) !TODO put into fv_io_register_restart + !No longer need to register restarts in fv_restart_mod with fms2_io implementation + ! The two calls are needed until everything uses fms2io + if (Atm(n)%neststruct%nested .and. n==this_grid) then + write(gnn,'(A4, I2.2)') "nest", Atm(n)%grid_number + call set_filename_appendix(gnn) + call fmsset_filename_appendix(gnn) + endif !3preN. Topography BCs for nest, including setup for blending diff --git a/tools/fv_surf_map.F90 b/tools/fv_surf_map.F90 index 6fcc1b263..84252bec9 100644 --- a/tools/fv_surf_map.F90 +++ b/tools/fv_surf_map.F90 @@ -20,9 +20,9 @@ !*********************************************************************** module fv_surf_map_mod - use fms_mod, only: file_exist, check_nml_error, & - open_namelist_file, close_file, stdlog, & + use fms_mod, only: check_nml_error, stdlog, & mpp_pe, mpp_root_pe, FATAL, error_mesg + use fms2_io_mod, only: file_exists use mpp_mod, only: get_unit, input_nml_file, mpp_error use mpp_domains_mod, only: mpp_update_domains, domain2d use constants_mod, only: grav, radius, pi=>pi_8 @@ -175,7 +175,7 @@ subroutine surfdrv(npx, npy, grid, agrid, area, dx, dy, dxa, dya, dxc, dyc, sin_ ! ! surface file must be in NetCDF format ! - if ( file_exist(surf_file) ) then + if ( file_exists(surf_file) ) then if (surf_format == "netcdf") then status = nf_open (surf_file, NF_NOWRITE, ncid) @@ -1561,18 +1561,8 @@ subroutine read_namelist if (namelist_read) return -#ifdef INTERNAL_FILE_NML read (input_nml_file, nml=surf_map_nml, iostat=io) ierr = check_nml_error(io,'surf_map_nml') -#else - unit = open_namelist_file ( ) - ierr=1 - do while (ierr /= 0) - read (unit, nml=surf_map_nml, iostat=io, end=10) - ierr = check_nml_error(io,'surf_map_nml') - enddo - 10 call close_file (unit) -#endif ! write version and namelist to log file diff --git a/tools/fv_treat_da_inc.F90 b/tools/fv_treat_da_inc.F90 index e42878741..7165b4a84 100644 --- a/tools/fv_treat_da_inc.F90 +++ b/tools/fv_treat_da_inc.F90 @@ -35,8 +35,7 @@ module fv_treat_da_inc_mod - use fms_mod, only: file_exist, read_data, & - field_exist, write_version_number + use fms2_io_mod, only: file_exists use mpp_mod, only: mpp_error, FATAL, NOTE, mpp_pe use mpp_domains_mod, only: mpp_get_tile_id, & domain2d, & @@ -138,7 +137,7 @@ subroutine read_da_inc(Atm, fv_domain, bd, npz_in, nq, u, v, q, delp, pt, is_in, fname = 'INPUT/'//Atm%flagstruct%res_latlon_dynamics - if( file_exist(fname) ) then + if( file_exists(fname) ) then call open_ncfile( fname, ncid ) ! open the file call get_ncdim1( ncid, 'lon', tsize(1) ) call get_ncdim1( ncid, 'lat', tsize(2) ) diff --git a/tools/test_cases.F90 b/tools/test_cases.F90 index 55c8e9e87..e162f536b 100644 --- a/tools/test_cases.F90 +++ b/tools/test_cases.F90 @@ -6819,19 +6819,9 @@ subroutine read_namelist_test_case_nml(nml_filename) bubble_do = .false. test_case = 11 ! (USGS terrain) -#ifdef INTERNAL_FILE_NML ! Read Test_Case namelist read (input_nml_file,test_case_nml,iostat=ios) ierr = check_nml_error(ios,'test_case_nml') -#else - f_unit = open_namelist_file(nml_filename) - - ! Read Test_Case namelist - rewind (f_unit) - read (f_unit,test_case_nml,iostat=ios) - ierr = check_nml_error(ios,'test_case_nml') - call close_file(f_unit) -#endif write(unit, nml=test_case_nml)