Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

sync dev/emc #16

Merged
merged 6 commits into from
Aug 15, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
29 changes: 16 additions & 13 deletions driver/fvGFS/atmosphere.F90
Original file line number Diff line number Diff line change
Expand Up @@ -324,10 +324,12 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area)
logical :: dycore_only = .false.
logical :: debug = .false.
logical :: sync = .false.
logical :: ignore_rst_cksum = .false.
integer, parameter :: maxhr = 4096
real, dimension(maxhr) :: fdiag = 0.
real :: fhmax=384.0, fhmaxhf=120.0, fhout=3.0, fhouthf=1.0,avg_max_length=3600.
namelist /atmos_model_nml/ blocksize, chksum_debug, dycore_only, debug, sync, fdiag, fhmax, fhmaxhf, fhout, fhouthf, ccpp_suite, avg_max_length
namelist /atmos_model_nml/ blocksize, chksum_debug, dycore_only, debug, sync, fdiag, fhmax, fhmaxhf, fhout, fhouthf, ccpp_suite, &
& avg_max_length, ignore_rst_cksum
! *DH 20210326

!For regional
Expand Down Expand Up @@ -357,7 +359,7 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area)
#ifdef MOVING_NEST
call fv_tracker_init(size(Atm))
if (mygrid .eq. 2) call allocate_tracker(mygrid, Atm(mygrid)%bd%isc, Atm(mygrid)%bd%iec, Atm(mygrid)%bd%jsc, Atm(mygrid)%bd%jec)
#endif
#endif

Atm(mygrid)%Time_init = Time_init

Expand Down Expand Up @@ -450,6 +452,15 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area)
!--- allocate pref
allocate(pref(npz+1,2), dum1d(npz+1))

! DH* 20210326
! First, read atmos_model_nml namelist section - this is a workaround to avoid
! unnecessary additional changes to the input namelists, in anticipation of the
! implementation of a generic interface for GFDL and CCPP fast physics soon
read(input_nml_file, nml=atmos_model_nml, iostat=io)
ierr = check_nml_error(io, 'atmos_model_nml')
!write(0,'(a)') "It's me, and my physics suite is '" // trim(ccpp_suite) // "'"
! *DH 20210326

call fv_restart(Atm(mygrid)%domain, Atm, seconds, days, cold_start, Atm(mygrid)%gridstruct%grid_type, mygrid)

fv_time = Time
Expand Down Expand Up @@ -491,15 +502,6 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area)

! Do CCPP fast physics initialization before call to adiabatic_init (since this calls fv_dynamics)

! DH* 20210326
! First, read atmos_model_nml namelist section - this is a workaround to avoid
! unnecessary additional changes to the input namelists, in anticipation of the
! implementation of a generic interface for GFDL and CCPP fast physics soon
read(input_nml_file, nml=atmos_model_nml, iostat=io)
ierr = check_nml_error(io, 'atmos_model_nml')
!write(0,'(a)') "It's me, and my physics suite is '" // trim(ccpp_suite) // "'"
! *DH 20210326

