diff --git a/mediator/esmFldsExchange_nems_mod.F90 b/mediator/esmFldsExchange_nems_mod.F90 index f93739618..54f3bad6c 100644 --- a/mediator/esmFldsExchange_nems_mod.F90 +++ b/mediator/esmFldsExchange_nems_mod.F90 @@ -38,6 +38,8 @@ subroutine esmFldsExchange_nems(gcomp, phase, rc) use esmFlds , only : addmap_from => med_fldList_addmap_from use esmFlds , only : addfld_aoflux => med_fldList_addfld_aoflux use esmFlds , only : addmap_aoflux => med_fldList_addmap_aoflux + use esmFlds , only : addfld_ocnalb => med_fldList_addfld_ocnalb + use esmFlds , only : addmap_ocnalb => med_fldList_addmap_ocnalb ! input/output parameters: type(ESMF_GridComp) :: gcomp @@ -172,6 +174,14 @@ subroutine esmFldsExchange_nems(gcomp, phase, rc) call addfld_from(compice, 'mean_sw_pen_to_ocn') end if + ! Advertise the ocean albedos. These are not sent to the ATM in UFS. + if (phase == 'advertise') then + call addfld_ocnalb('So_avsdr') + call addfld_ocnalb('So_avsdf') + call addfld_ocnalb('So_anidr') + call addfld_ocnalb('So_anidf') + end if + !===================================================================== ! FIELDS TO ATMOSPHERE (compatm) !===================================================================== @@ -306,7 +316,7 @@ subroutine esmFldsExchange_nems(gcomp, phase, rc) else if ( fldchk(is_local%wrap%FBexp(compatm) , 'Sw_z0', rc=rc) .and. & fldchk(is_local%wrap%FBImp(compwav,compwav), 'Sw_z0', rc=rc)) then - call addmap_from(compwav, 'Sw_z0', compatm, mapnstod_consf, 'one', 'unset') + call addmap_from(compwav, 'Sw_z0', compatm, mapbilnr_nstod, 'one', 'unset') call addmrg_to(compatm, 'Sw_z0', mrg_from=compwav, mrg_fld='Sw_z0', mrg_type='copy') end if end if @@ -453,13 +463,13 @@ subroutine esmFldsExchange_nems(gcomp, phase, rc) ! to ocn: evaporation water flux (custom merge in med_phases_prep_ocn) if (phase == 'advertise') then if (is_local%wrap%comp_present(compatm) .and. is_local%wrap%comp_present(compocn)) then - call addfld_from(compatm, 'Faxa_lat') + call addfld_from(compatm, 'Faxa_evap') call addfld_to(compocn, 'Faxa_evap') end if else if ( fldchk(is_local%wrap%FBexp(compocn) , 'Faxa_evap', rc=rc) .and. & - fldchk(is_local%wrap%FBImp(compatm,compatm), 'Faxa_lat' , rc=rc)) then - call addmap_from(compatm, 'Faxa_lat', compocn, mapconsf_aofrac, 'aofrac', 'unset') + fldchk(is_local%wrap%FBImp(compatm,compatm), 'Faxa_evap' , rc=rc)) then + call addmap_from(compatm, 'Faxa_evap', compocn, mapconsf_aofrac, 'aofrac', 'unset') end if end if else if (trim(coupling_mode) == 'nems_orig_data' .or. trim(coupling_mode) == 'nems_frac_aoflux') then @@ -698,7 +708,7 @@ subroutine esmFldsExchange_nems(gcomp, phase, rc) else if ( fldchk(is_local%wrap%FBexp(compwav) , trim(fldname), rc=rc) .and. & fldchk(is_local%wrap%FBImp(compatm,compatm), trim(fldname), rc=rc)) then - call addmap_from(compatm, trim(fldname), compwav, mapnstod_consf, 'one', 'unset') + call addmap_from(compatm, trim(fldname), compwav, mapbilnr_nstod, 'one', 'unset') call addmrg_to(compwav, trim(fldname), mrg_from=compatm, mrg_fld=trim(fldname), mrg_type='copy') end if end if diff --git a/mediator/med.F90 b/mediator/med.F90 index 56fcb7621..3efc94a6e 100644 --- a/mediator/med.F90 +++ b/mediator/med.F90 @@ -661,7 +661,7 @@ subroutine AdvertiseFields(gcomp, importState, exportState, clock, rc) use NUOPC , only : NUOPC_CompAttributeGet, NUOPC_CompAttributeSet, NUOPC_CompAttributeAdd use esmFlds, only : med_fldlist_init1, med_fld_GetFldInfo, med_fldList_entry_type use med_phases_history_mod, only : med_phases_history_init - use med_methods_mod , only : mediator_checkfornans + use med_methods_mod , only : mediator_checkfornans ! input/output variables type(ESMF_GridComp) :: gcomp @@ -921,7 +921,7 @@ subroutine AdvertiseFields(gcomp, importState, exportState, clock, rc) call NUOPC_CompAttributeGet(gcomp, name="check_for_nans", value=cvalue, isPresent=isPresent, isSet=isSet, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return if(isPresent .and. isSet) then - read(cvalue, *) mediator_checkfornans + read(cvalue, *) mediator_checkfornans else mediator_checkfornans = .false. endif @@ -1804,7 +1804,8 @@ subroutine DataInitialize(gcomp, rc) call esmFldsExchange_cesm(gcomp, phase='initialize', rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return else if (trim(coupling_mode(1:4)) == 'nems') then - call esmFldsExchange_nems(gcomp, phase='initialize', rc=rc) + call esmFldsExchange_nems(gcomp, phase='initialize', rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return else if (trim(coupling_mode) == 'hafs') then call esmFldsExchange_hafs(gcomp, phase='initialize', rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return @@ -1939,14 +1940,12 @@ subroutine DataInitialize(gcomp, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return !---------------------------------------------------------- - ! Initialize ocean albedos (this is needed for cesm and hafs) + ! Initialize ocean albedos !---------------------------------------------------------- - if (trim(coupling_mode(1:5)) /= 'nems_') then - if (is_local%wrap%comp_present(compocn) .or. is_local%wrap%comp_present(compatm)) then - call med_phases_ocnalb_run(gcomp, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - end if + if (is_local%wrap%comp_present(compocn) .or. is_local%wrap%comp_present(compatm)) then + call med_phases_ocnalb_run(gcomp, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return end if !--------------------------------------- diff --git a/mediator/med_io_mod.F90 b/mediator/med_io_mod.F90 index d55ebc724..265a5ddda 100644 --- a/mediator/med_io_mod.F90 +++ b/mediator/med_io_mod.F90 @@ -13,7 +13,7 @@ module med_io_mod use NUOPC , only : NUOPC_FieldDictionaryGetEntry use NUOPC , only : NUOPC_FieldDictionaryHasEntry use pio , only : file_desc_t, iosystem_desc_t - use med_internalstate_mod , only : logunit, med_id + use med_internalstate_mod , only : logunit, med_id, maintask use med_constants_mod , only : dbug_flag => med_constants_dbug_flag use med_methods_mod , only : FB_getFieldN => med_methods_FB_getFieldN use med_methods_mod , only : FB_getFldPtr => med_methods_FB_getFldPtr @@ -75,7 +75,7 @@ module med_io_mod character(*),parameter :: prefix = "med_io_" character(*),parameter :: modName = "(med_io_mod) " character(*),parameter :: version = "cmeps0" - + integer :: pio_iotype integer :: pio_ioformat type(iosystem_desc_t), pointer :: io_subsystem @@ -1737,7 +1737,10 @@ subroutine med_io_read_init_iodesc(FB, name1, pioid, iodesc, rc) deallocate(dof) deallocate(minIndexPTile, maxIndexPTile) - + else + if(maintask) write(logunit,'(a)') trim(subname)//' ERROR: '//trim(name1)//' is not present, aborting ' + call ESMF_LogWrite(trim(subname)//' ERROR: '//trim(name1)//' is not present, aborting ', ESMF_LOGMSG_ERROR) + rc = ESMF_FAILURE end if ! end if rcode check end subroutine med_io_read_init_iodesc diff --git a/mediator/med_map_mod.F90 b/mediator/med_map_mod.F90 index 3ab205bd6..54bcbb154 100644 --- a/mediator/med_map_mod.F90 +++ b/mediator/med_map_mod.F90 @@ -259,7 +259,8 @@ subroutine med_map_RouteHandles_initfrom_esmflds(gcomp, flds_scalar_name, llogun if (chkerr(rc,__LINE__,u_FILE_u)) return if (maintask) then write(logunit,'(a)') trim(subname)//' created field_NormOne for '& - //compname(n1)//'->'//compname(n2)//' with mapping '//trim(mapnames(mapindex)) + //trim(compname(n1))//'->'//trim(compname(n2))//' with mapping '& + //trim(mapnames(mapindex)) end if end if end do ! end of loop over map_indiex mappers diff --git a/mediator/med_phases_history_mod.F90 b/mediator/med_phases_history_mod.F90 index 5f150a4b7..7d59a7fea 100644 --- a/mediator/med_phases_history_mod.F90 +++ b/mediator/med_phases_history_mod.F90 @@ -25,7 +25,7 @@ module med_phases_history_mod use med_io_mod , only : med_io_write, med_io_wopen, med_io_enddef, med_io_close use perf_mod , only : t_startf, t_stopf use pio , only : file_desc_t - + implicit none private diff --git a/mediator/med_phases_ocnalb_mod.F90 b/mediator/med_phases_ocnalb_mod.F90 index a5ef002c7..31bd211f0 100644 --- a/mediator/med_phases_ocnalb_mod.F90 +++ b/mediator/med_phases_ocnalb_mod.F90 @@ -6,13 +6,11 @@ module med_phases_ocnalb_mod use med_utils_mod , only : chkerr => med_utils_chkerr use med_methods_mod , only : FB_diagnose => med_methods_FB_diagnose use med_methods_mod , only : State_GetScalar => med_methods_State_GetScalar - use med_internalstate_mod , only : mapconsf, mapnames, compatm, compocn + use med_internalstate_mod , only : mapconsf, mapnames, compatm, compocn, maintask use perf_mod , only : t_startf, t_stopf -#ifdef CESMCOUPLED use shr_orb_mod , only : shr_orb_cosz, shr_orb_decl use shr_orb_mod , only : shr_orb_params, SHR_ORB_UNDEF_INT, SHR_ORB_UNDEF_REAL use shr_log_mod , only : shr_log_unit -#endif implicit none private @@ -26,11 +24,10 @@ module med_phases_ocnalb_mod !-------------------------------------------------------------------------- ! Private interfaces !-------------------------------------------------------------------------- -#ifdef CESMCOUPLED + private med_phases_ocnalb_init private med_phases_ocnalb_orbital_update private med_phases_ocnalb_orbital_init -#endif !-------------------------------------------------------------------------- ! Private data @@ -47,25 +44,30 @@ module med_phases_ocnalb_mod logical :: created ! has memory been allocated here end type ocnalb_type - ! Conversion from degrees to radians character(*),parameter :: u_FILE_u = & __FILE__ -#ifdef CESMCOUPLED character(len=CL) :: orb_mode ! attribute - orbital mode integer :: orb_iyear ! attribute - orbital year integer :: orb_iyear_align ! attribute - associated with model year real(R8) :: orb_obliq ! attribute - obliquity in degrees real(R8) :: orb_mvelp ! attribute - moving vernal equinox longitude real(R8) :: orb_eccen ! attribute and update- orbital eccentricity -#endif + character(len=*) , parameter :: orb_fixed_year = 'fixed_year' character(len=*) , parameter :: orb_variable_year = 'variable_year' character(len=*) , parameter :: orb_fixed_parameters = 'fixed_parameters' + ! used, reused in module + logical :: flux_albav ! use average dif and dir albedos + logical :: use_nextswcday ! use the scalar field for next time (otherwise, will be set using clock) + logical :: use_min_albedo ! apply minimum value of albedo for direct vis, nir + real(R8) :: min_albedo ! minimum value of albedo for direct vis, nir + real(R8) :: albdif ! 60 deg reference albedo, diffuse + real(R8) :: albdir ! 60 deg reference albedo, direct !=============================================================================== contains !=============================================================================== -#ifdef CESMCOUPLED + subroutine med_phases_ocnalb_init(gcomp, ocnalb, rc) !----------------------------------------------------------------------- @@ -74,11 +76,12 @@ subroutine med_phases_ocnalb_init(gcomp, ocnalb, rc) ! All input field bundles are ASSUMED to be on the ocean grid !----------------------------------------------------------------------- - use ESMF , only : ESMF_LogWrite, ESMF_LOGMSG_INFO, ESMF_SUCCESS, ESMF_FAILURE - use ESMF , only : ESMF_VM, ESMF_VMGet, ESMF_Mesh, ESMF_MeshGet - use ESMF , only : ESMF_GridComp, ESMF_GridCompGet - use ESMF , only : ESMF_FieldBundleGet, ESMF_Field, ESMF_FieldGet - use ESMF , only : operator(==) + use ESMF , only : ESMF_LogWrite, ESMF_LOGMSG_INFO, ESMF_SUCCESS, ESMF_FAILURE + use ESMF , only : ESMF_VM, ESMF_VMGet, ESMF_Mesh, ESMF_MeshGet + use ESMF , only : ESMF_GridComp, ESMF_GridCompGet + use ESMF , only : ESMF_FieldBundleGet, ESMF_Field, ESMF_FieldGet + use NUOPC , only : NUOPC_CompAttributeGet + use ESMF , only : operator(==) ! Arguments type(ESMF_GridComp) :: gcomp @@ -97,7 +100,11 @@ subroutine med_phases_ocnalb_init(gcomp, ocnalb, rc) type(InternalState) :: is_local real(R8), pointer :: ownedElemCoords(:) character(len=CL) :: tempc1,tempc2 + character(len=CS) :: cvalue + logical :: use_min_ocnalb + logical :: isPresent, isSet integer :: fieldCount + character(CL) :: msg type(ESMF_Field), pointer :: fieldlist(:) character(*), parameter :: subname = '(med_phases_ocnalb_init) ' !----------------------------------------------------------------------- @@ -186,13 +193,65 @@ subroutine med_phases_ocnalb_init(gcomp, ocnalb, rc) call med_phases_ocnalb_orbital_init(gcomp, logunit, iam==0, rc) if (chkerr(rc,__LINE__,u_FILE_u)) return + ! Determine if reference albedos are used + flux_albav = .false. + call NUOPC_CompAttributeGet(gcomp, name='flux_albav', value=cvalue, isPresent=isPresent, isSet=isSet, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + if (isPresent .and. isSet) then + read(cvalue,*) flux_albav + end if + ! Set reference albedo values + call NUOPC_CompAttributeGet(gcomp, name="albdif", value=cvalue, isPresent=isPresent, isSet=isSet, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (isPresent .and. isSet) then + read(cvalue,*) albdif + else + albdif = 0.06_r8 + end if + call NUOPC_CompAttributeGet(gcomp, name="albdir", value=cvalue, isPresent=isPresent, isSet=isSet, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (isPresent .and. isSet) then + read(cvalue,*) albdir + else + albdir = 0.07_r8 + end if + ! Determine if direct albedo should have a minimum value + call NUOPC_CompAttributeGet(gcomp, name="ocean_albedo_limit", value=cvalue, isPresent=isPresent, isSet=isSet, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (isPresent .and. isSet) then + read(cvalue,*) min_albedo + use_min_albedo = .true. + else + min_albedo = 0.0_R8 + use_min_ocnalb = .false. + endif + ! Allow setting of albedo timestep using the clock instead of the atm's next timestep + use_nextswcday = .true. + call NUOPC_CompAttributeGet(gcomp, name="ScalarFieldIdxNextSwCday", isPresent=isPresent, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (.not. isPresent ) then + use_nextswcday = .false. + endif + + if (flux_albav) then + write(msg,'(2(A,f8.2))') trim(subname)//': mean albedos set: albdif = ',albdif,', albdir = ',albdir + call ESMF_LogWrite(trim(msg), ESMF_LOGMSG_INFO) + else + if (use_min_albedo) then + write(msg,'(A,f8.2)') trim(subname)//': min_albedo setting = ',min_albedo + call ESMF_LogWrite(trim(msg), ESMF_LOGMSG_INFO) + end if + end if + write(msg,'(A,l)') trim(subname)//': use_nextswcday setting is ',use_nextswcday + call ESMF_LogWrite(trim(msg), ESMF_LOGMSG_INFO) + if (dbug_flag > 5) then call ESMF_LogWrite(trim(subname)//": done", ESMF_LOGMSG_INFO) endif call t_stopf('MED:'//subname) end subroutine med_phases_ocnalb_init -#endif + !=============================================================================== subroutine med_phases_ocnalb_run(gcomp, rc) @@ -201,8 +260,10 @@ subroutine med_phases_ocnalb_run(gcomp, rc) ! Compute ocean albedos (on the ocean grid) !----------------------------------------------------------------------- + use NUOPC_Mediator, only : NUOPC_MediatorGet use ESMF , only : ESMF_GridComp, ESMF_GridCompGet, ESMF_TimeInterval use ESMF , only : ESMF_Clock, ESMF_ClockGet, ESMF_Time, ESMF_TimeGet + use ESMF , only : ESMF_ClockIsCreated, ESMF_ClockGetNextTime use ESMF , only : ESMF_VM, ESMF_VMGet use ESMF , only : ESMF_LogWrite, ESMF_LogFoundError use ESMF , only : ESMF_SUCCESS, ESMF_FAILURE, ESMF_LOGMSG_INFO @@ -211,11 +272,11 @@ subroutine med_phases_ocnalb_run(gcomp, rc) use ESMF , only : operator(+) use NUOPC , only : NUOPC_CompAttributeGet use med_constants_mod , only : shr_const_pi + use med_phases_history_mod, only : med_phases_history_write_med ! input/output variables type(ESMF_GridComp) :: gcomp integer, intent(out) :: rc -#ifdef CESMCOUPLED ! local variables type(ocnalb_type), save :: ocnalb type(ESMF_VM) :: vm @@ -224,12 +285,13 @@ subroutine med_phases_ocnalb_run(gcomp, rc) logical :: update_alb type(InternalState) :: is_local type(ESMF_Clock) :: clock + type(ESMF_Clock) :: dclock type(ESMF_Time) :: currTime + type(ESMF_Time) :: nextTime type(ESMF_TimeInterval) :: timeStep character(CL) :: cvalue character(CS) :: starttype ! config start type character(CL) :: runtype ! initial, continue, hybrid, branch - logical :: flux_albav ! flux avg option real(R8) :: nextsw_cday ! calendar day of next atm shortwave real(R8), pointer :: ofrac(:) real(R8), pointer :: ofrad(:) @@ -246,21 +308,13 @@ subroutine med_phases_ocnalb_run(gcomp, rc) real(R8) :: obliqr ! Earth orbit real(R8) :: delta ! Solar declination angle in radians real(R8) :: eccf ! Earth orbit eccentricity factor - real(R8), parameter :: albdif = 0.06_r8 ! 60 deg reference albedo, diffuse - real(R8), parameter :: albdir = 0.07_r8 ! 60 deg reference albedo, direct real(R8), parameter :: const_deg2rad = shr_const_pi/180.0_R8 ! deg to rads character(CL) :: msg logical :: first_call = .true. character(len=*) , parameter :: subname='(med_phases_ocnalb_run)' !--------------------------------------- -#endif - rc = ESMF_SUCCESS - -#ifndef CESMCOUPLED - RETURN ! the following code is not executed unless the model is CESM - -#else + rc = ESMF_SUCCESS ! Determine main task call ESMF_GridCompGet(gcomp, vm=vm, rc=rc) @@ -275,8 +329,7 @@ subroutine med_phases_ocnalb_run(gcomp, rc) ! Determine if ocnalb data type will be initialized - and if not return if (first_call) then - if ( ESMF_FieldBundleIsCreated(is_local%wrap%FBMed_aoflux_a, rc=rc) .and. & - ESMF_FieldBundleIsCreated(is_local%wrap%FBMed_aoflux_o, rc=rc)) then + if (ESMF_FieldBundleIsCreated(is_local%wrap%FBMed_ocnalb_o, rc=rc)) then ocnalb%created = .true. else ocnalb%created = .false. @@ -331,6 +384,26 @@ subroutine med_phases_ocnalb_run(gcomp, rc) call ESMF_TimeGet( currTime, dayOfYear_r8=nextsw_cday, rc=rc ) if (chkerr(rc,__LINE__,u_FILE_u)) return else + ! obtain nextsw_cday from atm if it is in the import state + if (use_nextswcday) then + call State_GetScalar(& + state=is_local%wrap%NstateImp(compatm), & + flds_scalar_name=is_local%wrap%flds_scalar_name, & + flds_scalar_num=is_local%wrap%flds_scalar_num, & + scalar_id=is_local%wrap%flds_scalar_index_nextsw_cday, & + scalar_value=nextsw_cday, rc=rc) + if (chkerr(rc,__LINE__,u_FILE_u)) return + else + call ESMF_TimeGet( currTime, dayOfYear_r8=nextsw_cday, rc=rc ) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + end if + end if + + first_call = .false. + + else + ! Note that med_methods_State_GetScalar includes a broadcast to all other pets + if (use_nextswcday) then call State_GetScalar(& state=is_local%wrap%NstateImp(compatm), & flds_scalar_name=is_local%wrap%flds_scalar_name, & @@ -338,27 +411,14 @@ subroutine med_phases_ocnalb_run(gcomp, rc) scalar_id=is_local%wrap%flds_scalar_index_nextsw_cday, & scalar_value=nextsw_cday, rc=rc) if (chkerr(rc,__LINE__,u_FILE_u)) return + else + call ESMF_ClockGetNextTime(clock, nextTime, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + call ESMF_TimeGet(nextTime, dayOfYear_r8=nextsw_cday, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return end if - - first_call = .false. - - else - - ! Note that med_methods_State_GetScalar includes a broadcast to all other pets - call State_GetScalar(& - state=is_local%wrap%NstateImp(compatm), & - flds_scalar_name=is_local%wrap%flds_scalar_name, & - flds_scalar_num=is_local%wrap%flds_scalar_num, & - scalar_id=is_local%wrap%flds_scalar_index_nextsw_cday, & - scalar_value=nextsw_cday, rc=rc) - if (chkerr(rc,__LINE__,u_FILE_u)) return - end if - call NUOPC_CompAttributeGet(gcomp, name='flux_albav', value=cvalue, rc=rc) - if (chkerr(rc,__LINE__,u_FILE_u)) return - read(cvalue,*) flux_albav - ! Get orbital values call med_phases_ocnalb_orbital_update(clock, logunit, iam==0, eccen, obliqr, lambm0, mvelpp, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return @@ -393,6 +453,9 @@ subroutine med_phases_ocnalb_run(gcomp, rc) ocnalb%anidr(n) = (.026_r8/(cosz**1.7_r8 + 0.065_r8)) + & (.150_r8*(cosz - 0.100_r8 ) * & (cosz - 0.500_r8 ) * (cosz - 1.000_r8 ) ) + if (use_min_albedo) then + ocnalb%anidr(n) = max (ocnalb%anidr(n), min_albedo) + end if ocnalb%avsdr(n) = ocnalb%anidr(n) ocnalb%anidf(n) = albdif ocnalb%avsdf(n) = albdif @@ -430,18 +493,25 @@ subroutine med_phases_ocnalb_run(gcomp, rc) ofrad(:) = ofrac(:) endif + if (ESMF_FieldBundleIsCreated(is_local%wrap%FBMed_ocnalb_o, rc=rc)) then + call NUOPC_MediatorGet(gcomp, driverClock=dClock, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (ESMF_ClockIsCreated(dclock)) then + call med_phases_history_write_med(gcomp, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + end if + end if + if (dbug_flag > 1) then call FB_diagnose(is_local%wrap%FBMed_ocnalb_o, string=trim(subname)//' FBMed_ocnalb_o', rc=rc) if (chkerr(rc,__LINE__,u_FILE_u)) return end if call t_stopf('MED:'//subname) -#endif - end subroutine med_phases_ocnalb_run !=============================================================================== -#ifdef CESMCOUPLED + subroutine med_phases_ocnalb_orbital_init(gcomp, logunit, maintask, rc) !---------------------------------------------------------- @@ -601,7 +671,6 @@ subroutine med_phases_ocnalb_orbital_update(clock, logunit, maintask, eccen, ob endif end subroutine med_phases_ocnalb_orbital_update -#endif !=============================================================================== diff --git a/mediator/med_phases_prep_ocn_mod.F90 b/mediator/med_phases_prep_ocn_mod.F90 index 373d92469..7a71f7e90 100644 --- a/mediator/med_phases_prep_ocn_mod.F90 +++ b/mediator/med_phases_prep_ocn_mod.F90 @@ -31,7 +31,7 @@ module med_phases_prep_ocn_mod public :: med_phases_prep_ocn_accum ! called from run sequence public :: med_phases_prep_ocn_avg ! called from run sequence - private :: med_phases_prep_ocn_custom_cesm + private :: med_phases_prep_ocn_custom private :: med_phases_prep_ocn_custom_nems character(*), parameter :: u_FILE_u = & @@ -217,10 +217,9 @@ subroutine med_phases_prep_ocn_accum(gcomp, rc) end if ! custom merges to ocean - if (trim(coupling_mode) == 'cesm') then - call med_phases_prep_ocn_custom_cesm(gcomp, rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - else if (trim(coupling_mode(1:5)) == 'nems_') then + call med_phases_prep_ocn_custom(gcomp, rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + if (trim(coupling_mode(1:5)) == 'nems_') then call med_phases_prep_ocn_custom_nems(gcomp, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return end if @@ -315,7 +314,7 @@ subroutine med_phases_prep_ocn_avg(gcomp, rc) end subroutine med_phases_prep_ocn_avg !----------------------------------------------------------------------------- - subroutine med_phases_prep_ocn_custom_cesm(gcomp, rc) + subroutine med_phases_prep_ocn_custom(gcomp, rc) !--------------------------------------- ! custom calculations for cesm @@ -372,7 +371,7 @@ subroutine med_phases_prep_ocn_custom_cesm(gcomp, rc) integer :: lsize real(R8) :: c1,c2,c3,c4 character(len=64), allocatable :: fldnames(:) - character(len=*), parameter :: subname='(med_phases_prep_ocn_custom_cesm)' + character(len=*), parameter :: subname='(med_phases_prep_ocn_custom)' !--------------------------------------- rc = ESMF_SUCCESS @@ -620,7 +619,7 @@ subroutine med_phases_prep_ocn_custom_cesm(gcomp, rc) end if call t_stopf('MED:'//subname) - end subroutine med_phases_prep_ocn_custom_cesm + end subroutine med_phases_prep_ocn_custom !----------------------------------------------------------------------------- subroutine med_phases_prep_ocn_custom_nems(gcomp, rc) @@ -643,7 +642,6 @@ subroutine med_phases_prep_ocn_custom_nems(gcomp, rc) real(R8), pointer :: ifrac(:) real(R8), pointer :: ofrac(:) integer :: lsize - real(R8) , parameter :: const_lhvap = 2.501e6_R8 ! latent heat of evaporation ~ J/kg character(len=*), parameter :: subname='(med_phases_prep_ocn_custom_nems)' !--------------------------------------- @@ -672,9 +670,9 @@ subroutine med_phases_prep_ocn_custom_nems(gcomp, rc) if (trim(coupling_mode) == 'nems_orig' .or. & trim(coupling_mode) == 'nems_frac' .or. & trim(coupling_mode) == 'nems_frac_aoflux_sbs') then - customwgt(:) = -ofrac(:) / const_lhvap + customwgt(:) = -ofrac(:) call med_merge_field(is_local%wrap%FBExp(compocn), 'Faxa_evap', & - FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_lat' , wgtA=customwgt, rc=rc) + FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_evap' , wgtA=customwgt, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return customwgt(:) = -ofrac(:) @@ -693,25 +691,6 @@ subroutine med_phases_prep_ocn_custom_nems(gcomp, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return end if - ! netsw_for_ocn = [downsw_from_atm*(1-ice_fraction)*(1-ocn_albedo)] + [pensw_from_ice*(ice_fraction)] - customwgt(:) = ofrac(:) * (1.0_R8 - 0.06_R8) - call med_merge_field(is_local%wrap%FBExp(compocn), 'Foxx_swnet_vdr', & - FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_swvdr' , wgtA=customwgt, & - FBinB=is_local%wrap%FBImp(compice,compocn), fnameB='Fioi_swpen_vdr', wgtB=ifrac, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - call med_merge_field(is_local%wrap%FBExp(compocn), 'Foxx_swnet_vdf', & - FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_swvdf' , wgtA=customwgt, & - FBinB=is_local%wrap%FBImp(compice,compocn), fnameB='Fioi_swpen_vdf', wgtB=ifrac, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - call med_merge_field(is_local%wrap%FBExp(compocn), 'Foxx_swnet_idr', & - FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_swndr' , wgtA=customwgt, & - FBinB=is_local%wrap%FBImp(compice,compocn), fnameB='Fioi_swpen_idr', wgtB=ifrac, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - call med_merge_field(is_local%wrap%FBExp(compocn), 'Foxx_swnet_idf', & - FBinA=is_local%wrap%FBImp(compatm,compocn), fnameA='Faxa_swndf' , wgtA=customwgt, & - FBinB=is_local%wrap%FBImp(compice,compocn), fnameB='Fioi_swpen_idf', wgtB=ifrac, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - deallocate(customwgt) if (dbug_flag > 20) then diff --git a/ufs/ccpp/data/MED_typedefs.F90 b/ufs/ccpp/data/MED_typedefs.F90 index 1b2ce51c5..e5e1b494f 100644 --- a/ufs/ccpp/data/MED_typedefs.F90 +++ b/ufs/ccpp/data/MED_typedefs.F90 @@ -4,9 +4,9 @@ module MED_typedefs !! \htmlinclude MED_typedefs.html !! use machine, only: kind_phys - use physcons, only: con_hvap, con_cp, con_rd, con_eps + use physcons, only: con_hvap, con_cp, con_rd, con_eps, con_rocp use physcons, only: con_epsm1, con_fvirt, con_g - use physcons, only: con_tice + use physcons, only: con_tice, karman implicit none @@ -68,7 +68,9 @@ module MED_typedefs real(kind=kind_phys), pointer :: fm10_water(:) => null() !< Monin-Obukhov similarity parameter for momentum at 10m over water real(kind=kind_phys), pointer :: prslki(:) => null() !< Exner function ratio bt midlayer and interface at 1st layer logical, pointer :: wet(:) => null() !< flag indicating presence of some ocean or lake surface area fraction - logical, pointer :: use_flake(:) => null() !< flag indicating lake points using flake model + integer, pointer :: use_lake_model(:)=>null() !< 0 for points that don't use a lake model, lkm for points that do + real (kind=kind_phys),pointer :: lake_t2m (:) => null() !< 2 meter temperature from CLM Lake model + real (kind=kind_phys),pointer :: lake_q2m (:) => null() !< 2 meter humidity from CLM Lake model real(kind=kind_phys), pointer :: wind(:) => null() !< wind speed at lowest model level (m/s) logical, pointer :: flag_iter(:) => null() !< flag for iteration real(kind=kind_phys), pointer :: qss_water(:) => null() !< surface air saturation specific humidity over water (kg/kg) @@ -172,7 +174,7 @@ module MED_typedefs integer :: sfc_z0_type !< surface roughness options over water logical :: thsfc_loc !< flag for reference pressure in theta calculation integer :: nstf_name(5) !< NSSTM flag: off/uncoupled/coupled=0/1/2 - integer :: lkm !< flag for flake model + integer :: lkm !< 0 = no lake model, 1 = lake model, 2 = lake & nsst on lake points logical :: first_time_step !< flag signaling first time step for time integration routine logical :: frac_grid !< flag for fractional grid logical :: cplwav2atm !< default no wav->atm coupling @@ -189,6 +191,16 @@ module MED_typedefs integer :: lsoil !< number of soil layers integer :: kice !< vertical loop extent for ice levels, start at 1 integer :: lsm_ruc !< flag for RUC land surface model + + ! Lake variables + logical :: frac_ice = .false. !< flag for fractional ice when fractional grid is not in use + logical :: use_lake2m = .false. !< use 2m T & Q calculated by the lake model + integer :: iopt_lake = 1 !< =1 flake, =2 clm lake + integer :: iopt_lake_flake = 1 + integer :: iopt_lake_clm = 2 + + logical :: diag_flux !< flag for flux method of 2-m diagnostics + logical :: diag_log !< flag for log 2-m diagnostics contains procedure :: init => control_initialize end type MED_control_type @@ -208,6 +220,8 @@ module MED_typedefs !! type MED_grid_type real(kind=kind_phys), pointer :: area(:) => null() !< area of the grid cell + real(kind=kind_phys), pointer :: xlat_d(:) => null() !< latitude in degrees + real(kind=kind_phys), pointer :: xlon_d(:) => null() !< longtitude in degrees contains procedure :: create => grid_create !< allocate array data end type MED_grid_type @@ -259,6 +273,7 @@ module MED_typedefs type MED_diag_type real(kind=kind_phys), pointer :: chh(:) => null() !< thermal exchange coefficient (kg m-2 s-1) real(kind=kind_phys), pointer :: cmm(:) => null() !< momentum exchange coefficient (m/s) + real(kind=kind_phys), pointer :: dpt2m(:) => null() !< 2-m dewpoint (K) contains procedure :: create => diag_create !< allocate array data end type MED_diag_type @@ -343,8 +358,12 @@ subroutine interstitial_create(interstitial, im) interstitial%prslki = clear_val allocate(interstitial%wet(im)) interstitial%wet = .false. - allocate(interstitial%use_flake(im)) - interstitial%use_flake = .false. + allocate(interstitial%use_lake_model(im)) + interstitial%use_lake_model = 0 + allocate(interstitial%lake_t2m(im)) + interstitial%lake_t2m=-9999 + allocate(interstitial%lake_q2m(im)) + interstitial%lake_q2m=-9999 allocate(interstitial%wind(im)) interstitial%wind = huge allocate(interstitial%flag_iter(im)) @@ -591,7 +610,9 @@ subroutine interstitial_phys_reset(interstitial) interstitial%tsurf_ice = huge interstitial%tsurf_land = huge interstitial%tsurf_water = huge - interstitial%use_flake = .false. + interstitial%use_lake_model = 0 + interstitial%lake_t2m = -9999 + interstitial%lake_q2m = -9999 interstitial%uustar_ice = huge interstitial%uustar_land = huge interstitial%uustar_water = huge @@ -636,6 +657,13 @@ subroutine control_initialize(model) model%lsoil = 4 model%kice = 2 model%lsm_ruc = 3 + model%frac_ice = .false. + model%use_lake2m = .false. + model%iopt_lake = 1 + model%iopt_lake_flake = 1 + model%iopt_lake_clm = 2 + model%diag_flux = .false. + model%diag_log = .false. end subroutine control_initialize @@ -658,6 +686,10 @@ subroutine grid_create(grid, im) allocate(grid%area(im)) grid%area = clear_val + allocate(grid%xlat_d(im)) + grid%xlat_d = clear_val + allocate(grid%xlon_d(im)) + grid%xlon_d = clear_val end subroutine grid_create @@ -745,6 +777,8 @@ subroutine diag_create(diag, im) diag%chh = clear_val allocate(diag%cmm(im)) diag%cmm = clear_val + allocate(diag%dpt2m(im)) + diag%dpt2m = clear_val end subroutine diag_create diff --git a/ufs/ccpp/data/MED_typedefs.meta b/ufs/ccpp/data/MED_typedefs.meta index 8177ae5ca..271110e9c 100644 --- a/ufs/ccpp/data/MED_typedefs.meta +++ b/ufs/ccpp/data/MED_typedefs.meta @@ -202,12 +202,28 @@ units = flag dimensions = (horizontal_loop_extent) type = logical -[use_flake] - standard_name = flag_for_using_flake - long_name = flag indicating lake points using flake model +[lake_t2m] + standard_name = temperature_at_2m_from_clm_lake + long_name = temperature at 2m from clm lake + units = K + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + active = (control_for_lake_model_selection == 2) +[lake_q2m] + standard_name = specific_humidity_at_2m_from_clm_lake + long_name = specific humidity at 2m from clm lake + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + active = (control_for_lake_model_selection == 2) +[use_lake_model] + standard_name = flag_for_using_lake_model + long_name = flag indicating lake points using a lake model units = flag dimensions = (horizontal_loop_extent) - type = logical + type = integer [wind] standard_name = wind_speed_at_lowest_model_layer long_name = wind speed at lowest model level @@ -817,9 +833,33 @@ units = flag dimensions = () type = integer +[iopt_lake] + standard_name = control_for_lake_model_selection + long_name = control for lake model selection + units = 1 + dimensions = () + type = integer +[iopt_lake_flake] + standard_name = flake_model_control_selection_value + long_name = value that indicates flake model in the control for lake model selection + units = 1 + dimensions = () + type = integer +[iopt_lake_clm] + standard_name = clm_lake_model_control_selection_value + long_name = value that indicates clm lake model in the control for lake model selection + units = 1 + dimensions = () + type = integer +[use_lake2m] + standard_name = use_2m_diagnostics_calculated_by_lake_model + long_name = model 2m diagnostics use the temperature and humidity calculated by the lake model + units = flag + dimensions = () + type = integer [lkm] - standard_name = control_for_lake_surface_scheme - long_name = flag for lake surface model + standard_name = control_for_lake_model_execution_method + long_name = control for lake model execution: 0=no lake, 1=lake, 2=lake+nsst units = flag dimensions = () type = integer @@ -835,6 +875,12 @@ units = flag dimensions = () type = logical +[frac_ice] + standard_name = flag_for_fractional_ice_when_fractional_landmask_is_disabled + long_name = flag for fractional ice when fractional landmask is disabled + units = flag + dimensions = () + type = logical [cplwav2atm] standard_name = flag_for_one_way_ocean_wave_coupling_to_atmosphere long_name = flag controlling ocean wave coupling to the atmosphere (default off) @@ -924,6 +970,18 @@ units = flag dimensions = () type = integer +[diag_flux] + standard_name = flag_for_flux_method_in_2m_diagnostics + long_name = flag for flux method in 2-m diagnostics + units = flag + dimensions = () + type = logical +[diag_log] + standard_name = flag_for_log_method_in_2m_diagnostics + long_name = flag for log method in 2-m diagnostics + units = flag + dimensions = () + type = logical ######################################################################## [ccpp-table-properties] @@ -964,6 +1022,20 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys +[xlat_d] + standard_name = latitude_in_degree + long_name = latitude in degree north + units = degree_north + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys +[xlon_d] + standard_name = longitude_in_degree + long_name = longitude in degree east + units = degree_east + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys ######################################################################## [ccpp-table-properties] @@ -1228,6 +1300,13 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys +[dpt2m] + standard_name = dewpoint_temperature_at_2m + long_name = 2 meter dewpoint temperature + units = K + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys ######################################################################## [ccpp-table-properties] @@ -1343,3 +1422,17 @@ dimensions = () type = real kind = kind_phys +[con_rocp] + standard_name = ratio_of_gas_constant_dry_air_to_specific_heat_of_dry_air_at_constant_pressure + long_name = (rd/cp) + units = none + dimensions = () + type = real + kind = kind_phys +[karman] + standard_name = von_karman_constant + long_name = von karman constant + units = none + dimensions = () + type = real +