diff --git a/.gitignore b/.gitignore index ac5fcf245..ce022e6e3 100644 --- a/.gitignore +++ b/.gitignore @@ -2,3 +2,4 @@ *.mod *.a *.pyc +*.i90 diff --git a/.gitmodules b/.gitmodules index 4760351ce..22e3da842 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,7 +1,7 @@ [submodule "atmos_cubed_sphere"] path = atmos_cubed_sphere - url = https://github.com/NOAA-EMC/GFDL_atmos_cubed_sphere - branch = dev/emc + url = https://github.com/NOAA-GSL/GFDL_atmos_cubed_sphere + branch = RRFS_cloud [submodule "ccpp/framework"] path = ccpp/framework url = https://github.com/NOAA-GSL/ccpp-framework @@ -9,4 +9,4 @@ [submodule "ccpp/physics"] path = ccpp/physics url = https://github.com/NOAA-GSL/ccpp-physics - branch = gsl/develop + branch = RRFS_cloud diff --git a/CMakeLists.txt b/CMakeLists.txt index df1fc2344..4b0947f55 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -165,6 +165,7 @@ add_dependencies(fv3atm ccppdriver ccppdata ccppphys ccpp) target_link_libraries(fv3atm PUBLIC ccppdriver ccppdata ccppphys ccpp) target_include_directories(fv3atm PRIVATE ${CMAKE_CURRENT_BINARY_DIR}/stochastic_physics) +target_include_directories(fv3atm PRIVATE ${CMAKE_CURRENT_BINARY_DIR}/../stochastic_physics) target_compile_definitions(fv3atm PRIVATE "${_fv3atm_defs_private}") target_link_libraries(fv3atm PUBLIC fv3dycore @@ -192,7 +193,7 @@ endif() ### Install ############################################################################### install( - TARGETS fv3atm fv3dycore io ${CCPP_LIBRARIES} cpl stochastic_physics stochastic_physics_wrapper + TARGETS fv3atm fv3dycore io ${CCPP_LIBRARIES} cpl stochastic_physics_wrapper EXPORT fv3atm-config LIBRARY DESTINATION lib ARCHIVE DESTINATION lib) diff --git a/atmos_cubed_sphere b/atmos_cubed_sphere index 306ff3137..72656741e 160000 --- a/atmos_cubed_sphere +++ b/atmos_cubed_sphere @@ -1 +1 @@ -Subproject commit 306ff31371e74694e5d9f4a57584295c7122b9ac +Subproject commit 72656741e386f5813ee7e908ee11590453aecc13 diff --git a/atmos_model.F90 b/atmos_model.F90 index 3730da692..04ecae56b 100644 --- a/atmos_model.F90 +++ b/atmos_model.F90 @@ -869,6 +869,7 @@ end subroutine update_atmos_model_state ! subroutine atmos_model_end (Atmos) + use get_stochy_pattern_mod, only: write_stoch_restart_atm type (atmos_data_type), intent(inout) :: Atmos !---local variables integer :: idx, seconds, ierr @@ -878,12 +879,12 @@ subroutine atmos_model_end (Atmos) call atmosphere_end (Atmos % Time, Atmos%grid, restart_endfcst) - call stochastic_physics_wrapper_end(GFS_control) - if(restart_endfcst) then call FV3GFS_restart_write (GFS_data, GFS_restart_var, Atm_block, & GFS_control, Atmos%domain) + call write_stoch_restart_atm('RESTART/atm_stoch.res.nc') endif + call stochastic_physics_wrapper_end(GFS_control) ! Fast physics (from dynamics) are finalized in atmosphere_end above; ! standard/slow physics (from IPD) are finalized in CCPP_step 'finalize'. diff --git a/ccpp/CMakeLists.txt b/ccpp/CMakeLists.txt index 2f4ab9da4..04cb44ddf 100644 --- a/ccpp/CMakeLists.txt +++ b/ccpp/CMakeLists.txt @@ -220,18 +220,6 @@ else (SIONLIB) message (STATUS "Disable SIONlib support") endif (SIONLIB) -#------------------------------------------------------------------------------ -# Set Intel MKL flags for preprocessor, compiler and linker (if defined) -if(MKL_DIR) - set (MKL_INC "-m64 -I${MKL_DIR}/include") - set (MKL_LIB "-L${MKL_DIR}/lib -Wl,-rpath,${MKL_DIR}/lib -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl") - set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} ${MKL_INC} ${MKL_LIB}") - ADD_DEFINITIONS(-DMKL) - message (STATUS "Enable Intel MKL support") -else(MKL_DIR) - message (STATUS "Disable Intel MKL support") -endif(MKL_DIR) - #------------------------------------------------------------------------------ # Set netCDF flags for preprocessor, compiler and linker (if defined) # Legacy settings for old make build diff --git a/ccpp/config/ccpp_prebuild_config.py b/ccpp/config/ccpp_prebuild_config.py index efc4cb6ec..200af941e 100755 --- a/ccpp/config/ccpp_prebuild_config.py +++ b/ccpp/config/ccpp_prebuild_config.py @@ -159,6 +159,7 @@ 'ccpp/physics/physics/mp_thompson_pre.F90', 'ccpp/physics/physics/mp_thompson.F90', 'ccpp/physics/physics/mp_thompson_post.F90', + 'ccpp/physics/physics/mp_nsslg.F90', 'ccpp/physics/physics/ozphys.f', 'ccpp/physics/physics/ozphys_2015.f', 'ccpp/physics/physics/precpd.f', @@ -301,6 +302,13 @@ 'rime_factor', ], }, + 'mp_nsslg' : { + 'mp_nsslg_run' : [ + 'effective_radius_of_stratiform_cloud_liquid_water_particle_in_um', + 'effective_radius_of_stratiform_cloud_ice_particle_in_um', + 'effective_radius_of_stratiform_cloud_snow_particle_in_um', + ], + }, 'rrtmgp_sw_rte' : { 'rrtmgp_sw_rte_run' : [ 'components_of_surface_downward_shortwave_fluxes', diff --git a/ccpp/data/GFS_typedefs.F90 b/ccpp/data/GFS_typedefs.F90 index 72e7ca565..d0882fad2 100644 --- a/ccpp/data/GFS_typedefs.F90 +++ b/ccpp/data/GFS_typedefs.F90 @@ -714,6 +714,8 @@ module GFS_typedefs integer :: idcor_con = 0 !< choice for decorrelation-length: Use constant value integer :: idcor_hogan = 1 !< choice for decorrelation-length: (https://rmets.onlinelibrary.wiley.com/doi/full/10.1002/qj.647) integer :: idcor_oreopoulos = 2 !< choice for decorrelation-length: (10.5194/acp-12-9097-2012) + integer :: imp_physics_nssl2m = 17 !< choice of NSSL microphysics scheme with background CCN + integer :: imp_physics_nssl2mccn = 18 !< choice of NSSL microphysics scheme with predicted CCN !--- Z-C microphysical parameters real(kind=kind_phys) :: psautco(2) !< [in] auto conversion coeff from ice to snow real(kind=kind_phys) :: prautco(2) !< [in] auto conversion coeff from cloud to rain @@ -758,6 +760,13 @@ module GFS_typedefs real(kind=kind_phys) :: shoc_parm(5) !< critical pressure in Pa for tke dissipation in shoc integer :: ncnd !< number of cloud condensate types + !--- NSSL microphysics params + real(kind=kind_phys) :: nssl_cccn !< CCN concentration (m-3) + real(kind=kind_phys) :: nssl_alphah !< graupel shape parameter + real(kind=kind_phys) :: nssl_alphahl !< hail shape parameter + logical :: nssl_hail_on !< NSSL flag to deactivate the hail category + logical :: nssl_invertccn !< NSSL flag to treat CCN as activated (true) or unactivated (false) + !--- Thompson's microphysical parameters logical :: ltaerosol !< flag for aerosol version logical :: lradar !< flag for radar reflectivity @@ -1070,11 +1079,15 @@ module GFS_typedefs !--- stochastic physics control parameters logical :: do_sppt + logical :: pert_clds + logical :: pert_radtend + logical :: pert_mp logical :: use_zmtnblck logical :: do_shum logical :: do_skeb integer :: skeb_npass integer :: lndp_type + real(kind=kind_phys) :: sppt_amp ! pjp cloud perturbations integer :: n_var_lndp logical :: lndp_each_step ! flag to indicate that land perturbations are applied at every time step, ! otherwise they are applied only after gcycle is run @@ -1096,12 +1109,18 @@ module GFS_typedefs integer :: ntrw !< tracer index for rain water integer :: ntsw !< tracer index for snow water integer :: ntgl !< tracer index for graupel + integer :: nthl !< tracer index for hail integer :: ntclamt !< tracer index for cloud amount integer :: ntlnc !< tracer index for liquid number concentration integer :: ntinc !< tracer index for ice number concentration integer :: ntrnc !< tracer index for rain number concentration integer :: ntsnc !< tracer index for snow number concentration integer :: ntgnc !< tracer index for graupel number concentration + integer :: nthnc !< tracer index for hail number concentration + integer :: ntccn !< tracer index for CCN + integer :: ntccna !< tracer index for activated CCN + integer :: ntgv !< tracer index for graupel particle volume + integer :: nthv !< tracer index for hail particle volume integer :: ntke !< tracer index for kinetic energy integer :: nto !< tracer index for oxygen ion integer :: nto2 !< tracer index for oxygen @@ -1287,7 +1306,7 @@ module GFS_typedefs real (kind=kind_phys), pointer :: acvt (:) => null() !< arrays used by cnvc90 top (cnvc90.f) !--- Stochastic physics properties calculated in physics_driver - real (kind=kind_phys), pointer :: dtdtr (:,:) => null() !< temperature change due to radiative heating per time step (K) + real (kind=kind_phys), pointer :: dtdtnp (:,:) => null() !< temperature change from physics that should not be perturbed with SPPT (k) real (kind=kind_phys), pointer :: dtotprcp (:) => null() !< change in totprcp (diag_type) real (kind=kind_phys), pointer :: dcnvprcp (:) => null() !< change in cnvprcp (diag_type) real (kind=kind_phys), pointer :: drain_cpl (:) => null() !< change in rain_cpl (coupling_type) @@ -1536,6 +1555,7 @@ module GFS_typedefs real (kind=kind_phys), pointer :: skebv_wts(:,:) => null() !< 10 meter v wind speed real (kind=kind_phys), pointer :: sppt_wts(:,:) => null() !< real (kind=kind_phys), pointer :: shum_wts(:,:) => null() !< + real (kind=kind_phys), pointer :: sfc_wts(:,:) => null() !< real (kind=kind_phys), pointer :: zmtnblck(:) => null() ! null() !< u momentum change due to physics real (kind=kind_phys), pointer :: dv3dt (:,:,:) => null() !< v momentum change due to physics @@ -1757,7 +1777,6 @@ module GFS_typedefs real (kind=kind_phys), pointer :: drain(:) => null() !< real (kind=kind_phys), pointer :: drain_in_m_sm1(:) => null() !< real (kind=kind_phys), pointer :: dtdt(:,:) => null() !< - real (kind=kind_phys), pointer :: dtdtc(:,:) => null() !< real (kind=kind_phys), pointer :: dtsfc1(:) => null() !< real (kind=kind_phys), pointer :: dtzm(:) => null() !< real (kind=kind_phys), pointer :: dt_mf(:,:) => null() !< @@ -1882,6 +1901,7 @@ module GFS_typedefs real (kind=kind_phys), pointer :: oc(:) => null() !< real (kind=kind_phys), pointer :: olyr(:,:) => null() !< logical , pointer :: otspt(:,:) => null() !< + logical , pointer :: otsptflag(:) => null() !< integer :: oz_coeffp5 !< logical :: phys_hydrostatic !< real (kind=kind_phys), pointer :: plvl(:,:) => null() !< @@ -3084,6 +3104,12 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & logical :: mg_do_ice_gmao = .false. !< set .true. to turn on gmao ice formulation logical :: mg_do_liq_liu = .true. !< set .true. to turn on liu liquid treatment + !--- NSSL microphysics params + real(kind=kind_phys) :: nssl_cccn = 0.6e9 !< CCN concentration (m-3) + real(kind=kind_phys) :: nssl_alphah = 0.0 !< graupel shape parameter + real(kind=kind_phys) :: nssl_alphahl = 1.0 !< hail shape parameter + logical :: nssl_hail_on = .true. !< NSSL flag to deactivate the hail category + logical :: nssl_invertccn = .false. !< NSSL flag to treat CCN as activated (true) or unactivated (false) !--- Thompson microphysical parameters logical :: ltaerosol = .false. !< flag for aerosol version @@ -3385,6 +3411,9 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & !--- stochastic physics control parameters logical :: do_sppt = .false. + logical :: pert_mp = .false. + logical :: pert_clds = .false. + logical :: pert_radtend = .true. logical :: use_zmtnblck = .false. logical :: do_shum = .false. logical :: do_skeb = .false. @@ -3428,7 +3457,8 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & mg_ncnst, mg_ninst, mg_ngnst, sed_supersat, do_sb_physics, & mg_alf, mg_qcmin, mg_do_ice_gmao, mg_do_liq_liu, & ltaerosol, lradar, nsradar_reset, lrefres, ttendlim, & - lgfdlmprad, & + lgfdlmprad, nssl_cccn, nssl_alphah, nssl_alphahl, & + nssl_invertccn, & !--- max hourly avg_max_length, & !--- land/surface model control @@ -3467,6 +3497,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & dlqf, rbcr, shoc_parm, psauras, prauras, wminras, & do_sppt, do_shum, do_skeb, & lndp_type, n_var_lndp, lndp_each_step, & + pert_mp,pert_clds,pert_radtend, & !--- Rayleigh friction prslrd0, ral_ts, ldiag_ugwp, do_ugwp, do_tofd, & ! --- Ferrier-Aligo @@ -3836,6 +3867,13 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & Model%tcr = tcr Model%tcrf = 1.0/(tcr-tf) +!-- NSSL microphysics params + Model%nssl_cccn = nssl_cccn + Model%nssl_alphah = nssl_alphah + Model%nssl_alphahl = nssl_alphahl + Model%nssl_hail_on = nssl_hail_on + Model%nssl_invertccn = nssl_invertccn + !--- Thompson MP parameters Model%ltaerosol = ltaerosol Model%lradar = lradar @@ -4147,6 +4185,9 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & ! physics that are parsed here and then compared in init_stochastic_physics ! to the stochastic physics namelist parametersto ensure consistency. Model%do_sppt = do_sppt + Model%pert_mp = pert_mp + Model%pert_clds = pert_clds + Model%pert_radtend = pert_radtend Model%use_zmtnblck = use_zmtnblck Model%do_shum = do_shum Model%do_skeb = do_skeb @@ -4208,12 +4249,18 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & Model%ntrw = get_tracer_index(Model%tracer_names, 'rainwat', Model%me, Model%master, Model%debug) Model%ntsw = get_tracer_index(Model%tracer_names, 'snowwat', Model%me, Model%master, Model%debug) Model%ntgl = get_tracer_index(Model%tracer_names, 'graupel', Model%me, Model%master, Model%debug) + Model%nthl = get_tracer_index(Model%tracer_names, 'hailwat', Model%me, Model%master, Model%debug) Model%ntclamt = get_tracer_index(Model%tracer_names, 'cld_amt', Model%me, Model%master, Model%debug) Model%ntlnc = get_tracer_index(Model%tracer_names, 'water_nc', Model%me, Model%master, Model%debug) Model%ntinc = get_tracer_index(Model%tracer_names, 'ice_nc', Model%me, Model%master, Model%debug) Model%ntrnc = get_tracer_index(Model%tracer_names, 'rain_nc', Model%me, Model%master, Model%debug) Model%ntsnc = get_tracer_index(Model%tracer_names, 'snow_nc', Model%me, Model%master, Model%debug) Model%ntgnc = get_tracer_index(Model%tracer_names, 'graupel_nc', Model%me, Model%master, Model%debug) + Model%nthnc = get_tracer_index(Model%tracer_names, 'hail_nc', Model%me, Model%master, Model%debug) + Model%ntccn = get_tracer_index(Model%tracer_names, 'ccn_nc', Model%me, Model%master, Model%debug) + Model%ntccna = get_tracer_index(Model%tracer_names, 'ccna_nc', Model%me, Model%master, Model%debug) + Model%ntgv = get_tracer_index(Model%tracer_names, 'graupel_vol',Model%me, Model%master, Model%debug) + Model%nthv = get_tracer_index(Model%tracer_names, 'hail_vol', Model%me, Model%master, Model%debug) Model%ntke = get_tracer_index(Model%tracer_names, 'sgs_tke', Model%me, Model%master, Model%debug) Model%nqrimef = get_tracer_index(Model%tracer_names, 'q_rimef', Model%me, Model%master, Model%debug) Model%ntwa = get_tracer_index(Model%tracer_names, 'liq_aero', Model%me, Model%master, Model%debug) @@ -4237,8 +4284,78 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & endif endif + + IF ( Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn ) THEN !{ + ! turn off tracer for CCCN if physics = 17 + IF ( Model%imp_physics == Model%imp_physics_nssl2m ) THEN + if (Model%me == Model%master) write(*,*) 'NSSL micro: CCN is OFF' + IF ( Model%ntccn > 1 ) THEN + IF (Model%me == Model%master) then + write(*,*) 'NSSL micro: error! CCN is OFF but ntccn > 1. Should remove ccn_nc from field_table' + write(0,*) 'NSSL micro: error! CCN is OFF but ntccn > 1. Should remove ccn_nc from field_table' + ENDIF + ENDIF + Model%ntccn = -99 + Model%ntccna = -99 + ELSEIF ( Model%ntccn < 1 ) THEN + if (Model%me == Model%master) then + write(*,*) 'NSSL micro: error! CCN is ON but ntccn < 1. Must set ccn_nc in field_table' + write(0,*) 'NSSL micro: error! CCN is ON but ntccn < 1. Must set ccn_nc in field_table' + ENDIF + stop + ELSE + if (Model%me == Model%master) then + write(*,*) 'NSSL micro: CCN is ON' + ENDIF + IF ( Model%ntccna > 1 .and. Model%me == Model%master ) THEN + write(*,*) 'NSSL micro: CCNA is ON' + ENDIF + ENDIF + +! IF ( nssl_hail_on .eqv. .false. ) THEN + if (Model%me == Model%master) then + write(*,*) 'Model%nthl = ',Model%nthl + ENDIF + IF ( ( Model%nthl < 1 ) ) THEN ! check if hail is in the field_table. If not, set flag so the microphysics knows. + if (Model%me == Model%master) then + write(*,*) 'NSSL micro: hail is OFF' + ENDIF + nssl_hail_on = .false. + Model%nssl_hail_on = .false. + ! pretend that hail exists so that bad arrays are not passed to microphysics + Model%nthl = Max( 1, Model%ntgl ) + Model%nthv = Max( 1, Model%ntgv ) + Model%nthnc = Max( 1, Model%ntgnc ) + ELSE + nssl_hail_on = .true. + Model%nssl_hail_on = .true. + if (Model%me == Model%master) then + write(*,*) 'NSSL micro: hail is ON' + ENDIF + IF ( Model%nthv < 1 .or. Model%nthnc < 1 ) THEN + if (Model%me == Model%master) write(0,*) 'missing needed tracers for NSSL hail! nthl > 1 but either volume or number is not in field_table' + stop + ENDIF + ENDIF + + Model%nssl_hail_on = nssl_hail_on + + IF ( Model%ntgl < 1 .or. Model%ntgv < 1 .or. Model%ntgnc < 1 .or. & + Model%ntsw < 1 .or. Model%ntsnc < 1 .or. & + Model%ntrw < 1 .or. Model%ntrnc < 1 .or. & + Model%ntiw < 1 .or. Model%ntinc < 1 .or. & + Model%ntcw < 1 .or. Model%ntlnc < 1 & + ) THEN + if (Model%me == Model%master) write(0,*) 'missing needed tracers for NSSL!' + stop + ENDIF + + ENDIF !} + ! -- setup aerosol scavenging factors - allocate(Model%fscav(Model%ntchm)) + n = max(Model%ntrac, Model%ntchm) + allocate(Model%fscav(n)) + Model%fscav = -9999.0 if (Model%ntchm > 0) then ! -- initialize to default Model%fscav = 0.6_kind_phys @@ -4647,6 +4764,54 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & Model%nseffr = 3 if (Model%me == Model%master) print *,' Using wsm6 microphysics' + elseif (Model%imp_physics == Model%imp_physics_nssl2m) then !NSSL microphysics + Model%npdf3d = 0 + Model%num_p3d = 3 + Model%num_p2d = 1 + Model%pdfcld = .false. + Model%shcnvcw = .false. + IF ( Model%nssl_hail_on ) THEN + Model%ncnd = 6 + ELSE + Model%ncnd = 5 + ENDIF + Model%nleffr = 1 + Model%nieffr = 2 + Model%nseffr = 3 + Model%nreffr = 4 + Model%lradar = .true. + if (.not. Model%effr_in) then + Model%effr_in = .true. + ENDIF + if (Model%me == Model%master) print *,' Using NSSL double moment', & + ' microphysics', & + ' lradar =',Model%lradar,' ttendlim =',Model%ttendlim, & + Model%num_p3d,Model%num_p2d + + elseif (Model%imp_physics == Model%imp_physics_nssl2mccn) then !NSSL microphysics + Model%npdf3d = 0 + Model%num_p3d = 3 + Model%num_p2d = 1 + Model%pdfcld = .false. + Model%shcnvcw = .false. + IF ( Model%nssl_hail_on ) THEN + Model%ncnd = 6 + ELSE + Model%ncnd = 5 + ENDIF + Model%nleffr = 1 + Model%nieffr = 2 + Model%nseffr = 3 + Model%lradar = .true. + if (.not. Model%effr_in) then + Model%effr_in = .true. + ENDIF + + if (Model%me == Model%master) print *,' Using NSSL double moment', & + ' microphysics with CCN', & + ' lradar =',Model%lradar,' ttendlim =',Model%ttendlim, & + Model%num_p3d,Model%num_p2d + elseif (Model%imp_physics == Model%imp_physics_thompson) then !Thompson microphysics Model%npdf3d = 0 Model%num_p3d = 3 @@ -4967,6 +5132,14 @@ subroutine control_print(Model) print *, ' ttendlim : ', Model%ttendlim print *, ' ' endif + if (Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn) then + print *, ' NSSL microphysical parameters' + print *, ' CCCN : ', Model%nssl_cccn + print *, ' graupel shape : ', Model%nssl_alphah + print *, ' hail shape : ', Model%nssl_alphahl + print *, ' hail flag : ', Model%nssl_hail_on + print *, ' lradar : ', Model%lradar + endif if (Model%imp_physics == Model%imp_physics_mg) then print *, ' M-G microphysical parameters' print *, ' fprcp : ', Model%fprcp @@ -5166,6 +5339,9 @@ subroutine control_print(Model) print *, ' ' print *, 'stochastic physics' print *, ' do_sppt : ', Model%do_sppt + print *, ' pert_mp : ', Model%pert_mp + print *, ' pert_clds : ', Model%pert_clds + print *, ' pert_radtend : ', Model%pert_radtend print *, ' do_shum : ', Model%do_shum print *, ' do_skeb : ', Model%do_skeb print *, ' lndp_type : ', Model%lndp_type @@ -5206,12 +5382,18 @@ subroutine control_print(Model) print *, ' ntrw : ', Model%ntrw print *, ' ntsw : ', Model%ntsw print *, ' ntgl : ', Model%ntgl + print *, ' nthl : ', Model%nthl print *, ' ntclamt : ', Model%ntclamt print *, ' ntlnc : ', Model%ntlnc print *, ' ntinc : ', Model%ntinc print *, ' ntrnc : ', Model%ntrnc print *, ' ntsnc : ', Model%ntsnc print *, ' ntgnc : ', Model%ntgnc + print *, ' nthnc : ', Model%nthnc + print *, ' ntccn : ', Model%ntccn + print *, ' ntccna : ', Model%ntccna + print *, ' ntgv : ', Model%ntgv + print *, ' nthv : ', Model%nthv print *, ' ntke : ', Model%ntke print *, ' nto : ', Model%nto print *, ' nto2 : ', Model%nto2 @@ -5408,10 +5590,10 @@ subroutine tbd_create (Tbd, IM, Model) endif if (Model%do_sppt .or. Model%ca_global) then - allocate (Tbd%dtdtr (IM,Model%levs)) + allocate (Tbd%dtdtnp (IM,Model%levs)) allocate (Tbd%dtotprcp (IM)) allocate (Tbd%dcnvprcp (IM)) - Tbd%dtdtr = clear_val + Tbd%dtdtnp = clear_val Tbd%dtotprcp = clear_val Tbd%dcnvprcp = clear_val endif @@ -5683,7 +5865,8 @@ subroutine diag_create (Diag, IM, Model) allocate (Diag%skebv_wts(IM,Model%levs)) allocate (Diag%sppt_wts(IM,Model%levs)) allocate (Diag%shum_wts(IM,Model%levs)) - allocate (Diag%zmtnblck(IM)) + allocate (Diag%sfc_wts(IM,Model%n_var_lndp)) + allocate (Diag%zmtnblck(IM)) allocate (Diag%ca1 (IM)) allocate (Diag%ca2 (IM)) allocate (Diag%ca3 (IM)) @@ -5707,6 +5890,9 @@ subroutine diag_create (Diag, IM, Model) allocate (Diag%dt3dt (IM,Model%levs,11)) if (Model%qdiag3d) then allocate (Diag%dq3dt (IM,Model%levs,13)) + allocate (Diag%upd_mf (IM,Model%levs)) + allocate (Diag%dwn_mf (IM,Model%levs)) + allocate (Diag%det_mf (IM,Model%levs)) else allocate (Diag%dq3dt (1,1,13)) endif @@ -6004,9 +6190,9 @@ subroutine diag_phys_zero (Diag, Model, linit, iauwindow_center) Diag%dt3dt = zero if (Model%qdiag3d) then Diag%dq3dt = zero -! Diag%upd_mf = zero -! Diag%dwn_mf = zero -! Diag%det_mf = zero + Diag%upd_mf = zero + Diag%dwn_mf = zero + Diag%det_mf = zero endif endif @@ -6207,6 +6393,7 @@ subroutine interstitial_create (Interstitial, IM, Model) integer :: iGas ! allocate (Interstitial%otspt (Model%ntracp1,2)) + allocate (Interstitial%otsptflag (Model%ntracp1)) ! Set up numbers of tracers for PBL, convection, etc: sets ! Interstitial%{nncl,nvdiff,mg3_as_mg2,nn,tracers_total,ntiwx,ntk,ntkev,otspt,nsamftrac,ncstrac,nscav} call interstitial_setup_tracers(Interstitial, Model) @@ -6267,7 +6454,6 @@ subroutine interstitial_create (Interstitial, IM, Model) allocate (Interstitial%dqsfc1 (IM)) allocate (Interstitial%drain (IM)) allocate (Interstitial%dtdt (IM,Model%levs)) - allocate (Interstitial%dtdtc (IM,Model%levs)) allocate (Interstitial%dtsfc1 (IM)) allocate (Interstitial%dt_mf (IM,Model%levs)) allocate (Interstitial%dtzm (IM)) @@ -6529,7 +6715,8 @@ subroutine interstitial_create (Interstitial, IM, Model) end if ! ! Allocate arrays that are conditional on physics choices - if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson) then + if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson .or. & + Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn) then allocate (Interstitial%graupelmp (IM)) allocate (Interstitial%icemp (IM)) allocate (Interstitial%rainmp (IM)) @@ -6645,6 +6832,7 @@ subroutine interstitial_setup_tracers(Interstitial, Model) class(GFS_interstitial_type) :: Interstitial type(GFS_control_type), intent(in) :: Model integer :: n, tracers + integer :: ihail !first, initialize the values (in case the values don't get initialized within if statements below) Interstitial%nncl = Model%ncld @@ -6656,6 +6844,7 @@ subroutine interstitial_setup_tracers(Interstitial, Model) Interstitial%ntkev = 0 Interstitial%tracers_total = 0 Interstitial%otspt(:,:) = .true. + Interstitial%otsptflag(:) = .true. Interstitial%nsamftrac = 0 Interstitial%ncstrac = 0 Interstitial%nscav = Model%ntrac-Model%ncld+2 @@ -6671,6 +6860,21 @@ subroutine interstitial_setup_tracers(Interstitial, Model) endif if (Model%satmedmf) Interstitial%nvdiff = Interstitial%nvdiff + 1 Interstitial%nncl = 5 + elseif (Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn) then + IF ( Model%nssl_hail_on ) THEN + ihail = 1 + Interstitial%nvdiff = 17 ! Model%ntrac ! 17 + ELSE + ihail = 0 + Interstitial%nvdiff = 14 ! turn off hail q,N, and volume + ENDIF + ! write(*,*) 'NSSL: nvdiff, ntrac = ',Interstitial%nvdiff, Model%ntrac + if (Model%satmedmf) Interstitial%nvdiff = Interstitial%nvdiff + 1 + if (Model%imp_physics == Model%imp_physics_nssl2m) then + Interstitial%nncl = 4 + ihail + else + Interstitial%nncl = 5 + ihail + end if elseif (Model%imp_physics == Model%imp_physics_wsm6) then Interstitial%nvdiff = Model%ntrac -3 if (Model%satmedmf) Interstitial%nvdiff = Interstitial%nvdiff + 1 @@ -6717,6 +6921,16 @@ subroutine interstitial_setup_tracers(Interstitial, Model) Interstitial%ntiwx = 3 ! total ice or total condensate elseif (Model%imp_physics == Model%imp_physics_mg) then Interstitial%ntiwx = 3 + elseif (Model%imp_physics == Model%imp_physics_nssl2m .or. & + Model%imp_physics == Model%imp_physics_nssl2mccn) then +! Interstitial%ntqvx = 1 +! Interstitial%ntcwx = 2 + Interstitial%ntiwx = 3 + IF ( Model%imp_physics == Model%imp_physics_nssl2mccn ) THEN +! Interstitial%ntozx = 13 + 3*ihail + ELSE +! Interstitial%ntozx = 12 + 3*ihail + ENDIF else Interstitial%ntiwx = 0 endif @@ -6762,21 +6976,27 @@ subroutine interstitial_setup_tracers(Interstitial, Model) if (Model%cscnv .or. Model%satmedmf .or. Model%trans_trac ) then Interstitial%otspt(:,:) = .true. ! otspt is used only for cscnv Interstitial%otspt(1:3,:) = .false. ! this is for sp.hum, ice and liquid water + Interstitial%otsptflag(:) = .true. tracers = 2 do n=2,Model%ntrac if ( n /= Model%ntcw .and. n /= Model%ntiw .and. n /= Model%ntclamt .and. & n /= Model%ntrw .and. n /= Model%ntsw .and. n /= Model%ntrnc .and. & - n /= Model%ntsnc .and. n /= Model%ntgl .and. n /= Model%ntgnc) then + n /= Model%ntsnc .and. n /= Model%ntgl .and. n /= Model%ntgnc .and. & + n /= Model%nthl .and. n /= Model%nthnc .and. n /= Model%ntgv .and. & + n /= Model%nthv .and. n /= Model%ntccn .and. n /= Model%ntccna ) then tracers = tracers + 1 if (Model%ntke == n ) then Interstitial%otspt(tracers+1,1) = .false. Interstitial%ntk = tracers endif - if (Model%ntlnc == n .or. Model%ntinc == n .or. Model%ntrnc == n .or. Model%ntsnc == n .or. Model%ntgnc == n) & + if (Model%ntlnc == n .or. Model%ntinc == n .or. Model%ntrnc == n .or. Model%ntsnc == n & + .or. Model%ntgnc == n .or. Model%nthnc == n) & ! if (ntlnc == n .or. ntinc == n .or. ntrnc == n .or. ntsnc == n .or.& ! ntrw == n .or. ntsw == n .or. ntgl == n) & Interstitial%otspt(tracers+1,1) = .false. if (Interstitial%trans_aero .and. Model%ntchs == n) Interstitial%itc = tracers + else + Interstitial%otsptflag(n) = .false. endif enddo Interstitial%tracers_total = tracers - 2 @@ -6961,7 +7181,6 @@ subroutine interstitial_phys_reset (Interstitial, Model) Interstitial%drain = clear_val Interstitial%dt_mf = clear_val Interstitial%dtdt = clear_val - Interstitial%dtdtc = clear_val Interstitial%dtsfc1 = clear_val Interstitial%dtzm = clear_val Interstitial%dudt = clear_val @@ -7143,7 +7362,8 @@ subroutine interstitial_phys_reset (Interstitial, Model) end if ! ! Reset fields that are conditional on physics choices - if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson) then + if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson .or. & + Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn) then Interstitial%graupelmp = clear_val Interstitial%icemp = clear_val Interstitial%rainmp = clear_val @@ -7313,7 +7533,6 @@ subroutine interstitial_print(Interstitial, Model, mpirank, omprank, blkno) write (0,*) 'sum(Interstitial%dqsfc1 ) = ', sum(Interstitial%dqsfc1 ) write (0,*) 'sum(Interstitial%drain ) = ', sum(Interstitial%drain ) write (0,*) 'sum(Interstitial%dtdt ) = ', sum(Interstitial%dtdt ) - write (0,*) 'sum(Interstitial%dtdtc ) = ', sum(Interstitial%dtdtc ) write (0,*) 'sum(Interstitial%dtsfc1 ) = ', sum(Interstitial%dtsfc1 ) write (0,*) 'sum(Interstitial%dtzm ) = ', sum(Interstitial%dtzm ) write (0,*) 'sum(Interstitial%dt_mf ) = ', sum(Interstitial%dt_mf ) @@ -7530,8 +7749,9 @@ subroutine interstitial_print(Interstitial, Model, mpirank, omprank, blkno) end if ! ! Print arrays that are conditional on physics choices - if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson) then - write (0,*) 'Interstitial_print: values specific to GFDL/Thompson microphysics' + if (Model%imp_physics == Model%imp_physics_gfdl .or. Model%imp_physics == Model%imp_physics_thompson .or. & + Model%imp_physics == Model%imp_physics_nssl2m .or. Model%imp_physics == Model%imp_physics_nssl2mccn) then + write (0,*) 'Interstitial_print: values specific to GFDL/Thompson/NSSL microphysics' write (0,*) 'sum(Interstitial%graupelmp ) = ', sum(Interstitial%graupelmp ) write (0,*) 'sum(Interstitial%icemp ) = ', sum(Interstitial%icemp ) write (0,*) 'sum(Interstitial%rainmp ) = ', sum(Interstitial%rainmp ) diff --git a/ccpp/data/GFS_typedefs.meta b/ccpp/data/GFS_typedefs.meta index 3b1019119..f6eeac570 100644 --- a/ccpp/data/GFS_typedefs.meta +++ b/ccpp/data/GFS_typedefs.meta @@ -195,6 +195,13 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys +[qgrs(:,:,index_for_hail)] + standard_name = hail_mixing_ratio + long_name = ratio of mass of hail to mass of dry air plus vapor (without condensates) + units = kg kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys [qgrs(:,:,index_for_ozone)] standard_name = ozone_mixing_ratio long_name = ozone mixing ratio @@ -254,6 +261,41 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys +[qgrs(:,:,index_for_hail_number_concentration)] + standard_name = hail_number_concentration + long_name = number concentration of hail + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[qgrs(:,:,index_for_cloud_condensation_nuclei_number_concentration)] + standard_name = cloud_condensation_nuclei_number_concentration + long_name = number concentration of cloud condensation nuclei + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[qgrs(:,:,index_for_activated_cloud_condensation_nuclei_number_concentration)] + standard_name = activated_cloud_condensation_nuclei_number_concentration + long_name = number concentration of activated cloud condensation nuclei + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[qgrs(:,:,index_for_graupel_volume)] + standard_name = graupel_volume + long_name = graupel particle volume + units = m3 kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[qgrs(:,:,index_for_hail_volume)] + standard_name = hail_volume + long_name = hail particle volume + units = m3 kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys [qgrs(:,:,index_for_turbulent_kinetic_energy)] standard_name = turbulent_kinetic_energy long_name = turbulent kinetic energy @@ -383,6 +425,13 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys +[gq0(:,:,index_for_hail)] + standard_name = hail_mixing_ratio_updated_by_physics + long_name = ratio of mass of hail to mass of dry air plus vapor (without condensates) updated by physics + units = kg kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys [gq0(:,:,index_for_mass_weighted_rime_factor)] standard_name = mass_weighted_rime_factor_updated_by_physics long_name = mass weighted rime factor updated by physics @@ -442,6 +491,41 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys +[gq0(:,:,index_for_hail_number_concentration)] + standard_name = hail_number_concentration_updated_by_physics + long_name = number concentration of hail updated by physics + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[gq0(:,:,index_for_cloud_condensation_nuclei_number_concentration)] + standard_name = cloud_condensation_nuclei_number_concentration_updated_by_physics + long_name = number concentration of cloud condensation nuclei updated by physics + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[gq0(:,:,index_for_activated_cloud_condensation_nuclei_number_concentration)] + standard_name = activated_cloud_condensation_nuclei_number_concentration_updated_by_physics + long_name = number concentration of cloud condensation nuclei updated by physics + units = kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[gq0(:,:,index_for_graupel_volume)] + standard_name = graupel_volume_updated_by_physics + long_name = graupel volume updated by physics + units = m3 kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys +[gq0(:,:,index_for_hail_volume)] + standard_name = hail_volume_updated_by_physics + long_name = hail volume updated by physics + units = m3 kg-1 + dimensions = (horizontal_dimension,vertical_dimension) + type = real + kind = kind_phys [gq0(:,:,index_for_cloud_amount)] standard_name = cloud_fraction_updated_by_physics long_name = cloud fraction updated by physics @@ -2960,6 +3044,18 @@ units = km dimensions = () type = real +[imp_physics_nssl2m] + standard_name = flag_for_nssl2m_microphysics_scheme + long_name = choice of NSSL 2-moment microphysics scheme + units = flag + dimensions = () + type = integer +[imp_physics_nssl2mccn] + standard_name = flag_for_nssl2mccn_microphysics_scheme + long_name = choice of NSSL 2-moment microphysics scheme with CCN + units = flag + dimensions = () + type = integer [psautco] standard_name = coefficient_from_cloud_ice_to_snow long_name = auto conversion coeff from ice to snow @@ -3132,6 +3228,39 @@ dimensions = () type = character kind = len=16 +[nssl_cccn] + standard_name = nssl_ccn_concentration + long_name = CCN concentration + units = m-3 + dimensions = () + type = real + kind = kind_phys +[nssl_alphah] + standard_name = nssl_alpha_graupel + long_name = graupel PSD shape parameter in NSSL micro + units = none + dimensions = () + type = real + kind = kind_phys +[nssl_alphahl] + standard_name = nssl_alpha_hail + long_name = hail PSD shape parameter in NSSL micro + units = none + dimensions = () + type = real + kind = kind_phys +[nssl_hail_on] + standard_name = nssl_hail_on + long_name = hail activation flag in NSSL micro + units = none + dimensions = () + type = logical +[nssl_invertccn] + standard_name = nssl_invertccn + long_name = flag to invert CCN in NSSL micro + units = none + dimensions = () + type = logical [tf] standard_name = frozen_cloud_threshold_temperature long_name = threshold temperature below which all cloud is ice @@ -4273,6 +4402,31 @@ units = flag dimensions = () type = logical +[pert_mp] + standard_name = flag_for_stochastic_microphysics_perturbations + long_name = flag for stochastic microphysics perturbations + units = flag + dimensions = () + type = logical +[pert_clds] + standard_name = flag_for_stochastic_cloud_fraction_perturbations + long_name = flag for stochastic cloud fraction perturbations + units = flag + dimensions = () + type = logical +[sppt_amp] + standard_name = total_ampltiude_of_sppt_perturbation + long_name = toal ampltidue of stochastic sppt perturbation + units = none + dimensions = () + type = real + kind = kind_phys +[pert_radtend] + standard_name = flag_for_stochastic_radiative_heating_perturbations + long_name = flag for stochastic radiative heating perturbations + units = flag + dimensions = () + type = logical [use_zmtnblck] standard_name = flag_for_mountain_blocking long_name = flag for mountain blocking @@ -4371,6 +4525,12 @@ units = index dimensions = () type = integer +[nthl] + standard_name = index_for_hail + long_name = tracer index for hail + units = index + dimensions = () + type = integer [ntclamt] standard_name = index_for_cloud_amount long_name = tracer index for cloud amount integer @@ -4407,6 +4567,36 @@ units = index dimensions = () type = integer +[nthnc] + standard_name = index_for_hail_number_concentration + long_name = tracer index for hail number concentration + units = index + dimensions = () + type = integer +[ntccn] + standard_name = index_for_cloud_condensation_nuclei_number_concentration + long_name = tracer index for cloud condensation nuclei number concentration + units = index + dimensions = () + type = integer +[ntccna] + standard_name = index_for_activated_cloud_condensation_nuclei_number_concentration + long_name = tracer index for activated cloud condensation nuclei number concentration + units = index + dimensions = () + type = integer +[ntgv] + standard_name = index_for_graupel_volume + long_name = tracer index for graupel particle volume + units = index + dimensions = () + type = integer +[nthv] + standard_name = index_for_hail_volume + long_name = tracer index for hail particle volume + units = index + dimensions = () + type = integer [ntke] standard_name = index_for_turbulent_kinetic_energy long_name = tracer index for turbulent kinetic energy @@ -5392,10 +5582,10 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys -[dtdtr] - standard_name = tendency_of_air_temperature_due_to_radiative_heating_on_physics_time_step - long_name = temp. change due to radiative heating per time step - units = K +[dtdtnp] + standard_name = tendency_of_air_temperature_to_withold_from_sppt + long_name = temp. change from physics that should not be perturbed by sppt + units = K s-1 dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys @@ -6510,6 +6700,13 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys +[sfc_wts] + standard_name = weights_for_stochastic_surface_physics_perturbation_flipped + long_name = weights for stochastic surface physics perturbation, flipped + units = none + dimensions = (horizontal_loop_extent,number_of_land_surface_variables_perturbed) + type = real + kind = kind_phys [zmtnblck] standard_name = level_of_dividing_streamline long_name = level of the dividing streamline @@ -8091,13 +8288,6 @@ dimensions = (horizontal_loop_extent,vertical_dimension) type = real kind = kind_phys -[dtdtc] - standard_name = tendency_of_air_temperature_due_to_radiative_heating_assuming_clear_sky - long_name = clear sky radiative (shortwave + longwave) heating rate at current time - units = K s-1 - dimensions = (horizontal_loop_extent,vertical_dimension) - type = real - kind = kind_phys [dtsfc1] standard_name = instantaneous_surface_upward_sensible_heat_flux long_name = surface upward sensible heat flux @@ -9065,6 +9255,12 @@ units = flag dimensions = (number_of_tracers_plus_one,2) type = logical +[otsptflag] + standard_name = flag_convective_tracer_transport_interstitial + long_name = flag for interstitial tracer transport + units = flag + dimensions = (number_of_tracers_plus_one) + type = logical [oz_coeffp5] standard_name = number_of_coefficients_in_ozone_forcing_data_plus_five long_name = number of coefficients in ozone forcing data plus five diff --git a/ccpp/driver/GFS_diagnostics.F90 b/ccpp/driver/GFS_diagnostics.F90 index 3f9bafde7..fea122446 100644 --- a/ccpp/driver/GFS_diagnostics.F90 +++ b/ccpp/driver/GFS_diagnostics.F90 @@ -2033,6 +2033,28 @@ subroutine GFS_externaldiag_populate (ExtDiag, Model, Statein, Stateout, Sfcprop ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%shum_wts(:,:) enddo + idx = idx + 1 + ExtDiag(idx)%axes = 2 + ExtDiag(idx)%name = 'sfc_wts1' + ExtDiag(idx)%desc = 'perturbation amplitude' + ExtDiag(idx)%unit = 'none' + ExtDiag(idx)%mod_name = 'gfs_phys' + allocate (ExtDiag(idx)%data(nblks)) + do nb = 1,nblks + ExtDiag(idx)%data(nb)%var2 => IntDiag(nb)%sfc_wts(:,1) + enddo + + idx = idx + 1 + ExtDiag(idx)%axes = 2 + ExtDiag(idx)%name = 'sfc_wts2' + ExtDiag(idx)%desc = 'perturbation amplitude' + ExtDiag(idx)%unit = 'none' + ExtDiag(idx)%mod_name = 'gfs_phys' + allocate (ExtDiag(idx)%data(nblks)) + do nb = 1,nblks + ExtDiag(idx)%data(nb)%var2 => IntDiag(nb)%sfc_wts(:,2) + enddo + idx = idx + 1 ExtDiag(idx)%axes = 2 ExtDiag(idx)%name = 'ca1' @@ -2749,7 +2771,6 @@ subroutine GFS_externaldiag_populate (ExtDiag, Model, Statein, Stateout, Sfcprop do nb = 1,nblks ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%dq3dt(:,:,9) enddo - end if if_qdiag3d idx = idx + 1 ExtDiag(idx)%axes = 3 @@ -2799,6 +2820,44 @@ subroutine GFS_externaldiag_populate (ExtDiag, Model, Statein, Stateout, Sfcprop ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%dq3dt(:,:,13) enddo + idx = idx + 1 + ExtDiag(idx)%axes = 3 + ExtDiag(idx)%name = 'upd_mf' + ExtDiag(idx)%desc = 'updraft convective mass flux' + ExtDiag(idx)%unit = 'kg m-1 s-3' + ExtDiag(idx)%mod_name = 'gfs_phys' + ExtDiag(idx)%time_avg = .TRUE. + allocate (ExtDiag(idx)%data(nblks)) + do nb = 1,nblks + ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%upd_mf(:,:) + enddo + + idx = idx + 1 + ExtDiag(idx)%axes = 3 + ExtDiag(idx)%name = 'dwn_mf' + ExtDiag(idx)%desc = 'downdraft convective mass flux' + ExtDiag(idx)%unit = 'kg m-1 s-3' + ExtDiag(idx)%mod_name = 'gfs_phys' + ExtDiag(idx)%time_avg = .TRUE. + allocate (ExtDiag(idx)%data(nblks)) + do nb = 1,nblks + ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%dwn_mf(:,:) + enddo + + idx = idx + 1 + ExtDiag(idx)%axes = 3 + ExtDiag(idx)%name = 'det_mf' + ExtDiag(idx)%desc = 'detrainment convective mass flux' + ExtDiag(idx)%unit = 'kg m-1 s-3' + ExtDiag(idx)%mod_name = 'gfs_phys' + ExtDiag(idx)%time_avg = .TRUE. + allocate (ExtDiag(idx)%data(nblks)) + do nb = 1,nblks + ExtDiag(idx)%data(nb)%var3 => IntDiag(nb)%det_mf(:,:) + enddo + + end if if_qdiag3d + end if if_ldiag3d !rab diff --git a/ccpp/physics b/ccpp/physics index 1488fa603..47df2d4b8 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit 1488fa603a48706926d4b8a33608146b0e4ad458 +Subproject commit 47df2d4b83cadd635b8b0cf13a6a76f20250e0b7 diff --git a/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn.xml b/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn.xml index d15325ed6..155809729 100644 --- a/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn.xml +++ b/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn.xml @@ -40,7 +40,7 @@ - sfc_diff + mynnsfc_wrapper GFS_surface_loop_control_part1 sfc_ocean lsm_noah diff --git a/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn_lam3km.xml b/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn_lam3km.xml index 5b943bc00..351a80017 100644 --- a/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn_lam3km.xml +++ b/ccpp/suites/suite_FV3_GFS_v15_thompson_mynn_lam3km.xml @@ -39,7 +39,7 @@ - sfc_diff + mynnsfc_wrapper GFS_surface_loop_control_part1 sfc_nst_pre sfc_nst @@ -65,7 +65,6 @@ get_phi_fv3 GFS_suite_interstitial_3 GFS_suite_interstitial_4 - cnvc90 GFS_MP_generic_pre mp_thompson_pre mp_thompson diff --git a/ccpp/suites/suite_FV3_GFS_v16_ras.xml b/ccpp/suites/suite_FV3_GFS_v16_ras.xml new file mode 100644 index 000000000..4a7fc2b27 --- /dev/null +++ b/ccpp/suites/suite_FV3_GFS_v16_ras.xml @@ -0,0 +1,94 @@ + + + + + + + fv_sat_adj + + + + + GFS_time_vary_pre + GFS_rrtmg_setup + GFS_rad_time_vary + GFS_phys_time_vary + + + + + GFS_suite_interstitial_rad_reset + GFS_rrtmg_pre + rrtmg_sw_pre + rrtmg_sw + rrtmg_sw_post + rrtmg_lw_pre + rrtmg_lw + rrtmg_lw_post + GFS_rrtmg_post + + + + + GFS_suite_interstitial_phys_reset + GFS_suite_stateout_reset + get_prs_fv3 + GFS_suite_interstitial_1 + GFS_surface_generic_pre + GFS_surface_composites_pre + dcyc2t3 + GFS_surface_composites_inter + GFS_suite_interstitial_2 + + + + sfc_diff + GFS_surface_loop_control_part1 + sfc_nst_pre + sfc_nst + sfc_nst_post + lsm_noah + sfc_sice + GFS_surface_loop_control_part2 + + + + GFS_surface_composites_post + sfc_diag + sfc_diag_post + GFS_surface_generic_post + GFS_PBL_generic_pre + satmedmfvdifq + GFS_PBL_generic_post + GFS_GWD_generic_pre + cires_ugwp + cires_ugwp_post + GFS_GWD_generic_post + rayleigh_damp + GFS_suite_stateout_update + ozphys_2015 + h2ophys + get_phi_fv3 + GFS_suite_interstitial_3 + GFS_DCNV_generic_pre + rascnv + GFS_DCNV_generic_post + GFS_SCNV_generic_pre + samfshalcnv + GFS_SCNV_generic_post + GFS_suite_interstitial_4 + cnvc90 + GFS_MP_generic_pre + gfdl_cloud_microphys + GFS_MP_generic_post + maximum_hourly_diagnostics + + + + + GFS_stochastics + phys_tend + + + + diff --git a/ccpp/suites/suite_HMTW_2020_mp2.xml b/ccpp/suites/suite_HMTW_2020_mp2.xml new file mode 100644 index 000000000..aeba3d232 --- /dev/null +++ b/ccpp/suites/suite_HMTW_2020_mp2.xml @@ -0,0 +1,87 @@ + + + + + + + GFS_time_vary_pre + GFS_rrtmg_setup + GFS_rad_time_vary + GFS_phys_time_vary + + + + + + GFS_suite_interstitial_rad_reset + GFS_rrtmg_pre + rrtmg_sw_pre + rrtmg_sw + rrtmg_sw_post + rrtmg_lw_pre + rrtmg_lw + rrtmg_lw_post + GFS_rrtmg_post + + + + + GFS_suite_interstitial_phys_reset + GFS_suite_stateout_reset + get_prs_fv3 + GFS_suite_interstitial_1 + GFS_surface_generic_pre + GFS_surface_composites_pre + dcyc2t3 + GFS_surface_composites_inter + GFS_suite_interstitial_2 + + + + sfc_diff + GFS_surface_loop_control_part1 + sfc_nst_pre + sfc_nst + sfc_nst_post + lsm_noah + sfc_sice + GFS_surface_loop_control_part2 + + + + GFS_surface_composites_post + sfc_diag + sfc_diag_post + GFS_surface_generic_post + GFS_PBL_generic_pre + hedmf + GFS_PBL_generic_post + GFS_GWD_generic_pre + cires_ugwp + cires_ugwp_post + GFS_GWD_generic_post + rayleigh_damp + GFS_suite_stateout_update + ozphys + GFS_DCNV_generic_pre + get_phi_fv3 + GFS_suite_interstitial_3 + GFS_DCNV_generic_post + GFS_SCNV_generic_pre + GFS_SCNV_generic_post + GFS_suite_interstitial_4 + cnvc90 + GFS_MP_generic_pre + mp_nssl_2mom + GFS_MP_generic_post + maximum_hourly_diagnostics + + + + + GFS_stochastics + + + + diff --git a/ccpp/suites/suite_nssl_mp_no_nsst.xml b/ccpp/suites/suite_nssl_mp_no_nsst.xml new file mode 100644 index 000000000..e702ea0d2 --- /dev/null +++ b/ccpp/suites/suite_nssl_mp_no_nsst.xml @@ -0,0 +1,78 @@ + + + + + + + GFS_time_vary_pre + GFS_rrtmg_setup + GFS_rad_time_vary + GFS_phys_time_vary + + + + + GFS_suite_interstitial_rad_reset + GFS_rrtmg_pre + rrtmg_sw_pre + rrtmg_sw + rrtmg_sw_post + rrtmg_lw_pre + rrtmg_lw + rrtmg_lw_post + GFS_rrtmg_post + + + + + GFS_suite_interstitial_phys_reset + GFS_suite_stateout_reset + get_prs_fv3 + GFS_suite_interstitial_1 + GFS_surface_generic_pre + GFS_surface_composites_pre + dcyc2t3 + GFS_surface_composites_inter + GFS_suite_interstitial_2 + + + + sfc_diff + GFS_surface_loop_control_part1 + sfc_ocean + lsm_noah + sfc_sice + GFS_surface_loop_control_part2 + + + + GFS_surface_composites_post + sfc_diag + sfc_diag_post + GFS_surface_generic_post + GFS_PBL_generic_pre + hedmf + GFS_PBL_generic_post + GFS_GWD_generic_pre + cires_ugwp_post + GFS_GWD_generic_post + rayleigh_damp + GFS_suite_stateout_update + ozphys_2015 + h2ophys + get_phi_fv3 + GFS_suite_interstitial_3 + GFS_suite_interstitial_4 + GFS_MP_generic_pre + mp_nsslg + GFS_MP_generic_post + maximum_hourly_diagnostics + + + + + GFS_stochastics + + + + diff --git a/coupler_main.F90 b/coupler_main.F90 deleted file mode 100644 index 3f36fb550..000000000 --- a/coupler_main.F90 +++ /dev/null @@ -1,506 +0,0 @@ -!*********************************************************************** -!* GNU General Public License * -!* This file is a part of fvGFS. * -!* * -!* fvGFS is free software; you can redistribute it and/or modify it * -!* and are expected to follow the terms of the GNU General Public * -!* License as published by the Free Software Foundation; either * -!* version 2 of the License, or (at your option) any later version. * -!* * -!* fvGFS is distributed in the hope that it will be useful, but * -!* WITHOUT ANY WARRANTY; without even the implied warranty of * -!* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * -!* General Public License for more details. * -!* * -!* For the full text of the GNU General Public License, * -!* write to: Free Software Foundation, Inc., * -!* 675 Mass Ave, Cambridge, MA 02139, USA. * -!* or see: http://www.gnu.org/licenses/gpl.html * -!*********************************************************************** - -program coupler_main - -!----------------------------------------------------------------------- -! -! program that couples component models for the atmosphere, -! ocean (amip), land, and sea-ice using the exchange module -! -!----------------------------------------------------------------------- - -use time_manager_mod, only: time_type, set_calendar_type, set_time, & - set_date, days_in_month, month_name, & - operator(+), operator (<), operator (>), & - operator (/=), operator (/), operator (==),& - operator (*), THIRTY_DAY_MONTHS, JULIAN, & - NOLEAP, NO_CALENDAR, date_to_string, & - get_date - -use atmos_model_mod, only: atmos_model_init, atmos_model_end, & - update_atmos_model_dynamics, & - update_atmos_radiation_physics, & - update_atmos_model_state, & - atmos_data_type, atmos_model_restart - -use constants_mod, only: constants_init -#ifdef INTERNAL_FILE_NML -use mpp_mod, only: input_nml_file -#else -use fms_mod, only: open_namelist_file -#endif -use fms_mod, only: file_exist, check_nml_error, & - error_mesg, fms_init, fms_end, close_file, & - write_version_number, uppercase - -use mpp_mod, only: mpp_init, mpp_pe, mpp_root_pe, mpp_npes, mpp_get_current_pelist, & - mpp_set_current_pelist, stdlog, mpp_error, NOTE, FATAL, WARNING -use mpp_mod, only: mpp_clock_id, mpp_clock_begin, mpp_clock_end, mpp_sync - -use mpp_io_mod, only: mpp_open, mpp_close, & - MPP_NATIVE, MPP_RDONLY, MPP_DELETE - -use mpp_domains_mod, only: mpp_get_global_domain, mpp_global_field, CORNER -use memutils_mod, only: print_memuse_stats -use sat_vapor_pres_mod,only: sat_vapor_pres_init - -use diag_manager_mod, only: diag_manager_init, diag_manager_end, & - get_base_date, diag_manager_set_time_end - -use data_override_mod, only: data_override_init - - -implicit none - -!----------------------------------------------------------------------- - -character(len=128) :: version = '$Id: coupler_main.F90,v 19.0.4.1.2.3 2014/09/09 23:51:59 Rusty.Benson Exp $' -character(len=128) :: tag = '$Name: ulm_201505 $' - -!----------------------------------------------------------------------- -!---- model defined-types ---- - - type (atmos_data_type) :: Atm - -!----------------------------------------------------------------------- -! ----- coupled model time ----- - - type (time_type) :: Time_atmos, Time_init, Time_end, & - Time_step_atmos, Time_step_ocean, & - Time_restart, Time_step_restart - integer :: num_cpld_calls, num_atmos_calls, nc, na, ret - -! ----- coupled model initial date ----- - - integer :: date_init(6) - integer :: calendar_type = -99 - -! ----- timing flags ----- - - integer :: initClock, mainClock, termClock - integer, parameter :: timing_level = 1 - -! ----- namelist ----- - integer, dimension(6) :: current_date = (/ 0, 0, 0, 0, 0, 0 /) - character(len=17) :: calendar = ' ' - logical :: force_date_from_namelist = .false. ! override restart values for date - integer :: months=0, days=0, hours=0, minutes=0, seconds=0 - integer :: dt_atmos = 0 - integer :: dt_ocean = 0 - integer :: restart_days = 0 - integer :: restart_secs = 0 - integer :: atmos_nthreads = 1 - logical :: memuse_verbose = .false. - logical :: use_hyper_thread = .false. - logical :: debug_affinity = .false. - integer :: ncores_per_node = 0 - - namelist /coupler_nml/ current_date, calendar, force_date_from_namelist, & - months, days, hours, minutes, seconds, & - dt_atmos, dt_ocean, atmos_nthreads, memuse_verbose, & - use_hyper_thread, ncores_per_node, debug_affinity, & - restart_secs, restart_days - -! ----- local variables ----- - character(len=32) :: timestamp - logical :: intrm_rst - -!####################################################################### - - call fms_init() - call mpp_init() - initClock = mpp_clock_id( 'Initialization' ) - call mpp_clock_begin (initClock) !nesting problem - - call fms_init - call constants_init - call sat_vapor_pres_init - - call coupler_init - call print_memuse_stats('after coupler init') - - call mpp_set_current_pelist() - call mpp_clock_end (initClock) !end initialization - mainClock = mpp_clock_id( 'Main loop' ) - termClock = mpp_clock_id( 'Termination' ) - call mpp_clock_begin(mainClock) !begin main loop - - do nc = 1, num_cpld_calls - - Time_atmos = Time_atmos + Time_step_atmos - - call update_atmos_model_dynamics (Atm) - - call update_atmos_radiation_physics (Atm) - - call update_atmos_model_state (Atm) - -!--- intermediate restart - if (intrm_rst) then - if ((nc /= num_cpld_calls) .and. (Time_atmos == Time_restart)) then - timestamp = date_to_string (Time_restart) - call atmos_model_restart(Atm, timestamp) - call coupler_res(timestamp) - Time_restart = Time_restart + Time_step_restart - endif - endif - - call print_memuse_stats('after full step') - - enddo - -!----------------------------------------------------------------------- - -#ifdef AVEC_TIMERS - call avec_timers_output -#endif - call mpp_set_current_pelist() - call mpp_clock_end(mainClock) - call mpp_clock_begin(termClock) - - call coupler_end - call mpp_set_current_pelist() - call mpp_clock_end(termClock) - - call fms_end - -!----------------------------------------------------------------------- - - stop - -contains - -!####################################################################### - - subroutine coupler_init - -!----------------------------------------------------------------------- -! initialize all defined exchange grids and all boundary maps -!----------------------------------------------------------------------- - integer :: total_days, total_seconds, unit, ierr, io - integer :: n, gnlon, gnlat - integer :: date(6), flags - type (time_type) :: Run_length - character(len=9) :: month - logical :: use_namelist - - logical, allocatable, dimension(:,:) :: mask - real, allocatable, dimension(:,:) :: glon_bnd, glat_bnd - integer :: omp_get_thread_num, get_cpu_affinity, base_cpu -!----------------------------------------------------------------------- -!----- initialization timing identifiers ---- - -!----- read namelist ------- -!----- for backwards compatibilty read from file coupler.nml ----- - -#ifdef INTERNAL_FILE_NML - read(input_nml_file, nml=coupler_nml, iostat=io) - ierr = check_nml_error(io, 'coupler_nml') -#else - if (file_exist('input.nml')) then - unit = open_namelist_file () - else - call error_mesg ('program coupler', & - 'namelist file input.nml does not exist', FATAL) - endif - - ierr=1 - do while (ierr /= 0) - read (unit, nml=coupler_nml, iostat=io, end=10) - ierr = check_nml_error (io, 'coupler_nml') - enddo -10 call close_file (unit) -#endif - -!----- write namelist to logfile ----- - call write_version_number (version, tag) - if (mpp_pe() == mpp_root_pe()) write(stdlog(),nml=coupler_nml) - -!----- allocate and set the pelist (to the global pelist) ----- - allocate( Atm%pelist (mpp_npes()) ) - call mpp_get_current_pelist(Atm%pelist) - -!----- read restart file ----- - - if (file_exist('INPUT/coupler.res')) then - call mpp_open( unit, 'INPUT/coupler.res', action=MPP_RDONLY ) - read (unit,*,err=999) calendar_type - read (unit,*) date_init - read (unit,*) date - goto 998 !back to fortran-4 - ! read old-style coupler.res - 999 call mpp_close (unit) - call mpp_open (unit, 'INPUT/coupler.res', action=MPP_RDONLY, form=MPP_NATIVE) - read (unit) calendar_type - read (unit) date - 998 call mpp_close(unit) - else - force_date_from_namelist = .true. - endif - -!----- use namelist value (either no restart or override flag on) --- - - if ( force_date_from_namelist ) then - - if ( sum(current_date) <= 0 ) then - call error_mesg ('program coupler', & - 'no namelist value for current_date', FATAL) - else - date = current_date - endif - -!----- override calendar type with namelist value ----- - - select case( uppercase(trim(calendar)) ) - case( 'JULIAN' ) - calendar_type = JULIAN - case( 'NOLEAP' ) - calendar_type = NOLEAP - case( 'THIRTY_DAY' ) - calendar_type = THIRTY_DAY_MONTHS - case( 'NO_CALENDAR' ) - calendar_type = NO_CALENDAR - case default - call mpp_error ( FATAL, 'COUPLER_MAIN: coupler_nml entry calendar must '// & - 'be one of JULIAN|NOLEAP|THIRTY_DAY|NO_CALENDAR.' ) - end select - - endif - -!$ base_cpu = get_cpu_affinity() -!$ call omp_set_num_threads(atmos_nthreads) -!$OMP PARALLEL NUM_THREADS(atmos_nthreads) -!$ if(omp_get_thread_num() < atmos_nthreads/2 .OR. (.not. use_hyper_thread)) then -!$ call set_cpu_affinity(base_cpu + omp_get_thread_num()) -!$ else -!$ call set_cpu_affinity(base_cpu + omp_get_thread_num() + & -!$ ncores_per_node - atmos_nthreads/2) -!$ endif -!$ if (debug_affinity) then -!$ write(6,*) mpp_pe()," atmos ",get_cpu_affinity(), base_cpu, omp_get_thread_num() -!$ call flush(6) -!$ endif -!$OMP END PARALLEL - - call set_calendar_type (calendar_type) - -!----- write current/initial date actually used to logfile file ----- - - if ( mpp_pe() == mpp_root_pe() ) then - write (stdlog(),16) date(1),trim(month_name(date(2))),date(3:6) - endif - - 16 format (' current date used = ',i4,1x,a,2i3,2(':',i2.2),' gmt') - -!----------------------------------------------------------------------- -!------ initialize diagnostics manager ------ - - call diag_manager_init (TIME_INIT=date) - -!----- always override initial/base date with diag_manager value ----- - - call get_base_date ( date_init(1), date_init(2), date_init(3), & - date_init(4), date_init(5), date_init(6) ) - -!----- use current date if no base date ------ - - if ( date_init(1) == 0 ) date_init = date - -!----- set initial and current time types ------ - - Time_init = set_date (date_init(1), date_init(2), date_init(3), & - date_init(4), date_init(5), date_init(6)) - - Time_atmos = set_date (date(1), date(2), date(3), & - date(4), date(5), date(6)) - -!----------------------------------------------------------------------- -!----- compute the ending time (compute days in each month first) ----- -! -! (NOTE: if run length in months then starting day must be <= 28) - - if ( months > 0 .and. date(3) > 28 ) & - call error_mesg ('program coupler', & - 'if run length in months then starting day must be <= 28', FATAL) - - Time_end = Time_atmos - total_days = 0 - do n = 1, months - total_days = total_days + days_in_month(Time_end) - Time_end = Time_atmos + set_time (0,total_days) - enddo - - total_days = total_days + days - total_seconds = hours*3600 + minutes*60 + seconds - Run_length = set_time (total_seconds,total_days) - Time_end = Time_atmos + Run_length - - !Need to pass Time_end into diag_manager for multiple thread case. - call diag_manager_set_time_end(Time_end) - - -!----------------------------------------------------------------------- -!----- write time stamps (for start time and end time) ------ - - call mpp_open( unit, 'time_stamp.out', nohdrs=.TRUE. ) - - month = month_name(date(2)) - if ( mpp_pe() == mpp_root_pe() ) write (unit,20) date, month(1:3) - - call get_date (Time_end, date(1), date(2), date(3), & - date(4), date(5), date(6)) - month = month_name(date(2)) - if ( mpp_pe() == mpp_root_pe() ) write (unit,20) date, month(1:3) - - call mpp_close (unit) - - 20 format (6i4,2x,a3) - -!----------------------------------------------------------------------- -!----- compute the time steps ------ - -Time_step_atmos = set_time (dt_atmos,0) -Time_step_ocean = set_time (dt_ocean,0) -num_cpld_calls = Run_length / Time_step_ocean -num_atmos_calls = Time_step_ocean / Time_step_atmos -Time_step_restart = set_time (restart_secs, restart_days) -Time_restart = Time_atmos + Time_step_restart -intrm_rst = .false. -if (restart_days > 0 .or. restart_secs > 0) intrm_rst = .true. - -!----------------------------------------------------------------------- -!------------------- some error checks --------------------------------- - -!----- initial time cannot be greater than current time ------- - - if ( Time_init > Time_atmos ) call error_mesg ('program coupler', & - 'initial time is greater than current time', FATAL) - -!----- make sure run length is a multiple of ocean time step ------ - - if ( num_cpld_calls * Time_step_ocean /= Run_length ) & - call error_mesg ('program coupler', & - 'run length must be multiple of ocean time step', FATAL) - -! ---- make sure cpld time step is a multiple of atmos time step ---- - - if ( num_atmos_calls * Time_step_atmos /= Time_step_ocean ) & - call error_mesg ('program coupler', & - 'atmos time step is not a multiple of the ocean time step', FATAL) - -!------ initialize component models ------ - - call atmos_model_init (Atm, Time_init, Time_atmos, Time_step_atmos) - - call print_memuse_stats('after atmos model init') - - call mpp_get_global_domain(Atm%Domain, xsize=gnlon, ysize=gnlat) - allocate ( glon_bnd(gnlon+1,gnlat+1), glat_bnd(gnlon+1,gnlat+1) ) - call mpp_global_field(Atm%Domain, Atm%lon_bnd, glon_bnd, position=CORNER) - call mpp_global_field(Atm%Domain, Atm%lat_bnd, glat_bnd, position=CORNER) - - call data_override_init ( ) ! Atm_domain_in = Atm%domain, & - ! Ice_domain_in = Ice%domain, & - ! Land_domain_in = Land%domain ) - -!----------------------------------------------------------------------- -!---- open and close dummy file in restart dir to check if dir exists -- - - if (mpp_pe() == 0 ) then - call mpp_open( unit, 'RESTART/file' ) - call mpp_close(unit, MPP_DELETE) - endif - -!----------------------------------------------------------------------- - - end subroutine coupler_init - -!####################################################################### - subroutine coupler_res(timestamp) - character(len=32), intent(in) :: timestamp - - integer :: unit, date(6) - -!----- compute current date ------ - - call get_date (Time_atmos, date(1), date(2), date(3), & - date(4), date(5), date(6)) - -!----- write restart file ------ - - if (mpp_pe() == mpp_root_pe())then - call mpp_open( unit, 'RESTART/'//trim(timestamp)//'.coupler.res', nohdrs=.TRUE. ) - write( unit, '(i6,8x,a)' )calendar_type, & - '(Calendar: no_calendar=0, thirty_day_months=1, julian=2, gregorian=3, noleap=4)' - - write( unit, '(6i6,8x,a)' )date_init, & - 'Model start time: year, month, day, hour, minute, second' - write( unit, '(6i6,8x,a)' )date, & - 'Current model time: year, month, day, hour, minute, second' - call mpp_close(unit) - endif - end subroutine coupler_res - -!####################################################################### - - subroutine coupler_end - - integer :: unit, date(6) -!----------------------------------------------------------------------- - - call atmos_model_end (Atm) - -!----- compute current date ------ - - call get_date (Time_atmos, date(1), date(2), date(3), & - date(4), date(5), date(6)) - -!----- check time versus expected ending time ---- - - if (Time_atmos /= Time_end) call error_mesg ('program coupler', & - 'final time does not match expected ending time', WARNING) - -!----- write restart file ------ - - call mpp_open( unit, 'RESTART/coupler.res', nohdrs=.TRUE. ) - if (mpp_pe() == mpp_root_pe())then - write( unit, '(i6,8x,a)' )calendar_type, & - '(Calendar: no_calendar=0, thirty_day_months=1, julian=2, gregorian=3, noleap=4)' - - write( unit, '(6i6,8x,a)' )date_init, & - 'Model start time: year, month, day, hour, minute, second' - write( unit, '(6i6,8x,a)' )date, & - 'Current model time: year, month, day, hour, minute, second' - endif - call mpp_close(unit) - -!----- final output of diagnostic fields ---- - - call diag_manager_end (Time_atmos) - -!----------------------------------------------------------------------- - - end subroutine coupler_end - -!####################################################################### - -end program coupler_main - diff --git a/io/inline_post_stub.F90 b/io/inline_post_stub.F90 index 35c400ba0..f33c78d6e 100644 --- a/io/inline_post_stub.F90 +++ b/io/inline_post_stub.F90 @@ -34,23 +34,20 @@ subroutine inline_post_run(wrt_int_state,mypei,mpicomp,lead_write, & integer,intent(in) :: mynfmin integer,intent(in) :: mynfsec ! - print *,'in stub post_run_gfs - not supported on this machine, return' + print *,'in stub inline_post_run - not supported on this machine, return' ! end subroutine inline_post_run ! !----------------------------------------------------------------------- ! - subroutine inline_post_getattr(wrt_int_state, fldbundle) -! - use esmf + subroutine inline_post_getattr(wrt_int_state) ! implicit none ! type(wrt_internal_state),intent(inout) :: wrt_int_state - type(ESMF_FieldBundle), intent(in) :: fldbundle ! ! - print *,'in stub post_getattr_gfs - not supported on this machine, return' + print *,'in stub inline_post_getattr - not supported on this machine, return' ! end subroutine inline_post_getattr diff --git a/io/post_gfs_stub.F90 b/io/post_gfs_stub.F90 deleted file mode 100644 index 3e61e31a1..000000000 --- a/io/post_gfs_stub.F90 +++ /dev/null @@ -1,58 +0,0 @@ -!----------------------------------------------------------------------- -!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -!----------------------------------------------------------------------- -! -module post_gfs - - use module_fv3_io_def, only : wrttasks_per_group,filename_base - use write_internal_state, only : wrt_internal_state - - implicit none - - public post_run_gfs, post_getattr_gfs - - contains - - subroutine post_run_gfs(wrt_int_state,mypei,mpicomp,lead_write, & - mynfhr,mynfmin,mynfsec) -! -! revision history: -! Jul 2019 J. Wang create interface to run inline post for FV3 -! -! -!----------------------------------------------------------------------- -! - implicit none -! -!----------------------------------------------------------------------- -! - type(wrt_internal_state),intent(in) :: wrt_int_state - integer,intent(in) :: mypei - integer,intent(in) :: mpicomp - integer,intent(in) :: lead_write - integer,intent(in) :: mynfhr - integer,intent(in) :: mynfmin - integer,intent(in) :: mynfsec -! - print *,'in stub post_run_gfs - not supported on this machine, return' -! - end subroutine post_run_gfs -! -!----------------------------------------------------------------------- -! - subroutine post_getattr_gfs(wrt_int_state, fldbundle) -! - use esmf -! - implicit none -! - type(wrt_internal_state),intent(inout) :: wrt_int_state - type(ESMF_FieldBundle), intent(in) :: fldbundle -! -! - print *,'in stub post_getattr_gfs - not supported on this machine, return' -! - end subroutine post_getattr_gfs - - - end module post_gfs diff --git a/io/post_regional.F90 b/io/post_regional.F90 index 953461a1f..eecd5c10f 100644 --- a/io/post_regional.F90 +++ b/io/post_regional.F90 @@ -250,8 +250,8 @@ subroutine post_getattr_regional(wrt_int_state) fldbundle = wrt_int_state%wrtFB(nfb) ! set grid spec: - if(mype==0) print*,'in post_getattr_lam,output_grid=',trim(output_grid),'nfb=',nfb - if(mype==0) print*,'in post_getattr_lam, lon1=',lon1,lon2,lat1,lat2,dlon,dlat +! if(mype==0) print*,'in post_getattr_lam,output_grid=',trim(output_grid),'nfb=',nfb +! if(mype==0) print*,'in post_getattr_lam, lon1=',lon1,lon2,lat1,lat2,dlon,dlat gdsdegr = 1000000. if(trim(output_grid) == 'regional_latlon') then @@ -274,8 +274,8 @@ subroutine post_getattr_regional(wrt_int_state) dxval = dlon*gdsdegr dyval = dlat*gdsdegr - if(mype==0) print*,'lonstart,latstart,dyval,dxval', & - lonstart,lonlast,latstart,latlast,dyval,dxval +! if(mype==0) print*,'lonstart,latstart,dyval,dxval', & +! lonstart,lonlast,latstart,latlast,dyval,dxval else if(trim(output_grid) == 'lambert_conformal') then MAPTYPE=1 @@ -347,8 +347,8 @@ subroutine post_getattr_regional(wrt_int_state) dyval = spval endif - if(mype==0) print*,'rotated latlon,lonstart,latstart,cenlon,cenlat,dyval,dxval', & - lonstart_r,lonlast_r,latstart_r,latlast_r,cenlon,cenlat,dyval,dxval +! if(mype==0) print*,'rotated latlon,lonstart,latstart,cenlon,cenlat,dyval,dxval', & +! lonstart_r,lonlast_r,latstart_r,latlast_r,cenlon,cenlat,dyval,dxval endif ! look at the field bundle attributes @@ -557,8 +557,8 @@ subroutine set_postvars_regional(wrt_int_state,mpicomp,setvar_atmfile, & tsrfc = tprec tmaxmin = tprec td3d = tprec - if(mype==0)print*,'MP_PHYSICS= ',imp_physics,'nbdl=',nbdl, 'tprec=',tprec,'tclod=',tclod, & - 'dtp=',dtp,'tmaxmin=',tmaxmin,'jsta=',jsta,jend,im,jm +! if(mype==0)print*,'MP_PHYSICS= ',imp_physics,'nbdl=',nbdl, 'tprec=',tprec,'tclod=',tclod, & +! 'dtp=',dtp,'tmaxmin=',tmaxmin,'jsta=',jsta,jend,im,jm ! !$omp parallel do default(shared),private(i,j) @@ -803,7 +803,7 @@ subroutine set_postvars_regional(wrt_int_state,mpicomp,setvar_atmfile, & idat(4) = wrt_int_state%fdate(4) idat(5) = wrt_int_state%fdate(5) ! - if(mype==0) print *,'idat=',idat,'sdat=',sdat,'ihrst=',ihrst +! if(mype==0) print *,'idat=',idat,'sdat=',sdat,'ihrst=',ihrst ! CALL W3DIFDAT(JDATE,IDATE,0,RINC) ! ! if(mype==0)print *,' rinc=',rinc @@ -925,14 +925,14 @@ subroutine set_postvars_regional(wrt_int_state,mpicomp,setvar_atmfile, & if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) return ! bail out - if(mype==0) print *,'in setvar, allocate fcstField,ibdl=',ibdl,'count=',ncount_field,'wrtFBname=',trim(wrtFBName) +! if(mype==0) print *,'in setvar, allocate fcstField,ibdl=',ibdl,'count=',ncount_field,'wrtFBname=',trim(wrtFBName) allocate(fcstField(ncount_field)) call ESMF_FieldBundleGet(wrt_int_state%wrtFB(ibdl), & fieldList=fcstField, itemorderflag=ESMF_ITEMORDER_ADDORDER, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) return ! bail out - if(mype==0) print *,'in setvar, read field, ibdl=',ibdl, 'nfield=',ncount_field +! if(mype==0) print *,'in setvar, read field, ibdl=',ibdl, 'nfield=',ncount_field do n=1, ncount_field ! call ESMF_FieldGet(fcstField(n),typekind=typekind, name=fieldname, & diff --git a/module_fcst_grid_comp.F90 b/module_fcst_grid_comp.F90 index 9c7b64bf1..776a89bf4 100644 --- a/module_fcst_grid_comp.F90 +++ b/module_fcst_grid_comp.F90 @@ -73,6 +73,7 @@ module module_fcst_grid_comp quilting, calendar_type, cpl, & cplprint_flag, force_date_from_configure, & num_restart_interval, frestart, restart_endfcst + use get_stochy_pattern_mod, only: write_stoch_restart_atm ! !----------------------------------------------------------------------- ! @@ -827,6 +828,7 @@ subroutine fcst_run_phase_2(fcst_comp, importState, exportState,clock,rc) atm_int_state%Time_restart = atm_int_state%Time_atstart + restart_inctime timestamp = date_to_string (atm_int_state%Time_restart) call atmos_model_restart(atm_int_state%Atm, timestamp) + call write_stoch_restart_atm('RESTART/'//trim(timestamp)//'.atm_stoch.res.nc') call wrt_atmres_timestamp(atm_int_state,timestamp) endif diff --git a/stochastic_physics/stochastic_physics_wrapper.F90 b/stochastic_physics/stochastic_physics_wrapper.F90 index 479270a9f..0d0b29164 100644 --- a/stochastic_physics/stochastic_physics_wrapper.F90 +++ b/stochastic_physics/stochastic_physics_wrapper.F90 @@ -95,9 +95,15 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr) initalize_stochastic_physics: if (GFS_Control%kdt==0) then if (GFS_Control%do_sppt .OR. GFS_Control%do_shum .OR. GFS_Control%do_skeb .OR. (GFS_Control%lndp_type .GT. 0) ) then + allocate(xlat(1:Atm_block%nblks,maxval(GFS_Control%blksz))) + allocate(xlon(1:Atm_block%nblks,maxval(GFS_Control%blksz))) + do nb=1,Atm_block%nblks + xlat(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlat(:) + xlon(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlon(:) + end do ! Initialize stochastic physics - call init_stochastic_physics(GFS_Control%levs, GFS_Control%blksz, GFS_Control%dtp, & - GFS_Control%input_nml_file, GFS_Control%fn_nml, GFS_Control%nlunit, GFS_Control%do_sppt, GFS_Control%do_shum, & + call init_stochastic_physics(GFS_Control%levs, GFS_Control%blksz, GFS_Control%dtp, GFS_Control%sppt_amp, & + GFS_Control%input_nml_file, GFS_Control%fn_nml, GFS_Control%nlunit, xlon, xlat, GFS_Control%do_sppt, GFS_Control%do_shum, & GFS_Control%do_skeb, GFS_Control%lndp_type, GFS_Control%n_var_lndp, GFS_Control%use_zmtnblck, GFS_Control%skeb_npass, & GFS_Control%lndp_var_list, GFS_Control%lndp_prt_list, & GFS_Control%ak, GFS_Control%bk, nthreads, GFS_Control%master, GFS_Control%communicator, ierr) @@ -106,8 +112,6 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr) return endif end if - allocate(xlat(1:Atm_block%nblks,maxval(GFS_Control%blksz))) - allocate(xlon(1:Atm_block%nblks,maxval(GFS_Control%blksz))) if (GFS_Control%do_sppt) then allocate(sppt_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs)) end if @@ -143,17 +147,13 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr) allocate(zorll(1:Atm_block%nblks,maxval(GFS_Control%blksz))) endif - do nb=1,Atm_block%nblks - xlat(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlat(:) - xlon(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlon(:) - end do if ( GFS_Control%lndp_type .EQ. 1 ) then ! this scheme sets perts once allocate(sfc_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%n_var_lndp)) - call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%phour, GFS_Control%blksz, xlat=xlat, xlon=xlon, & + call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%fhour, GFS_Control%blksz, & sppt_wts=sppt_wts, shum_wts=shum_wts, skebu_wts=skebu_wts, skebv_wts=skebv_wts, sfc_wts=sfc_wts, & nthreads=nthreads) - ! Copy contiguous data back; no need to copy xlat/xlon, these are intent(in) - just deallocate + ! Copy contiguous data back do nb=1,Atm_block%nblks GFS_Data(nb)%Coupling%sfc_wts(:,:) = sfc_wts(nb,1:GFS_Control%blksz(nb),:) end do @@ -170,12 +170,11 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr) endif else initalize_stochastic_physics - if (GFS_Control%do_sppt .OR. GFS_Control%do_shum .OR. GFS_Control%do_skeb .OR. (GFS_Control%lndp_type .EQ. 2) ) then - call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%phour, GFS_Control%blksz, xlat=xlat, xlon=xlon, & + call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%fhour, GFS_Control%blksz, & sppt_wts=sppt_wts, shum_wts=shum_wts, skebu_wts=skebu_wts, skebv_wts=skebv_wts, sfc_wts=sfc_wts, & nthreads=nthreads) - ! Copy contiguous data back; no need to copy xlat/xlon, these are intent(in) - just deallocate + ! Copy contiguous data back if (GFS_Control%do_sppt) then do nb=1,Atm_block%nblks GFS_Data(nb)%Coupling%sppt_wts(:,:) = sppt_wts(nb,1:GFS_Control%blksz(nb),:) @@ -239,6 +238,7 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr) GFS_Control%n_var_lndp, GFS_Control%lndp_var_list, GFS_Control%lndp_prt_list, & sfc_wts, xlon, xlat, stype, GFS_Control%pores, GFS_Control%resid,param_update_flag, & smc, slc, stc, vfrac, alvsf, alnsf, alvwf, alnwf, facsf, facwf, snoalb, semis, zorll, ierr) + if (ierr/=0) then write(6,*) 'call to GFS_apply_lndp failed' return