! For fast physics running over the entire domain, block
! and thread number are not used; set to safe values
cdata%blk_no = 1
Expand Down Expand Up @@ -954,10 +956,10 @@ end subroutine get_nth_domain_info
!! the "domain2d" variable associated with the coupling grid and the
!! decomposition for the current cubed-sphere tile.
!>@detail Coupling is done using the mass/temperature grid with no halos.
subroutine atmosphere_domain ( fv_domain, layout, regional, nested, &
subroutine atmosphere_domain ( fv_domain, rd_domain, layout, regional, nested, &
moving_nest_parent, is_moving_nest, &
ngrids_atmos, mygrid_atmos, pelist )
type(domain2d), intent(out) :: fv_domain
type(domain2d), intent(out) :: fv_domain, rd_domain
integer, intent(out) :: layout(2)
logical, intent(out) :: regional
logical, intent(out) :: nested
Expand All @@ -970,6 +972,7 @@ subroutine atmosphere_domain ( fv_domain, layout, regional, nested, &
integer :: n

fv_domain = Atm(mygrid)%domain_for_coupler
rd_domain = Atm(mygrid)%domain_for_read
layout(1:2) = Atm(mygrid)%layout(1:2)
regional = Atm(mygrid)%flagstruct%regional
nested = ngrids > 1
Expand Down
2 changes: 2 additions & 0 deletions model/fv_arrays.F90
Original file line number Diff line number Diff line change
Expand Up @@ -516,6 +516,7 @@ module fv_arrays_mod
!-----------------------------------------------------------------------------------------------

logical :: reset_eta = .false.
logical :: ignore_rst_cksum = .false. !< enfore (.false.) or override (.true.) data integrity restart checksums
real :: p_fac = 0.05 !< Safety factor for minimum nonhydrostatic pressures, which
!< will be limited so the full pressure is no less than p_fac
!< times the hydrostatic pressure. This is only of concern in mid-top
Expand Down Expand Up @@ -1303,6 +1304,7 @@ module fv_arrays_mod
#if defined(SPMD)

type(domain2D) :: domain_for_coupler !< domain used in coupled model with halo = 1.
type(domain2D) :: domain_for_read !< domain used for reads to increase performance when io_layout=(1,1)

!global tile and tile_of_mosaic only have a meaning for the CURRENT pe
integer :: num_contact, npes_per_tile, global_tile, tile_of_mosaic, npes_this_grid
Expand Down
9 changes: 6 additions & 3 deletions model/fv_control.F90
Original file line number Diff line number Diff line change
Expand Up @@ -198,7 +198,7 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split,
real, intent(in) :: dt_atmos
integer, intent(OUT) :: this_grid
logical, allocatable, intent(OUT) :: grids_on_this_pe(:)
character(len=32), optional, intent(in) :: nml_filename_in ! alternate nml
character(len=32), optional, intent(in) :: nml_filename_in ! alternate nml
logical, optional, intent(in) :: skip_nml_read_in ! use previously loaded nml

integer, intent(INOUT) :: p_split
Expand Down Expand Up @@ -291,6 +291,7 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split,
real(kind=R_GRID) , pointer :: target_lon

logical , pointer :: reset_eta
logical , pointer :: ignore_rst_cksum
real , pointer :: p_fac
real , pointer :: a_imp
integer , pointer :: n_split
Expand Down Expand Up @@ -676,7 +677,8 @@ subroutine fv_control_init(Atm, dt_atmos, this_grid, grids_on_this_pe, p_split,
Atm(this_grid)%flagstruct%grid_type,Atm(this_grid)%neststruct%nested, &
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)
Atm(this_grid)%domain_for_coupler,Atm(this_grid)%domain_for_read,Atm(this_grid)%num_contact, &
Atm(this_grid)%pelist)
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)
Expand Down Expand Up @@ -857,6 +859,7 @@ subroutine set_namelist_pointers(Atm)
regional_bcs_from_gsi => Atm%flagstruct%regional_bcs_from_gsi
write_restart_with_bcs => Atm%flagstruct%write_restart_with_bcs
reset_eta => Atm%flagstruct%reset_eta
ignore_rst_cksum => Atm%flagstruct%ignore_rst_cksum
p_fac => Atm%flagstruct%p_fac
a_imp => Atm%flagstruct%a_imp
n_split => Atm%flagstruct%n_split
Expand Down Expand Up @@ -1075,7 +1078,7 @@ subroutine read_namelist_fv_core_nml(Atm)
write_coarse_restart_files,&
write_coarse_diagnostics,&
write_only_coarse_intermediate_restarts, &
write_coarse_agrid_vel_rst, write_coarse_dgrid_vel_rst
write_coarse_agrid_vel_rst, write_coarse_dgrid_vel_rst, ignore_rst_cksum


! Read FVCORE namelist
Expand Down
2 changes: 1 addition & 1 deletion model/fv_fill.F90
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ module fv_fill_mod
! </table>

use mpp_domains_mod, only: mpp_update_domains, domain2D
use platform_mod, only: kind_phys => r8_kind
use GFS_typedefs, only: kind_phys

implicit none
public fillz
Expand Down
59 changes: 26 additions & 33 deletions tools/external_ic.F90
Original file line number Diff line number Diff line change
Expand Up @@ -206,10 +206,9 @@ module external_ic_mod

contains

subroutine get_external_ic( Atm, fv_domain, cold_start )
subroutine get_external_ic( Atm, cold_start )

type(fv_atmos_type), intent(inout), target :: Atm
type(domain2d), intent(inout) :: fv_domain
logical, intent(IN) :: cold_start
real:: alpha = 0.
real rdg
Expand Down Expand Up @@ -260,14 +259,14 @@ subroutine get_external_ic( Atm, fv_domain, cold_start )
enddo
enddo

call mpp_update_domains( f0, fv_domain )
call mpp_update_domains( f0, Atm%domain )
if ( Atm%gridstruct%cubed_sphere .and. (.not. Atm%gridstruct%bounded_domain))then
call fill_corners(f0, Atm%npx, Atm%npy, YDir)
endif

! Read in cubed_sphere terrain
if ( Atm%flagstruct%mountain ) then
call get_cubed_sphere_terrain(Atm, fv_domain)
call get_cubed_sphere_terrain(Atm)
else
if (.not. Atm%neststruct%nested) Atm%phis = 0. !TODO: Not sure about this line --- lmh 30 may 18
endif
Expand All @@ -276,32 +275,32 @@ subroutine get_external_ic( Atm, fv_domain, cold_start )
if ( Atm%flagstruct%ncep_ic ) then
nq = 1
call timing_on('NCEP_IC')
call get_ncep_ic( Atm, fv_domain, nq )
call get_ncep_ic( Atm, nq )
call timing_off('NCEP_IC')
#ifdef FV_TRACERS
if (.not. cold_start) then
call fv_io_read_tracers( fv_domain, Atm )
call fv_io_read_tracers( Atm )
if(is_master()) write(*,*) 'All tracers except sphum replaced by FV IC'
endif
#endif
elseif ( Atm%flagstruct%nggps_ic ) then
call timing_on('NGGPS_IC')
call get_nggps_ic( Atm, fv_domain )
call get_nggps_ic( Atm )
call timing_off('NGGPS_IC')
elseif ( Atm%flagstruct%hrrrv3_ic ) then
call timing_on('HRRR_IC')
call get_hrrr_ic( Atm, fv_domain )
call get_hrrr_ic( Atm )
call timing_off('HRRR_IC')
elseif ( Atm%flagstruct%ecmwf_ic ) then
if( is_master() ) write(*,*) 'Calling get_ecmwf_ic'
call timing_on('ECMWF_IC')
call get_ecmwf_ic( Atm, fv_domain )
call get_ecmwf_ic( Atm )
call timing_off('ECMWF_IC')
else
! The following is to read in legacy lat-lon FV core restart file
! is Atm%q defined in all cases?
nq = size(Atm%q,4)
call get_fv_ic( Atm, fv_domain, nq )
call get_fv_ic( Atm, nq )
endif

call prt_maxmin('PS', Atm%ps, is, ie, js, je, ng, 1, 0.01)
Expand Down Expand Up @@ -368,9 +367,8 @@ end subroutine get_external_ic


!------------------------------------------------------------------
subroutine get_cubed_sphere_terrain( Atm, fv_domain )
subroutine get_cubed_sphere_terrain( Atm )
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
Expand All @@ -393,13 +391,13 @@ subroutine get_cubed_sphere_terrain( Atm, fv_domain )
jed = Atm%bd%jed
ng = Atm%bd%ng

tile_id = mpp_get_tile_id( fv_domain )
tile_id = mpp_get_tile_id( Atm%domain )

fname = 'INPUT/fv_core.res.nc'
call mpp_error(NOTE, 'external_ic: looking for '//fname)


if( open_file(Fv_core, fname, "read", fv_domain, is_restart=.true.) ) then
if( open_file(Fv_core, fname, "read", Atm%domain_for_read, is_restart=.true.) ) then
call read_data(Fv_core, 'phis', Atm%phis(is:ie,js:je))
call close_file(Fv_core)
else
Expand Down Expand Up @@ -428,7 +426,7 @@ end subroutine get_cubed_sphere_terrain
!>@brief The subroutine 'get_nggps_ic' reads in data after it has been preprocessed with
!! NCEP/EMC orography maker and 'global_chgres', and has been horiztontally
!! interpolated to the current cubed-sphere grid
subroutine get_nggps_ic (Atm, fv_domain)
subroutine get_nggps_ic (Atm)

!>variables read in from 'gfs_ctrl.nc'
!> VCOORD - level information
Expand Down Expand Up @@ -456,7 +454,6 @@ subroutine get_nggps_ic (Atm, fv_domain)
#endif

type(fv_atmos_type), intent(inout) :: Atm
type(domain2d), intent(inout) :: fv_domain
! local:
real, dimension(:), allocatable:: ak, bk
real, dimension(:,:), allocatable:: wk2, ps, oro_g
Expand Down Expand Up @@ -585,7 +582,7 @@ subroutine get_nggps_ic (Atm, fv_domain)

!--- read in surface temperature (k) and land-frac
! surface skin temperature
if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain_for_read, is_restart=.true., dont_add_res_to_filename=.true.) ) then
naxis_dims = get_variable_num_dimensions(SFC_restart, 'tsea')
allocate (dim_names_alloc(naxis_dims))
call get_variable_dimension_names(SFC_restart, 'tsea', dim_names_alloc)
Expand All @@ -605,7 +602,7 @@ subroutine get_nggps_ic (Atm, fv_domain)
dim_names_2d(2) = "lon"

! terrain surface height -- (needs to be transformed into phis = zs*grav)
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain_for_read, 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
Expand Down Expand Up @@ -915,7 +912,7 @@ subroutine read_gfs_ic()
dim_names_3d4(1) = "levp"

! surface pressure (Pa)
if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain_for_read, 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)
Expand Down Expand Up @@ -945,7 +942,7 @@ subroutine read_gfs_ic()
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
Expand All @@ -964,7 +961,7 @@ subroutine read_gfs_ic()
end subroutine get_nggps_ic
!------------------------------------------------------------------
!------------------------------------------------------------------
subroutine get_hrrr_ic (Atm, fv_domain)
subroutine get_hrrr_ic (Atm)
! read in data after it has been preprocessed with
! NCEP/EMC orography maker
!
Expand All @@ -990,7 +987,6 @@ subroutine get_hrrr_ic (Atm, fv_domain)


type(fv_atmos_type), intent(inout) :: Atm
type(domain2d), intent(inout) :: fv_domain
! local:
real, dimension(:), allocatable:: ak, bk
real, dimension(:,:), allocatable:: wk2, ps, oro_g
Expand Down Expand Up @@ -1104,7 +1100,7 @@ subroutine get_hrrr_ic (Atm, fv_domain)

!--- read in surface temperature (k) and land-frac
! surface skin temperature
if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(SFC_restart, fn_sfc_ics, "read", Atm%domain_for_read, is_restart=.true., dont_add_res_to_filename=.true.) ) then
call get_variable_dimension_names(SFC_restart, 'tsea', dim_names_2d)
call register_axis(SFC_restart, dim_names_2d(2), "y")
call register_axis(SFC_restart, dim_names_2d(1), "x")
Expand All @@ -1121,7 +1117,7 @@ subroutine get_hrrr_ic (Atm, fv_domain)
dim_names_2d(2) = "lon"

! terrain surface height -- (needs to be transformed into phis = zs*grav)
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain_for_read, 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
Expand Down Expand Up @@ -1164,7 +1160,7 @@ subroutine get_hrrr_ic (Atm, fv_domain)
dim_names_3d4(1) = "levp"

! edge pressure (Pa)
if( open_file(HRRR_restart, fn_hrr_ics, "read", Atm%domain,is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(HRRR_restart, fn_hrr_ics, "read", Atm%domain_for_read, 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)
Expand Down Expand Up @@ -1360,9 +1356,8 @@ end subroutine get_hrrr_ic
!------------------------------------------------------------------
!------------------------------------------------------------------
!>@brief The subroutine 'get_ncep_ic' reads in the specified NCEP analysis or reanalysis dataset
subroutine get_ncep_ic( Atm, fv_domain, nq )
subroutine get_ncep_ic( Atm, nq )
type(fv_atmos_type), intent(inout) :: Atm
type(domain2d), intent(inout) :: fv_domain
integer, intent(in):: nq
! local:
#ifdef HIWPP_ETA
Expand Down Expand Up @@ -1818,7 +1813,7 @@ end subroutine get_ncep_ic
!>@brief The subroutine 'get_ecmwf_ic' reads in initial conditions from ECMWF analyses
!! (EXPERIMENTAL: contact Jan-Huey Chen jan-huey.chen@noaa.gov for support)
!>@authors Jan-Huey Chen, Xi Chen, Shian-Jiann Lin
subroutine get_ecmwf_ic( Atm, fv_domain )
subroutine get_ecmwf_ic( Atm )

#ifdef __PGI
use GFS_restart, only : GFS_restart_type
Expand All @@ -1827,7 +1822,6 @@ subroutine get_ecmwf_ic( Atm, fv_domain )
#endif

type(fv_atmos_type), intent(inout) :: Atm
type(domain2d), intent(inout) :: fv_domain
! local:
real :: ak_ec(138), bk_ec(138)
data ak_ec/ 0.000000, 2.000365, 3.102241, 4.666084, 6.827977, 9.746966, &
Expand Down Expand Up @@ -2062,7 +2056,7 @@ subroutine get_ecmwf_ic( Atm, fv_domain )
dim_names_3d4(1) = "levp"

!! Read in model terrain from oro_data.tile?.nc
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(ORO_restart, fn_oro_ics, "read", Atm%domain_for_read, 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
Expand All @@ -2088,7 +2082,7 @@ subroutine get_ecmwf_ic( Atm, fv_domain )
allocate (ps_gfs(is:ie,js:je))
allocate (zh_gfs(is:ie,js:je,levp_gfs+1))

if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain, is_restart=.true., dont_add_res_to_filename=.true.) ) then
if( open_file(GFS_restart, fn_gfs_ics, "read", Atm%domain_for_read, 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, "levp", size(zh_gfs,3))
Expand Down Expand Up @@ -2622,9 +2616,8 @@ subroutine get_ecmwf_ic( Atm, fv_domain )
end subroutine get_ecmwf_ic
!------------------------------------------------------------------
!------------------------------------------------------------------
subroutine get_fv_ic( Atm, fv_domain, nq )
subroutine get_fv_ic( Atm, nq )
type(fv_atmos_type), intent(inout) :: Atm
type(domain2d), intent(inout) :: fv_domain
integer, intent(in):: nq

type(FmsNetcdfFile_t) :: Latlon_dyn, Latlon_tra
Expand Down
Loading