From d58a72780827ef3dd0a5df1d677836c3f9f369da Mon Sep 17 00:00:00 2001 From: Soren Rasmussen Date: Wed, 1 May 2024 14:19:25 -0600 Subject: [PATCH 1/2] Adding w3emc library as a module within the project and replacing calls to exterior w3emc dependency with local subroutines from w3emc module. Removing sp dependency. --- CMakeLists.txt | 4 +- physics/GWD/cires_tauamf_data.F90 | 2 +- physics/GWD/cires_ugwpv1_module.F90 | 10 +- .../GFS_phys_time_vary.scm.F90 | 2 + .../UFS_SCM_NEPTUNE/GFS_time_vary_pre.scm.F90 | 4 +- .../UFS_SCM_NEPTUNE/iccninterp.F90 | 2 + .../Interstitials/UFS_SCM_NEPTUNE/sfcsub.F | 2 + physics/MP/Morrison_Gettelman/aerinterp.F90 | 1 + physics/Radiation/radiation_astronomy.f | 16 +- .../SFC_Layer/UFS/module_nst_water_prop.f90 | 2 +- physics/photochem/h2ointerp.f90 | 2 +- physics/tools/w3emc.F90 | 343 ++++++++++++++++++ 12 files changed, 370 insertions(+), 20 deletions(-) create mode 100644 physics/tools/w3emc.F90 diff --git a/CMakeLists.txt b/CMakeLists.txt index c56070123..6112967dd 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -183,9 +183,7 @@ set_target_properties(ccpp_physics PROPERTIES VERSION ${PROJECT_VERSION} target_include_directories(ccpp_physics PUBLIC $) -target_link_libraries(ccpp_physics PUBLIC w3emc::w3emc_d - sp::sp_d - NetCDF::NetCDF_Fortran) +target_link_libraries(ccpp_physics PUBLIC NetCDF::NetCDF_Fortran) # Define where to install the library install(TARGETS ccpp_physics diff --git a/physics/GWD/cires_tauamf_data.F90 b/physics/GWD/cires_tauamf_data.F90 index 4f12b2ec1..f5fbf256d 100644 --- a/physics/GWD/cires_tauamf_data.F90 +++ b/physics/GWD/cires_tauamf_data.F90 @@ -1,5 +1,5 @@ module cires_tauamf_data - + use w3emc, only: w3doxdat, w3kind, w3movdat use machine, only: kind_phys !........................................................................................... ! tabulated GW-sources: GRACILE/Ern et al., 2018 and/or Resolved GWs from C384-Annual run diff --git a/physics/GWD/cires_ugwpv1_module.F90 b/physics/GWD/cires_ugwpv1_module.F90 index 9c3fa24ee..809abc84d 100644 --- a/physics/GWD/cires_ugwpv1_module.F90 +++ b/physics/GWD/cires_ugwpv1_module.F90 @@ -11,10 +11,11 @@ module cires_ugwpv1_module !................................................................................... ! ! - use machine, only : kind_phys - use ugwp_common, only : arad, pi, pi2, hpscale, rhp, rhp2, rh4, rhp4, khp, hpskm - use ugwp_wmsdis_init, only : ilaunch, nslope, lhmet, lzmax, lzmin, lzstar - use ugwp_wmsdis_init, only : tau_min, tamp_mpa + use machine, only : kind_phys + use ugwp_common, only : arad, pi, pi2, hpscale, rhp, rhp2, rh4, rhp4, khp, hpskm + use ugwp_wmsdis_init, only : ilaunch, nslope, lhmet, lzmax, lzmin, lzstar + use ugwp_wmsdis_init, only : tau_min, tamp_mpa + use w3emc, only : iw3jdn implicit none logical :: module_is_initialized @@ -436,7 +437,6 @@ subroutine calendar_ugwp(yr, mm, dd, ddd_ugwp) integer, intent(in) :: yr, mm, dd integer :: ddd_ugwp - integer :: iw3jdn integer :: jd1, jddd jd1 = iw3jdn(yr,1,1) jddd = iw3jdn(yr,mm,dd) diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.scm.F90 b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.scm.F90 index 6b1b29af1..c8742000f 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.scm.F90 +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.scm.F90 @@ -7,6 +7,8 @@ !> @{ module GFS_phys_time_vary + use w3emc, only: w3doxdat, w3kind, w3movdat + use machine, only : kind_phys, kind_dbl_prec, kind_sngl_prec use mersenne_twister, only: random_setseed, random_number diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.scm.F90 b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.scm.F90 index 3293e09e4..4f708c964 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.scm.F90 +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.scm.F90 @@ -3,6 +3,8 @@ module GFS_time_vary_pre + use w3emc, only: iw3jdn, w3difdat, w3kind + use funcphys, only: gfuncphys implicit none @@ -93,7 +95,7 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr, & real(kind=kind_phys), parameter :: con_hr = 3600.0_kind_phys real(kind=kind_dbl_prec) :: rinc8(5) - integer :: iw3jdn + integer :: w3kindreal,w3kindint integer :: jd0, jd1 real :: fjd diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/iccninterp.F90 b/physics/Interstitials/UFS_SCM_NEPTUNE/iccninterp.F90 index dd752d9b8..fa2a2738e 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/iccninterp.F90 +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/iccninterp.F90 @@ -7,6 +7,8 @@ !! IN and CCN data. module iccninterp + use w3emc, only: w3doxdat, w3kind, w3movdat + implicit none private diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F b/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F index 369a94358..bfd677800 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F @@ -7069,6 +7069,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & &, tile_num_ch, i_index, j_index) ! use machine , only : kind_io8,kind_io4, kind_dbl_prec + use w3emc, only : w3movdat implicit none character(len=*), intent(in) :: tile_num_ch integer, intent(in) :: i_index(len), j_index(len) @@ -8574,6 +8575,7 @@ subroutine fixrda(lugb,fngrib,kpds5,slmask, & &, outlat, outlon, me) use machine , only : kind_io8,kind_dbl_prec,kind_sngl_prec use sfccyc_module, only : mdata + use w3emc, only : w3movdat implicit none integer nrepmx,nvalid,imo,iyr,idy,jret,ihr,nrept,lskip,lugi, & & lgrib,j,ndata,i,inttyp,jmax,imax,ijmax,ij,jday,len,iret, & diff --git a/physics/MP/Morrison_Gettelman/aerinterp.F90 b/physics/MP/Morrison_Gettelman/aerinterp.F90 index 174a1a1a1..d7c8f9a87 100644 --- a/physics/MP/Morrison_Gettelman/aerinterp.F90 +++ b/physics/MP/Morrison_Gettelman/aerinterp.F90 @@ -6,6 +6,7 @@ !! This module contain subroutines of reading and interpolating !! aerosol data for MG microphysics. module aerinterp + use w3emc, only: w3doxdat, w3kind, w3movdat implicit none diff --git a/physics/Radiation/radiation_astronomy.f b/physics/Radiation/radiation_astronomy.f index b25c89a8c..a21040938 100644 --- a/physics/Radiation/radiation_astronomy.f +++ b/physics/Radiation/radiation_astronomy.f @@ -75,9 +75,9 @@ !> \brief This module sets up astronomical quantities for solar radiation !! calculations. !! -!! Operational GFS selection for Solar constant value -!! (namelist control parameter - \b ISOL = 2) -!! \n ISOL=0: presribed value = 1366 \f$W m^{-2}\f$ (old) +!! Operational GFS selection for Solar constant value +!! (namelist control parameter - \b ISOL = 2) +!! \n ISOL=0: presribed value = 1366 \f$W m^{-2}\f$ (old) !! \n ISOL=10: prescibed value = 1361 \f$W m^{-2}\f$ (new) !! \n ISOL=1: NOAA old yearly solar constant table with 11-year cycle (range: 1944-2006) !! \n ISOL=2: NOAA new yearly solar constant table with 11-year cycle (range: 1850-2019) @@ -86,10 +86,11 @@ !! \version NCEP-Radiation_astronomy v5.2 Jan 2013 !> This module sets up astronomy quantities for solar radiation calculations. - module module_radiation_astronomy + module module_radiation_astronomy ! - use machine, only : kind_phys + use machine, only : kind_phys use module_iounitdef, only : NIRADSF + use w3emc, only : iw3jdn, w3fs26 ! implicit none ! @@ -196,7 +197,7 @@ subroutine sol_init & degrad = 180.0/con_pi tpi = 2.0 * con_pi hpi = 0.5 * con_pi - pid12 = con_pi/f12 + pid12 = con_pi/f12 ! --- initialization isolflg = isolar @@ -400,7 +401,6 @@ subroutine sol_update & real (kind=kind_phys) :: fjd, fjd1, dlt, r1, alp integer :: jd, jd1, iyear, imon, iday, ihr, imin, isec - integer :: iw3jdn integer :: i, iyr, iyr1, iyr2, jyr, nn, nswr, icy1, icy2, icy logical :: file_exist @@ -901,7 +901,7 @@ subroutine coszmn & coszdg(i) = coszen(i) * rstp if (istsun(i) > 0 .and. coszen(i) /= 0.0_kind_phys) then coszen(i) = coszen(i) / istsun(i) - endif + endif enddo ! return diff --git a/physics/SFC_Layer/UFS/module_nst_water_prop.f90 b/physics/SFC_Layer/UFS/module_nst_water_prop.f90 index 858659e90..edb798c37 100644 --- a/physics/SFC_Layer/UFS/module_nst_water_prop.f90 +++ b/physics/SFC_Layer/UFS/module_nst_water_prop.f90 @@ -8,6 +8,7 @@ module module_nst_water_prop use machine , only : kind_phys use module_nst_parameters , only : t0k, zero, one, half + use w3emc, only : iw3jdn implicit none ! @@ -544,7 +545,6 @@ subroutine compjd(jyr,jmnth,jday,jhr,jmn,jd,fjd) !$$$ ! integer :: jyr,jmnth,jday,jhr,jmn,jd - integer :: iw3jdn real (kind=kind_phys) fjd jd=iw3jdn(jyr,jmnth,jday) if(jhr.lt.12) then diff --git a/physics/photochem/h2ointerp.f90 b/physics/photochem/h2ointerp.f90 index f5a1f36c6..711c1d54a 100644 --- a/physics/photochem/h2ointerp.f90 +++ b/physics/photochem/h2ointerp.f90 @@ -6,7 +6,7 @@ !> This module contains subroutines of reading and interpolating !! h2o coefficients. module h2ointerp - + use w3emc, only: w3doxdat, w3kind, w3movdat implicit none private diff --git a/physics/tools/w3emc.F90 b/physics/tools/w3emc.F90 new file mode 100644 index 000000000..58bcf23a3 --- /dev/null +++ b/physics/tools/w3emc.F90 @@ -0,0 +1,343 @@ +module w3emc + use machine, only: kind_sngl_prec, kind_dbl_prec + implicit none + public :: w3fs26 + + interface w3difdat + module procedure :: w3difdat32 + module procedure :: w3difdat64 + end interface w3difdat + + interface w3movdat + module procedure :: w3movdat32 + module procedure :: w3movdat64 + end interface w3movdat + + interface w3reddat + module procedure :: w3reddat32 + module procedure :: w3reddat64 + end interface w3reddat + +contains + + ! @brief Computes julian day number from year (4 digits), month, and day. + ! @author Ralph Jones @date 1987-03-29 + ! Computes julian day number from year (4 digits), month, + ! and day. iw3jdn is valid for years 1583 a.d. to 3300 a.d. + ! Julian day number can be used to compute day of week, day of + ! year, record numbers in an archive, replace day of century, + ! find the number of days between two dates. + ! @param[in] IYEAR Integer year (4 Digits) + ! @param[in] MONTH Integer month of year (1 - 12) + ! @param[in] IDAY Integer day of month (1 - 31) + ! @return IW3JDN Integer Julian day number + ! - Jan 1, 1960 is Julian day number 2436935 + ! - Jan 1, 1987 is Julian day number 2446797 + + ! @note Julian period was devised by joseph scaliger in 1582. + ! Julian day number #1 started on Jan. 1,4713 B.C. Three major + ! chronological cycles begin on the same day. A 28-year solar + ! cycle, a 19-year luner cycle, a 15-year indiction cycle, used + ! in ancient rome to regulate taxes. It will take 7980 years + ! to complete the period, the product of 28, 19, and 15. + ! scaliger named the period, date, and number after his father + ! Julius (not after the julian calendar). This seems to have + ! caused a lot of confusion in text books. Scaliger name is + ! spelled three different ways. Julian date and Julian day + ! number are interchanged. A Julian date is used by astronomers + ! to compute accurate time, it has a fraction. When truncated to + ! an integer it is called an Julian day number. This function + ! was in a letter to the editor of the communications of the acm + ! volume 11 / number 10 / october 1968. The Julian day number + ! can be converted to a year, month, day, day of week, day of + ! year by calling subroutine w3fs26. + ! @author Ralph Jones @date 1987-03-29 + function iw3jdn(iyear, month, iday) + integer, intent(in) :: iyear, month, iday + integer :: iw3jdn + iw3jdn = iday - 32075 & + + 1461 * (iyear + 4800 + (month - 14) / 12) / 4 & + + 367 * (month - 2 - (month -14) / 12 * 12) / 12 & + - 3 * ((iyear + 4900 + (month - 14) / 12) / 100) / 4 + end function iw3jdn + + + + !> @brief Return the real kind and integer kind used in w3 lib. + !> @author Jun Wang @date 2011-06-24 + !> This subprogram returns the real kind and the integer kind that the w3 lib + !> is compiled with. + !> @param[out] KINDREAL Kind of real number in w3 lib + !> @param[out] KINDINT Kind of integer number in w3 lib + !> + !> @author Jun Wang @date 2011-06-24 + subroutine w3kind(kindreal, kindint) + implicit none + integer,intent(out) :: kindreal,kindint + ! get real kind from a real number + kindreal=kind(1.0) + kindint=kind(1) + end subroutine w3kind + + + !> @brief Returns the integer day of week, the day + !> of year, and julian day given an NCEP absolute date and time. + !> @author Mark Iredell @date 1998-01-05 + + !> @param[in] IDAT Integer (8) NCEP absolute date and time + !> (year, month, day, time zone, hour, minute, second, millisecond) + !> @param[out] JDOW Integer day of week (1-7, where 1 is sunday) + !> @param[out] JDOY Integer day of year (1-366, where 1 is january 1) + !> @param[out] JDAY Integer julian day (day number from jan. 1,4713 b.c.) + !> + !> @author Mark Iredell @date 1998-01-05 + subroutine w3doxdat(idat, jdow, jdoy, jday) + integer :: idat(8), jdow, jdoy, jday + integer :: jy, jm, jd + + ! get julian day and then get day of week and day of year + jday=iw3jdn(idat(1), idat(2), idat(3)) + call w3fs26(jday, jy, jm, jd, jdow, jdoy) + end subroutine w3doxdat + + + !> @brief Return a time interval between two dates. + !> @author Mark Iredell @date 1998-01-05 + !> Returns the elapsed time interval from + !> an NCEP absolute date and time given in the second argument until + !> an NCEP absolute date and time given in the first argument. + !> The output time interval is in one of seven canonical forms + !> of the ncep relative time interval data structure. + !> @param[in] JDAT Integer (8) ncep absolute date and time + !> (year, month, day, time zone, hour, minute, second, millisecond) + !> @param[in] IDAT Integer (8) ncep absolute date and time + !> (year, month, day, time zone, hour, minute, second, millisecond) + !> @param[in] IT Integer relative time interval format type + !> (-1 for first reduced type (hours always positive), + !> 0 for second reduced type (hours can be negative), + !> 1 for days only, 2 for hours only, 3 for minutes only, + !> 4 for seconds only, 5 for milliseconds only) + !> @param[out] RINC Real (5) ncep relative time interval + !> (days, hours, minutes, seconds, milliseconds) + !> (time interval is positive if jdat is later than idat.) + !> + !> @author Mark Iredell @date 1998-01-05 + subroutine w3difdat32(jdat, idat, it, rinc) + integer :: it + integer :: jdat(8),idat(8) + real(kind_sngl_prec) :: rinc(5) + real(kind_sngl_prec) :: rinc1(5) + + ! difference the days and time and put into canonical form + rinc1(1)=iw3jdn(jdat(1),jdat(2),jdat(3)) - & + iw3jdn(idat(1),idat(2),idat(3)) + rinc1(2:5)=jdat(5:8)-idat(5:8) + call w3reddat(it,rinc1,rinc) + end subroutine w3difdat32 + + subroutine w3difdat64(jdat, idat, it, rinc) + integer :: it + integer :: jdat(8),idat(8) + real(kind_dbl_prec) :: rinc(5) + real(kind_dbl_prec) :: rinc1(5) + + ! difference the days and time and put into canonical form + rinc1(1)=iw3jdn(jdat(1),jdat(2),jdat(3)) - & + iw3jdn(idat(1),idat(2),idat(3)) + rinc1(2:5)=jdat(5:8)-idat(5:8) + call w3reddat(it,rinc1,rinc) + end subroutine w3difdat64 + + + + ! @brief Year, month, day from julian day number + ! @author Ralph Jones @date 1987-03-29 + subroutine w3fs26(jldayn, iyear, month, iday, idaywk, idayyr) + integer :: jldayn, iyear, month, iday, idaywk, idayyr + integer :: i, j, l + real :: n + l = jldayn + 68569 + n = 4 * l / 146097 + l = l - (146097 * n + 3) / 4 + i = 4000 * (l + 1) / 1461001 + l = l - 1461 * i / 4 + 31 + j = 80 * l / 2447 + iday = l - 2447 * j / 80 + l = j / 11 + month = j + 2 - 12 * l + iyear = 100 * (n - 49) + i + l + idaywk = mod((jldayn + 1),7) + 1 + idayyr = jldayn - & + (-31739 +1461 * (iyear+4799) / 4 - 3 * ((iyear+4899)/100)/4) + end subroutine w3fs26 + + !> @file + !> @brief Reduce a time interval to a canonical form. + !> @author Mark Iredell @date 1998-01-05 + !> ### Program History Log: + !> Date | Programmer | Comment + !> -----|------------|-------- + !> 1998-01-05 | Mark Iredell | Initial. + !> + !> @param[in] IT Relative time interval format type + !> - (-1 for first reduced type (hours always positive), + !> - 0 for second reduced type (hours can be negative), + !> - 1 for days only, 2 for hours only, 3 for minutes only, + !> - 4 for seconds only, 5 for milliseconds only) + !> @param[in] RINC NCEP relative time interval (days, hours, minutes, seconds, + !> milliseconds) + !> @param[out] DINC NCEP relative time interval (days, hours, minutes, + !> seconds, milliseconds) + !> + !> @author Mark Iredell @date 1998-01-05 + subroutine w3reddat32(it, rinc, dinc) + integer :: it + real(kind_sngl_prec) :: rinc(5), dinc(5) + ! parameters for number of units in a day + ! and number of milliseconds in a unit + ! and number of next smaller units in a unit, respectively + integer, dimension(5), parameter :: itd=(/1,24,1440,86400,86400000/), & + itm=itd(5)/itd + integer, dimension(4), parameter :: itn=itd(2:5)/itd(1:4) + integer, parameter :: np=16 + integer :: iinc(4), jinc(5), kinc(5) + integer :: ms + + ! first reduce to the first reduced form + iinc=floor(rinc(1:4)) + ! convert all positive fractional parts to milliseconds + ! and determine canonical milliseconds + jinc(5)=nint(dot_product(rinc(1:4)-iinc,real(itm(1:4)))+rinc(5)) + kinc(5)=modulo(jinc(5),itn(4)) + ! convert remainder to seconds and determine canonical seconds + jinc(4)=iinc(4)+(jinc(5)-kinc(5))/itn(4) + kinc(4)=modulo(jinc(4),itn(3)) + ! convert remainder to minutes and determine canonical minutes + jinc(3)=iinc(3)+(jinc(4)-kinc(4))/itn(3) + kinc(3)=modulo(jinc(3),itn(2)) + ! convert remainder to hours and determine canonical hours + jinc(2)=iinc(2)+(jinc(3)-kinc(3))/itn(2) + kinc(2)=modulo(jinc(2),itn(1)) + ! convert remainder to days and compute milliseconds of the day + kinc(1)=iinc(1)+(jinc(2)-kinc(2))/itn(1) + ms=dot_product(kinc(2:5),itm(2:5)) + + ! next reduce to either single value canonical form + ! or to one of the two reduced forms + if(it.ge.1.and.it.le.5) then + ! ensure that exact multiples of 1./np dinc(it)=real(kinc(1))*itd(it)+rp/np + else + ! the reduced form is done except the second reduced form is modified + ! for negative time intervals with fractional days + dinc=kinc + if(it.eq.0.and.kinc(1).lt.0.and.ms.gt.0) then + dinc(1)=dinc(1)+1 + dinc(2:5)=mod(ms-itm(1),itm(1:4))/itm(2:5) + endif + endif + end subroutine w3reddat32 + + subroutine w3reddat64(it, rinc, dinc) + integer :: it + real(kind_dbl_prec) :: rinc(5), dinc(5) + ! parameters for number of units in a day + ! and number of milliseconds in a unit + ! and number of next smaller units in a unit, respectively + integer, dimension(5), parameter :: itd=(/1,24,1440,86400,86400000/), & + itm=itd(5)/itd + integer, dimension(4), parameter :: itn=itd(2:5)/itd(1:4) + integer, parameter :: np=16 + integer :: iinc(4), jinc(5), kinc(5) + integer :: ms + + ! first reduce to the first reduced form + iinc=floor(rinc(1:4)) + ! convert all positive fractional parts to milliseconds + ! and determine canonical milliseconds + jinc(5)=nint(dot_product(rinc(1:4)-iinc,real(itm(1:4)))+rinc(5)) + kinc(5)=modulo(jinc(5),itn(4)) + ! convert remainder to seconds and determine canonical seconds + jinc(4)=iinc(4)+(jinc(5)-kinc(5))/itn(4) + kinc(4)=modulo(jinc(4),itn(3)) + ! convert remainder to minutes and determine canonical minutes + jinc(3)=iinc(3)+(jinc(4)-kinc(4))/itn(3) + kinc(3)=modulo(jinc(3),itn(2)) + ! convert remainder to hours and determine canonical hours + jinc(2)=iinc(2)+(jinc(3)-kinc(3))/itn(2) + kinc(2)=modulo(jinc(2),itn(1)) + ! convert remainder to days and compute milliseconds of the day + kinc(1)=iinc(1)+(jinc(2)-kinc(2))/itn(1) + ms=dot_product(kinc(2:5),itm(2:5)) + + ! next reduce to either single value canonical form + ! or to one of the two reduced forms + if(it.ge.1.and.it.le.5) then + ! ensure that exact multiples of 1./np dinc(it)=real(kinc(1))*itd(it)+rp/np + else + ! the reduced form is done except the second reduced form is modified + ! for negative time intervals with fractional days + dinc=kinc + if(it.eq.0.and.kinc(1).lt.0.and.ms.gt.0) then + dinc(1)=dinc(1)+1 + dinc(2:5)=mod(ms-itm(1),itm(1:4))/itm(2:5) + endif + endif + end subroutine w3reddat64 + + + !> @file + !> @brief Return a date from a time interval and date + !> @author Mark Iredell @date 1998-08-01 + !> This subprogram returns the date and time that is a given + !> NCEP relative time interval from an NCEP absolute date and time. + !> The output is in the NCEP absolute date and time data structure. + !> + !> ### Program History Log: + !> Date | Programmer | Comment + !> -----|------------|-------- + !> 1998-01-05 | Mark Iredell | Initial. + !> + !> @param[in] RINC NCEP relative time interval (days, hours, minutes, seconds + !> milliseconds) + !> @param[in] IDAT NCEP absolute date and time (year, month, day, time zone, + !> hour, minute, second, millisecond) + !> @param[in] JDAT NCEP absolute date and time (year, month, day, time zone, + !> hour, minute, second, millisecond) (jdat is later than idat if time + !> interval is positive.) + !> + !> @author Mark Iredell @date 1998-08-01 + subroutine w3movdat32(rinc, idat, jdat) + real(kind_sngl_prec) :: rinc(5) + integer :: idat(8), jdat(8), jdow, jdoy, jldayn + real(kind_sngl_prec) :: rinc1(5), rinc2(5) + + ! add the interval to the input time of day and put into reduced form + ! and then compute new date using julian day arithmetic. + rinc1(1)=rinc(1) + rinc1(2:5)=rinc(2:5)+idat(5:8) + call w3reddat(-1,rinc1,rinc2) + jldayn=iw3jdn(idat(1),idat(2),idat(3))+nint(rinc2(1)) + call w3fs26(jldayn,jdat(1),jdat(2),jdat(3),jdow,jdoy) + jdat(4)=idat(4) + jdat(5:8)=nint(rinc2(2:5)) + ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + end subroutine w3movdat32 + + subroutine w3movdat64(rinc, idat, jdat) + real(kind_dbl_prec) :: rinc(5) + integer :: idat(8), jdat(8), jdow, jdoy, jldayn + real(kind_dbl_prec) :: rinc1(5), rinc2(5) + + ! add the interval to the input time of day and put into reduced form + ! and then compute new date using julian day arithmetic. + rinc1(1)=rinc(1) + rinc1(2:5)=rinc(2:5)+idat(5:8) + call w3reddat(-1,rinc1,rinc2) + jldayn=iw3jdn(idat(1),idat(2),idat(3))+nint(rinc2(1)) + call w3fs26(jldayn,jdat(1),jdat(2),jdat(3),jdow,jdoy) + jdat(4)=idat(4) + jdat(5:8)=nint(rinc2(2:5)) + ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + end subroutine w3movdat64 + +end module w3emc From 40ac1999e8f9441cd5442a8851976d6f351baede Mon Sep 17 00:00:00 2001 From: Soren Rasmussen Date: Thu, 20 Jun 2024 12:46:41 -0600 Subject: [PATCH 2/2] Adding changes to build UFS with the new local w3emc module. UFS still needs to link in w3emc module for other functions it needs. SP dependency removed, added local module in sfcsub.F to build needed splat routine. --- CMakeLists.txt | 4 +- .../GFS_phys_time_vary.fv3.F90 | 4 +- .../UFS_SCM_NEPTUNE/GFS_time_vary_pre.fv3.F90 | 2 +- .../Interstitials/UFS_SCM_NEPTUNE/sfcsub.F | 468 +++++++++++++++--- 4 files changed, 402 insertions(+), 76 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 6112967dd..cb3b418c9 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -167,8 +167,10 @@ if(CMAKE_BUILD_TYPE STREQUAL "Release" AND ${CMAKE_Fortran_COMPILER_ID} STREQUAL endif() #------------------------------------------------------------------------------ +set(LOCAL_W3EMC + ${LOCAL_CURRENT_SOURCE_DIR}/physics/tools/w3emc.F90) -add_library(ccpp_physics STATIC ${SCHEMES} ${SCHEMES_OPENMP_OFF} ${SCHEMES_DYNAMICS} ${CAPS}) +add_library(ccpp_physics STATIC ${LOCAL_W3EMC} ${SCHEMES} ${SCHEMES_OPENMP_OFF} ${SCHEMES_DYNAMICS} ${CAPS}) # Generate list of Fortran modules from defined sources foreach(source_f90 ${CAPS}) get_filename_component(tmp_source_f90 ${source_f90} NAME) diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.fv3.F90 b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.fv3.F90 index 406887f00..85f9d8c55 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.fv3.F90 +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_phys_time_vary.fv3.F90 @@ -12,6 +12,8 @@ module GFS_phys_time_vary use machine, only : kind_phys, kind_dbl_prec, kind_sngl_prec + use w3emc, only : w3movdat + use mersenne_twister, only: random_setseed, random_number use module_ozphys, only: ty_ozphys @@ -458,7 +460,7 @@ subroutine GFS_phys_time_vary_init ( !$OMP shared(isbarren_table,isice_table,isurban_table) & !$omp shared(iswater_table,laim_table,sla_table,bexp_table) & !$omp shared(stc,smc,slc,tg3,snowxy,tsnoxy,snicexy,snliqxy) & -!$omp shared(zsnsoxy,stype,smcmax_table,smcwlt_table,zs,dzs) & +!$omp shared(zsnsoxy,stype,smcmax_table,smcwlt_table,zs,dzs) & !$omp shared(dwsat_table,dksat_table,psisat_table,smoiseq) & !$OMP shared(smcwtdxy,deeprechxy,rechxy,errmsg,errflg) & !$OMP private(vegtyp,masslai,masssai,snd,dzsno,dzsnso,isnow) & diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.fv3.F90 b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.fv3.F90 index 6a556a88d..2cb892825 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.fv3.F90 +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/GFS_time_vary_pre.fv3.F90 @@ -73,7 +73,7 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr, kdt, julian, yearlen, ipt, lprnt, lssav, lsswr, lslwr, solhr, errmsg, errflg) use machine, only: kind_phys, kind_dbl_prec, kind_sngl_prec - + use w3emc, only: w3difdat implicit none integer, intent(in) :: idate(:) diff --git a/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F b/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F index bfd677800..bdf58c33f 100644 --- a/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F +++ b/physics/Interstitials/UFS_SCM_NEPTUNE/sfcsub.F @@ -3,7 +3,7 @@ !>\defgroup mod_sfcsub GFS sfcsub Module -!!\ingroup mod_GFS_phys_time_vary +!!\ingroup mod_GFS_phys_time_vary !> @{ !! This module contains grib code for each parameter-used in subroutines sfccycle() !! and setrmsk(). @@ -16,7 +16,7 @@ module sfccyc_module ! integer kpdtsf,kpdwet,kpdsno,kpdzor,kpdais,kpdtg3,kpdplr,kpdgla, & kpdmxi,kpdscv,kpdsmc,kpdoro,kpdmsk,kpdstc,kpdacn,kpdveg, - & kpdvet,kpdsot,kpdsoc, + & kpdvet,kpdsot,kpdsoc, & kpdvmn,kpdvmx,kpdslp,kpdabs &, kpdsnd, kpdabs_0, kpdabs_1, kpdalb(4) parameter(kpdtsf=11, kpdwet=86, kpdsno=65, kpdzor=83, @@ -54,6 +54,325 @@ end function message end module sfccyc_module +! this file only uses the procedure splat from the sp library. +! the local_sp module has splat and the two functions from lapack_gen.F it calls +! Note: the second time splat calls ludcmp, the original call had an extra +! argument has now been removed + module local_sp + contains + +! lubksb modified from https://github.com/NOAA-EMC/NCEPLIBS-sp/blob/develop/src/lapack_gen.F +! @brief Two Numerical Recipes routines for matrix inversion From Numerical Recipes. +! +! ### Program History Log +! Date | Programmer | Comments +! -----|------------|--------- +! 2012-11-05 | E.Mirvis | separated this generic LU from the splat.F + +! Solves a system of linear equations, follows call to ludcmp(). +! +! @param A +! @param N +! @param NP +! @param INDX +! @param B + subroutine lubksb(A,N,NP,INDX,B) + real A(NP,NP),B(N) + integer INDX(N) + II=0 + DO 12 I=1,N + LL=INDX(I) + SUM=B(LL) + B(LL)=B(I) + IF (II.NE.0)THEN + DO 11 J=II,I-1 + SUM=SUM-A(I,J)*B(J) + 11 CONTINUE + ELSE IF (SUM.NE.0.) THEN + II=I + ENDIF + B(I)=SUM + 12 CONTINUE + DO 14 I=N,1,-1 + SUM=B(I) + IF(I.LT.N)THEN + DO 13 J=I+1,N + SUM=SUM-A(I,J)*B(J) + 13 CONTINUE + ENDIF + B(I)=SUM/A(I,I) + 14 CONTINUE + RETURN + end subroutine lubksb + +! ludcmp modified from https://github.com/NOAA-EMC/NCEPLIBS-sp/blob/develop/src/lapack_gen.F +! Replaces an NxN matrix a with the LU decomposition. +! +! @param A +! @param N +! @param NP +! @param INDX + subroutine ludcmp(A,N,NP,INDX) +! parameter (NMAX=400,TINY=1.0E-20) + parameter (TINY=1.0E-20) +!==EM==^^^ + + real A(NP,NP),VV(N),D +! real A(NP,NP),VV(NMAX),D +!==EM==^^^ + integer INDX(N) + D=1. + DO 12 I=1,N + AAMAX=0. + DO 11 J=1,N + IF (ABS(A(I,J)).GT.AAMAX) AAMAX=ABS(A(I,J)) + 11 CONTINUE + IF (AAMAX.EQ.0.) print *, 'SINGULAR MATRIX.' + VV(I)=1./AAMAX + 12 CONTINUE + DO 19 J=1,N + IF (J.GT.1) THEN + DO 14 I=1,J-1 + SUM=A(I,J) + IF (I.GT.1)THEN + DO 13 K=1,I-1 + SUM=SUM-A(I,K)*A(K,J) + 13 CONTINUE + A(I,J)=SUM + ENDIF + 14 CONTINUE + ENDIF + AAMAX=0. + DO 16 I=J,N + SUM=A(I,J) + IF (J.GT.1)THEN + DO 15 K=1,J-1 + SUM=SUM-A(I,K)*A(K,J) + 15 CONTINUE + A(I,J)=SUM + ENDIF + DUM=VV(I)*ABS(SUM) + IF (DUM.GE.AAMAX) THEN + IMAX=I + AAMAX=DUM + ENDIF + 16 CONTINUE + IF (J.NE.IMAX)THEN + DO 17 K=1,N + DUM=A(IMAX,K) + A(IMAX,K)=A(J,K) + A(J,K)=DUM + 17 CONTINUE + D=-D + VV(IMAX)=VV(J) + ENDIF + INDX(J)=IMAX + IF(J.NE.N)THEN + IF(A(J,J).EQ.0.)A(J,J)=TINY + DUM=1./A(J,J) + DO 18 I=J+1,N + A(I,J)=A(I,J)*DUM + 18 CONTINUE + ENDIF + 19 CONTINUE + IF(A(N,N).EQ.0.)A(N,N)=TINY + RETURN + end subroutine ludcmp + +! splat modified from https://github.com/NOAA-EMC/NCEPLIBS-sp/blob/develop/src/splat.F +! @brief Computes cosines of colatitude and Gaussian weights +! for sets of latitudes. +! +! ### Program History Log +! Date | Programmer | Comments +! -----|------------|--------- +! 96-02-20 | Iredell | Initial. +! 97-10-20 | Iredell | Adjust precision. +! 98-06-11 | Iredell | Generalize precision using FORTRAN 90 intrinsic. +! 1998-12-03 | Iredell | Generalize precision further. +! 1998-12-03 | Iredell | Uses AIX ESSL BLAS calls. +! 2009-12-27 | D. Stark | Updated to switch between ESSL calls on an AIX platform, and Numerical Recipies calls elsewise. +! 2010-12-30 | Slovacek | Update alignment so preprocessor does not cause compilation failure. +! 2012-09-01 | E. Mirvis & M.Iredell | Merging & debugging linux errors of _d and _8 using generic LU factorization. +! 2012-11-05 | E. Mirvis | Generic FFTPACK and LU lapack were removed. +! +! @author Iredell @date 96-02-20 + +! Computes cosines of colatitude and Gaussian weights +! for one of the following specific global sets of latitudes. +! - Gaussian latitudes (IDRT=4) +! - Equally-spaced latitudes including poles (IDRT=0) +! - Equally-spaced latitudes excluding poles (IDRT=256) +! +! The Gaussian latitudes are located at the zeroes of the +! Legendre polynomial of the given order. These latitudes +! are efficient for reversible transforms from spectral space. +! (About twice as many equally-spaced latitudes are needed.) +! The weights for the equally-spaced latitudes are based on +! Ellsaesser (JAM,1966). (No weight is given the pole point.) +! Note that when analyzing grid to spectral in latitude pairs, +! if an equator point exists, its weight should be halved. +! This version invokes the ibm essl matrix solver. +! +! @param[in] IDRT grid identifier +! - 4 for Gaussian grid +! - 0 for equally-spaced grid including poles +! - 256 for equally-spaced grid excluding poles +! @param[in] JMAX number of latitudes +! @param[out] SLAT sines of latitude +! @param[out] WLAT Gaussian weights +! +! @author Iredell @date 96-02-20 + subroutine splat(idrt,jmax,slat,wlat) + real slat(jmax),wlat(jmax) + integer,parameter:: kd=selected_real_kind(15,45) + real(kind=kd):: pk(jmax/2),pkm1(jmax/2),pkm2(jmax/2) + real(kind=kd):: slatd(jmax/2),sp,spmax,eps=10.*epsilon(sp) + parameter(jz=50) + real bz(jz) + data bz / 2.4048255577, 5.5200781103, + $ 8.6537279129, 11.7915344391, 14.9309177086, 18.0710639679, + $ 21.2116366299, 24.3524715308, 27.4934791320, 30.6346064684, + $ 33.7758202136, 36.9170983537, 40.0584257646, 43.1997917132, + $ 46.3411883717, 49.4826098974, 52.6240518411, 55.7655107550, + $ 58.9069839261, 62.0484691902, 65.1899648002, 68.3314693299, + $ 71.4729816036, 74.6145006437, 77.7560256304, 80.8975558711, + $ 84.0390907769, 87.1806298436, 90.3221726372, 93.4637187819, + $ 96.6052679510, 99.7468198587, 102.888374254, 106.029930916, + $ 109.171489649, 112.313050280, 115.454612653, 118.596176630, + $ 121.737742088, 124.879308913, 128.020877005, 131.162446275, + $ 134.304016638, 137.445588020, 140.587160352, 143.728733573, + $ 146.870307625, 150.011882457, 153.153458019, 156.295034268 / + real:: dlt,d1=1. + real awork((jmax+1)/2,((jmax+1)/2)),bwork(((jmax+1)/2)) + integer:: jhe,jho,j0=0 + integer ipvt((jmax+1)/2) + real, parameter :: pi=3.14159265358979, c=(1.-(2./pi)**2)*0.25 + + ! GAUSSIAN LATITUDES + if(idrt.eq.4) then + jh=jmax/2 + jhe=(jmax+1)/2 + r=1./sqrt((jmax+0.5)**2+c) + do j=1,min(jh,jz) + slatd(j)=cos(bz(j)*r) + enddo + do j=jz+1,jh + slatd(j)=cos((bz(jz)+(j-jz)*pi)*r) + enddo + spmax=1. + do while(spmax.gt.eps) + spmax=0. + do j=1,jh + pkm1(j)=1. + pk(j)=slatd(j) + enddo + do n=2,jmax + do j=1,jh + pkm2(j)=pkm1(j) + pkm1(j)=pk(j) + pk(j)=((2*n-1)*slatd(j)*pkm1(j)-(n-1)*pkm2(j))/n + enddo + enddo + do j=1,jh + sp=pk(j)*(1.-slatd(j)**2)/(jmax*(pkm1(j)-slatd(j)*pk(j))) + slatd(j)=slatd(j)-sp + spmax=max(spmax,abs(sp)) + enddo + enddo +CDIR$ IVDEP + do j=1,jh + slat(j)=slatd(j) + wlat(j)=(2.*(1.-slatd(j)**2))/(jmax*pkm1(j))**2 + slat(jmax+1-j)=-slat(j) + wlat(jmax+1-j)=wlat(j) + enddo + if(jhe.gt.jh) then + slat(jhe)=0. + wlat(jhe)=2./jmax**2 + do n=2,jmax,2 + wlat(jhe)=wlat(jhe)*n**2/(n-1)**2 + enddo + endif + + ! equally-spaced latitudes including poles + elseif(idrt.eq.0) then + jh=jmax/2 + jhe=(jmax+1)/2 + jho=jhe-1 + dlt=pi/(jmax-1) + slat(1)=1. + do j=2,jh + slat(j)=cos((j-1)*dlt) + enddo + do js=1,jho + do j=1,jho + awork(js,j)=cos(2*(js-1)*j*dlt) + enddo + enddo + do js=1,jho + bwork(js)=-d1/(4*(js-1)**2-1) + enddo + + call ludcmp(awork,jho,jhe,ipvt) + call lubksb(awork,jho,jhe,ipvt,bwork) + + wlat(1)=0. + do j=1,jho + wlat(j+1)=bwork(j) + enddo +CDIR$ IVDEP + do j=1,jh + print *, j, jmax, jmax+1-j + slat(jmax+1-j)=-slat(j) + wlat(jmax+1-j)=wlat(j) + enddo + if(jhe.gt.jh) then + slat(jhe)=0. + wlat(jhe)=2.*wlat(jhe) + endif + + ! EQUALLY-SPACED LATITUDES EXCLUDING POLES + elseif(idrt.eq.256) then + jh=jmax/2 + jhe=(jmax+1)/2 + jho=jhe + dlt=pi/jmax + slat(1)=1. + do j=1,jh + slat(j)=cos((j-0.5)*dlt) + enddo + do js=1,jho + do j=1,jho + awork(js,j)=cos(2*(js-1)*(j-0.5)*dlt) + enddo + enddo + do js=1,jho + bwork(js)=-d1/(4*(js-1)**2-1) + enddo + + call ludcmp(awork,jho,jhe,ipvt) + call lubksb(awork,jho,jhe,ipvt,bwork) + + wlat(1)=0. + do j=1,jho + wlat(j)=bwork(j) + enddo +CDIR$ IVDEP + do j=1,jh + slat(jmax+1-j)=-slat(j) + wlat(jmax+1-j)=wlat(j) + enddo + if(jhe.gt.jh) then + slat(jhe)=0. + wlat(jhe)=2.*wlat(jhe) + endif + endif + + return + end subroutine splat + end module local_sp + !>\ingroup mod_GFS_phys_time_vary !! This subroutine reads or interpolates surface climatology data in analysis !! and forecast mode. @@ -82,7 +401,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & &, vmnfcs,vmxfcs,slpfcs,absfcs & &, tsffcs,snofcs,zorfcs,albfcs,tg3fcs & &, cnpfcs,smcfcs,stcfcs,slifcs,aisfcs & - &, vegfcs,vetfcs,sotfcs,socfcs,alffcs & + &, vegfcs,vetfcs,sotfcs,socfcs,alffcs & &, cvfcs,cvbfcs,cvtfcs,me,nthrds,nlunit & &, sz_nml,input_nml_file & &, min_ice & @@ -126,17 +445,17 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & vetsmn,vetimx,vetimn,vetjmx,vetjmn, & & sotlmx,sotlmn,sotomx,sotomn,sotsmx, & & sotsmn,sotimx,sotimn,sotjmx,sotjmn, & - & soclmx,soclmn,socomx,socomn,socsmx, & - & socsmn,socimx,socimn,socjmx,socjmn, & + & soclmx,soclmn,socomx,socomn,socsmx, & + & socsmn,socimx,socimn,socjmx,socjmn, & & alslmx,alslmn,alsomx,alsomn,alssmx, & & alssmn,alsimx,alsimn,alsjmx,alsjmn, & & epstsf,epsalb,epssno,epswet,epszor, & & epsplr,epsoro,epssmc,epsscv,eptsfc, & & epstg3,epsais,epsacn,epsveg,epsvet, & - & epssot,epssoc,epsalf,qctsfs,qcsnos,qctsfi, & + & epssot,epssoc,epsalf,qctsfs,qcsnos,qctsfi, & & aislim,snwmin,snwmax,cplrl,cplrs, & & cvegl,czors,csnol,csnos,czorl,csots, & - & csotl,csocs,csocl,cvwgs,cvetl,cvets,calfs, & + & csotl,csocs,csocl,cvwgs,cvetl,cvets,calfs, & & fcalfl,fcalfs,ccvt,ccnp,ccv,ccvb, & & calbl,calfl,calbs,ctsfs,grboro, & & grbmsk,ctsfl,deltf,caisl,caiss, & @@ -145,7 +464,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & faiss,fsnol,bltmsk,falbs,cvegs,percrit, & & deltsfc,critp2,critp3,blnmsk,critp1, & & fcplrl,fcplrs,fczors,fvets,fsotl,fsots, & - & fsocl,fsocs, & + & fsocl,fsocs, & & fvetl,fplrs,fvegl,fvegs,fcsnol,fcsnos, & & fczorl,fcalbs,fctsfl,fctsfs,fcalbl, & & falfs,falfl,fh,crit,zsca,ztsfc,tem1,tem2 & @@ -171,7 +490,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & &, sihnew integer imsk,jmsk,ifp,irtscv,irtacn,irtais,irtsno,irtzor, & - & irtalb,irtsot,irtsoc,irtalf,j,irtvet,irtsmc,irtstc,irtveg,& + & irtalb,irtsot,irtsoc,irtalf,j,irtvet,irtsmc,irtstc,irtveg,& & irtwet,k,iprnt,kk,irttsf,iret,i,igrdbg,iy,im,id, & & icalbl,icalbs,icalfl,ictsfs,lugb,len,lsoil,ih, & & ictsfl,iczors,icplrl,icplrs,iczorl,icalfs,icsnol, & @@ -477,7 +796,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & ! character*500 fntsfc,fnwetc,fnsnoc,fnzorc,fnalbc,fnaisc & &, fnplrc,fntg3c,fnscvc,fnsmcc,fnstcc,fnacnc & - &, fnvegc,fnvetc,fnsotc,fnsocc & + &, fnvegc,fnvetc,fnsotc,fnsocc & &, fnvmnc,fnvmxc,fnslpc,fnabsc, fnalbc2 real (kind=kind_io8) tsfclm(len), wetclm(len), snoclm(len) & &, zorclm(len), albclm(len,4), aisclm(len) & @@ -501,7 +820,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & &, tg3anl(len), acnanl(len), cnpanl(len) & &, cvanl (len), cvbanl(len), cvtanl(len) & &, scvanl(len), tsfan2(len), veganl(len) & - &, vetanl(len), sotanl(len), socanl(len) & + &, vetanl(len), sotanl(len), socanl(len) & &, alfanl(len,2), slianl(len) & &, smcanl(len,lsoil), stcanl(len,lsoil) & &, sihanl(len), sicanl(len) & @@ -621,11 +940,11 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & namelist/namsfc/fnglac,fnmxic, & fntsfc,fnwetc,fnsnoc,fnzorc,fnalbc,fnaisc, & fnplrc,fntg3c,fnscvc,fnsmcc,fnstcc,fnacnc, - & fnvegc,fnvetc,fnsotc,fnsocc,fnalbc2, + & fnvegc,fnvetc,fnsotc,fnsocc,fnalbc2, & fnvmnc,fnvmxc,fnslpc,fnabsc, & fntsfa,fnweta,fnsnoa,fnzora,fnalba,fnaisa, & fnplra,fntg3a,fnscva,fnsmca,fnstca,fnacna, - & fnvega,fnveta,fnsota,fnsoca, + & fnvega,fnveta,fnsota,fnsoca, & fnvmna,fnvmxa,fnslpa,fnabsa, & fnmskh, & ldebug,lgchek,lqcbgs,critp1,critp2,critp3, @@ -634,7 +953,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & ftsfl,ftsfs,falbl,falbs,faisl,faiss,fsnol,fsnos, & fzorl,fzors,fplrl,fplrs,fsmcl,fsmcs, & fstcl,fstcs,fvegl,fvegs,fvetl,fvets,fsotl,fsots, - & fsocl,fsocs, + & fsocl,fsocs, & fctsfl,fctsfs,fcalbl,fcalbs,fcsnol,fcsnos, & fczorl,fczors,fcplrl,fcplrs,fcsmcl,fcsmcs, & fcstcl,fcstcs,fsalfl,fsalfs,fcalfl,flalfs, @@ -660,7 +979,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & data fnalbc2/'global_albedo4.1x1.grb'/ data fntsfc/'global_sstclim.2x2.grb'/ data fnsotc/'global_soiltype.1x1.grb'/ - data fnsocc/''/ + data fnsocc/''/ data fnvegc/'global_vegfrac.1x1.grb'/ data fnvetc/'global_vegtype.1x1.grb'/ data fnglac/'global_glacier.2x2.grb'/ @@ -697,7 +1016,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & data fnvega/' '/ data fnveta/' '/ data fnsota/' '/ - data fnsoca/' '/ + data fnsoca/' '/ !clu [+4l] add fn()a for vmn, vmx, abs, slp data fnvmna/' '/ data fnvmxa/' '/ @@ -719,7 +1038,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & data fplrl/0.0/, fplrs/0.0/ data fvetl/0.0/, fvets/99999.0/ data fsotl/0.0/, fsots/99999.0/ - data fsocl/0.0/, fsocs/99999.0/ + data fsocl/0.0/, fsocs/99999.0/ data fvegl/0.0/, fvegs/99999.0/ !cwu [+4l] add f()l and f()s for sih, sic and aislim, sihlim data fsihl/99999.0/, fsihs/99999.0/ @@ -770,7 +1089,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & ftsfl,ftsfs,falbl,falbs,faisl,faiss,fsnol,fsnos, & fzorl,fzors,fplrl,fplrs,fsmcl,fsmcs,falfl,falfs, & fstcl,fstcs,fvegl,fvegs,fvetl,fvets,fsotl,fsots, - & fsocl,fsocs, + & fsocl,fsocs, & fctsfl,fctsfs,fcalbl,fcalbs,fcsnol,fcsnos, & fczorl,fczors,fcplrl,fcplrs,fcsmcl,fcsmcs, & fcstcl,fcstcs,fcalfl,fcalfs, @@ -789,7 +1108,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & ctsfl, ctsfs, calbl, calfl, calbs, calfs, csmcs, & csnol, csnos, czorl, czors, cplrl, cplrs, cstcl, & cstcs, cvegl, cvwgs, cvetl, cvets, csotl, csots, - & csocl, csocs, + & csocl, csocs, & csmcl !cwu [+1l] add c()l and c()s for sih, sic &, csihl, csihs, csicl, csics @@ -1137,16 +1456,16 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & fntsfc,fnwetc,fnsnoc,fnzorc,fnalbc,fnaisc, & fntg3c,fnscvc,fnsmcc,fnstcc,fnacnc,fnvegc, - & fnvetc,fnsotc,fnsocc, + & fnvetc,fnsotc,fnsocc, & fnvmnc,fnvmxc,fnslpc,fnabsc, & tsfclm,tsfcl2,wetclm,snoclm,zorclm,albclm,aisclm, & tg3clm,cvclm ,cvbclm,cvtclm, & cnpclm,smcclm,stcclm,sliclm,scvclm,acnclm,vegclm, - & vetclm,sotclm,socclm,alfclm, + & vetclm,sotclm,socclm,alfclm, & vmnclm,vmxclm,slpclm,absclm, & kpdtsf,kpdwet,kpdsno,kpdzor,kpdalb,kpdais, & kpdtg3,kpdscv,kpdacn,kpdsmc,kpdstc,kpdveg, - & kpdvet,kpdsot,kpdsoc,kpdalf,tsfcl0, + & kpdvet,kpdsot,kpdsoc,kpdalf,tsfcl0, & kpdvmn,kpdvmx,kpdslp,kpdabs, & deltsfc, lanom &, imsk, jmsk, slmskh, rla, rlo, gausm, blnmsk, bltmsk,me @@ -1418,7 +1737,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call monitr('vegclm',vegclm,sliclm,snoclm,len) call monitr('vetclm',vetclm,sliclm,snoclm,len) call monitr('sotclm',sotclm,sliclm,snoclm,len) - call monitr('socclm',socclm,sliclm,snoclm,len) + call monitr('socclm',socclm,sliclm,snoclm,len) !cwu [+2l] add sih, sic call monitr('sihclm',sihclm,sliclm,snoclm,len) call monitr('sicclm',sicclm,sliclm,snoclm,len) @@ -1442,13 +1761,13 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call filanl(tsfanl,tsfan2,wetanl,snoanl,zoranl,albanl,aisanl, & tg3anl,cvanl ,cvbanl,cvtanl, & cnpanl,smcanl,stcanl,slianl,scvanl,veganl, - & vetanl,sotanl,socanl,alfanl, + & vetanl,sotanl,socanl,alfanl, & sihanl,sicanl, & vmnanl,vmxanl,slpanl,absanl, & tsfclm,tsfcl2,wetclm,snoclm,zorclm,albclm,aisclm, & tg3clm,cvclm ,cvbclm,cvtclm, & cnpclm,smcclm,stcclm,sliclm,scvclm,vegclm, - & vetclm,sotclm,socclm,alfclm, + & vetclm,sotclm,socclm,alfclm, & sihclm,sicclm, & vmnclm,vmxclm,slpclm,absclm, & len,lsoil) @@ -1478,20 +1797,20 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call analy(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & fntsfa,fnweta,fnsnoa,fnzora,fnalba,fnaisa, & fntg3a,fnscva,fnsmca,fnstca,fnacna,fnvega, - & fnveta,fnsota,fnsoca, - & fnvmna,fnvmxa,fnslpa,fnabsa, + & fnveta,fnsota,fnsoca, + & fnvmna,fnvmxa,fnslpa,fnabsa, & tsfanl,wetanl,snoanl,zoranl,albanl,aisanl, & tg3anl,cvanl ,cvbanl,cvtanl, & smcanl,stcanl,slianl,scvanl,acnanl,veganl, - & vetanl,sotanl,socanl,alfanl,tsfan0, + & vetanl,sotanl,socanl,alfanl,tsfan0, & vmnanl,vmxanl,slpanl,absanl, & kpdtsf,kpdwet,kpdsno,kpdsnd,kpdzor,kpdalb,kpdais, & kpdtg3,kpdscv,kpdacn,kpdsmc,kpdstc,kpdveg, - & kpdvet,kpdsot,kpdsoc,kpdalf, + & kpdvet,kpdsot,kpdsoc,kpdalf, & kpdvmn,kpdvmx,kpdslp,kpdabs, & irttsf,irtwet,irtsno,irtzor,irtalb,irtais, & irttg3,irtscv,irtacn,irtsmc,irtstc,irtveg, - & irtvet,irtsot,irtsoc,irtalf + & irtvet,irtsot,irtsoc,irtalf &, irtvmn,irtvmx,irtslp,irtabs, & imsk, jmsk, slmskh, rla, rlo, gausm, blnmsk, bltmsk &, me, lanom) @@ -1803,7 +2122,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call monitr('veganl',veganl,slianl,snoanl,len) call monitr('vetanl',vetanl,slianl,snoanl,len) call monitr('sotanl',sotanl,slianl,snoanl,len) - call monitr('socanl',socanl,slianl,snoanl,len) + call monitr('socanl',socanl,slianl,snoanl,len) !cwu [+2l] add sih, sic call monitr('sihanl',sihanl,slianl,snoanl,len) call monitr('sicanl',sicanl,slianl,snoanl,len) @@ -1835,7 +2154,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call filfcs(tsffcs,wetfcs,snofcs,zorfcs,albfcs, & tg3fcs,cvfcs ,cvbfcs,cvtfcs, & cnpfcs,smcfcs,stcfcs,slifcs,aisfcs, - & vegfcs,vetfcs,sotfcs,socfcs,alffcs, + & vegfcs,vetfcs,sotfcs,socfcs,alffcs, !cwu [+1l] add ()fcs for sih, sic & sihfcs,sicfcs, !clu [+1l] add ()fcs for vmn, vmx, slp, abs @@ -1843,7 +2162,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & tsfanl,wetanl,snoanl,zoranl,albanl, & tg3anl,cvanl ,cvbanl,cvtanl, & cnpanl,smcanl,stcanl,slianl,aisanl, - & veganl,vetanl,sotanl,socanl,alfanl, + & veganl,vetanl,sotanl,socanl,alfanl, !cwu [+1l] add ()anl for sih, sic & sihanl,sicanl, !clu [+1l] add ()anl for vmn, vmx, slp, abs @@ -2000,7 +2319,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & sotlmx,sotlmn,sotomx,sotomn,sotimx,sotimn, & sotjmx,sotjmn,sotsmx,sotsmn,epssot, & rla,rlo,len,kqcm,percrit,lgchek,me) - call qcmxmn('socf ',socfcs,slmskl,snofcs,icefl1, + call qcmxmn('socf ',socfcs,slmskl,snofcs,icefl1, & soclmx,soclmn,socomx,socomn,socimx,socimn, & socjmx,socjmn,socsmx,socsmn,epssoc, & rla,rlo,len,kqcm,percrit,lgchek,me) @@ -2123,19 +2442,19 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & tsffcs,wetfcs,snofcs,zorfcs,albfcs,aisfcs, & cvfcs ,cvbfcs,cvtfcs, & cnpfcs,smcfcs,stcfcs,slifcs,vegfcs, - & vetfcs,sotfcs,socfcs,alffcs, + & vetfcs,sotfcs,socfcs,alffcs, & sihanl,sicanl, & vmnanl,vmxanl,slpanl,absanl, & tsfanl,tsfan2,wetanl,snoanl,zoranl,albanl,aisanl, & cvanl ,cvbanl,cvtanl, & cnpanl,smcanl,stcanl,slianl,veganl, - & vetanl,sotanl,socanl,alfanl, + & vetanl,sotanl,socanl,alfanl, & ctsfl,calbl,caisl,csnol,csmcl,czorl,cstcl,cvegl, & ctsfs,calbs,caiss,csnos,csmcs,czors,cstcs,cvegs, - & ccv,ccvb,ccvt,ccnp,cvetl,cvets,csotl,csots,csocl,csocs, + & ccv,ccvb,ccvt,ccnp,cvetl,cvets,csotl,csots,csocl,csocs, & calfl,calfs, & csihl,csihs,csicl,csics, - & cvmnl,cvmns,cvmxl,cvmxs,cslpl,cslps,cabsl,cabss, + & cvmnl,cvmns,cvmxl,cvmxs,cslpl,cslps,cabsl,cabss, & irttsf,irtwet,irtsno,irtzor,irtalb,irtais, & irttg3,irtscv,irtacn,irtsmc,irtstc,irtveg, & irtvmn,irtvmx,irtslp,irtabs, @@ -2164,7 +2483,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & ! if(lprnt) print *,'tsfanl=',tsfanl(iprnt),' tsffcs=',tsffcs(iprnt) ! if(lprnt) print *,' slian=',slianl(iprnt),' slifn=',slifcs(iprnt) ! if(lprnt) print *,' stcan=',stcanl(iprnt,:) - + ! set tsfc to tsnow over snow ! call snosfc(snoanl,tsfanl,tsfsmx,len,me) @@ -2235,7 +2554,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & sotjmx,sotjmn,sotsmx,sotsmn,epssot, & rla,rlo,len,kqcm,percrit,lgchek,me) - call qcmxmn('socm ',socanl,slmskl,snoanl,icefl1, + call qcmxmn('socm ',socanl,slmskl,snoanl,icefl1, & soclmx,soclmn,socomx,socomn,socimx,socimn, & socjmx,socjmn,socsmx,socsmn,epssoc, & rla,rlo,len,kqcm,percrit,lgchek,me) @@ -2348,7 +2667,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & & zorfcsd(len), slifcsd(len), aisfcsd(len), & & cnpfcsd(len), vegfcsd(len), vetfcsd(len), & & sotfcsd(len), socfcsd(len),sihfcsd(len), & - & sicfcsd(len), & + & sicfcsd(len), & & vmnfcsd(len), vmxfcsd(len), slpfcsd(len), & & absfcsd(len)) allocate (smcfcsd(len,lsoil), stcfcsd(len,lsoil), & @@ -2366,7 +2685,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & vegfcsd(i) = veganl(i) - vegfcs(i) vetfcsd(i) = vetanl(i) - vetfcs(i) sotfcsd(i) = sotanl(i) - sotfcs(i) - socfcsd(i) = socanl(i) - socfcs(i) + socfcsd(i) = socanl(i) - socfcs(i) !clu [+2l] add sih, sic sihfcsd(i) = sihanl(i) - sihfcs(i) sicfcsd(i) = sicanl(i) - sicfcs(i) @@ -2420,7 +2739,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call monitr('vegdif',vegfcsd,slianl,snoanl,len) call monitr('vetdif',vetfcsd,slianl,snoanl,len) call monitr('sotdif',sotfcsd,slianl,snoanl,len) - call monitr('socdif',socfcsd,slianl,snoanl,len) + call monitr('socdif',socfcsd,slianl,snoanl,len) !cwu [+2l] add sih, sic call monitr('sihdif',sihfcsd,slianl,snoanl,len) call monitr('sicdif',sicfcsd,slianl,snoanl,len) @@ -2431,7 +2750,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & call monitr('absdif',absfcsd,slianl,snoanl,len) endif deallocate (tsffcsd, snofcsd, tg3fcsd, zorfcsd, slifcsd, & - & aisfcsd, cnpfcsd, vegfcsd, vetfcsd, sotfcsd,socfcsd, & + & aisfcsd, cnpfcsd, vegfcsd, vetfcsd, sotfcsd,socfcsd, & & sihfcsd, sicfcsd, vmnfcsd, vmxfcsd, slpfcsd, & & absfcsd) deallocate (smcfcsd, stcfcsd, albfcsd) @@ -2454,7 +2773,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & vegfcs(i) = veganl(i) vetfcs(i) = vetanl(i) sotfcs(i) = sotanl(i) - socfcs(i) = socanl(i) + socfcs(i) = socanl(i) !clu [+4l] add vmn, vmx, slp, abs vmnfcs(i) = vmnanl(i) vmxfcs(i) = vmxanl(i) @@ -2559,7 +2878,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & ! set liquid soil moisture to a flag value of 1.0 if (landice) then do i = 1, len - if (slifcs(i) == 1.0 .and. + if (slifcs(i) == 1.0 .and. & nint(vetfcs(i)) == veg_type_landice) then do k=1, lsoil slcfcs(i,k) = 1.0 @@ -2570,7 +2889,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & ! ! ensure the consistency between snwdph and sheleg ! - if(fsnol < 99999.) then + if(fsnol < 99999.) then if(me == 0) then print *,'dbgx -- scale snwdph from sheleg' endif @@ -2622,7 +2941,7 @@ subroutine sfccycle(lugb,len,lsoil,sig1t,deltsfc & ! if(lprnt) print *,' stcfcsend=',stcfcs(iprnt,:) ! if(lprnt) print *,' slifcsend=',slifcs(iprnt) return - end subroutine sfccycle + end subroutine sfccycle !>\ingroup mod_sfcsub !! This subroutine counts number of points for the four surface @@ -3194,7 +3513,7 @@ subroutine la2ga(regin,imxin,jmxin,rinlon,rinlat,rlon,rlat,inttyp,& if (inttyp /=1) allocate (ifill(num_threads)) ! -!$omp parallel do default(none) +!$omp parallel do default(none) !$omp+private(i1_t,i2_t,len_thread,it,i,ii,i1,i2) !$omp+private(j,j1,j2,jq,ix,jy,nx,kxs,kxt,kmami) !$omp+private(alamd,denom,rnume,aphi,x,y,wsum,wsumiv,sum1,sum2) @@ -3455,7 +3774,7 @@ subroutine la2ga(regin,imxin,jmxin,rinlon,rinlat,rlon,rlat,inttyp,& if (num_threads == 1) then print*,'no matching mask found ',i,i1,j1,ix,jx & - &, ' slmask=',slmask(i),' me=',me & + &, ' slmask=',slmask(i),' me=',me & &, ' outlon=',outlon(i),' outlat=',outlat(i) &, 'set to default value.' endif @@ -3696,7 +4015,7 @@ subroutine filanl(tsfanl,tsfan2,wetanl,snoanl,zoranl,albanl, & & cnpanl,smcanl,stcanl,slianl,scvanl,veganl, & & vetanl,sotanl,socanl,alfanl, & !socanl: soil color & sihanl,sicanl, & !cwu [+1l] add ()anl for sih, sic - & vmnanl,vmxanl,slpanl,absanl, & !clu [+1l] add ()anl for vmn, vmx, slp, abs + & vmnanl,vmxanl,slpanl,absanl, & !clu [+1l] add ()anl for vmn, vmx, slp, abs & tsfclm,tsfcl2,wetclm,snoclm,zorclm,albclm, & & aisclm, & & tg3clm,cvclm ,cvbclm,cvtclm, & @@ -4779,14 +5098,14 @@ subroutine merge(len,lsoil,iy,im,id,ih,fh,deltsfc, & & ctsfl,calbl,caisl,csnol,csmcl,czorl,cstcl,cvegl, & & ctsfs,calbs,caiss,csnos,csmcs,czors,cstcs,cvegs, & & ccv,ccvb,ccvt,ccnp,cvetl,cvets,csotl,csots, & - & csocl,csocs, & !csocl,csocs:soil color + & csocl,csocs, & !csocl,csocs:soil color & calfl,calfs, & & csihl,csihs,csicl,csics, & & cvmnl,cvmns,cvmxl,cvmxs,cslpl,cslps,cabsl,cabss, & & irttsf,irtwet,irtsno,irtzor,irtalb,irtais, & & irttg3,irtscv,irtacn,irtsmc,irtstc,irtveg, & & irtvmn,irtvmx,irtslp,irtabs, & - & irtvet,irtsot,irtsoc,irtalf, landice, me) + & irtvet,irtsot,irtsoc,irtalf, landice, me) use machine , only : kind_io8,kind_io4 use sfccyc_module, only : veg_type_landice, soil_type_landice, & & num_threads, zero, one,soil_color_landice @@ -4805,7 +5124,7 @@ subroutine merge(len,lsoil,iy,im,id,ih,fh,deltsfc, & & qvetl,rtsfl,calbs,caiss,ctsfs,czorl,cvegl, & & csnos,ccvb,ccvt,ccv,czors,cvegs,caisl,csnol, & & calbl,fh,ctsfl,ccnp,csots,calfl,csotl,cvetl, & - & csocl,csocs, & !csocl,csocs:soil color + & csocl,csocs, & !csocl,csocs:soil color & cvets,calfs,deltsfc, & & csihl,csihs,csicl,csics, & & rsihl,rsihs,rsicl,rsics, & @@ -5044,7 +5363,7 @@ subroutine merge(len,lsoil,iy,im,id,ih,fh,deltsfc, & qvegs = 1. - rvegs qvets = 1. - rvets qsots = 1. - rsots - qsocs = 1. - rsocs + qsocs = 1. - rsocs qsihs = 1. - rsihs qsics = 1. - rsics @@ -5111,7 +5430,7 @@ subroutine merge(len,lsoil,iy,im,id,ih,fh,deltsfc, & ! albanl(i) = albfcs(i)*ralbs + albanl(i)*qalbs aisanl(i) = aisfcs(i)*raiss + aisanl(i)*qaiss snoanl(i) = snofcs(i)*rsnos + snoanl(i)*qsnos - + zoranl(i) = zorfcs(i)*rzors + zoranl(i)*qzors veganl(i) = vegfcs(i)*rvegs + veganl(i)*qvegs sihanl(i) = sihfcs(i)*rsihs + sihanl(i)*qsihs @@ -5161,7 +5480,7 @@ subroutine merge(len,lsoil,iy,im,id,ih,fh,deltsfc, & if (nint(slianl(i)) == 1) then if (nint(vetanl(i)) == veg_type_landice) then sotanl(i) = soil_type_landice - socanl(i) = soil_color_landice + socanl(i) = soil_color_landice veganl(i) = 0.0 slpanl(i) = 9.0 vmnanl(i) = 0.0 @@ -5225,7 +5544,7 @@ end subroutine merge !>\ingroup mod_sfcsub subroutine newice(slianl,slifcs,tsfanl,tsffcs,len,lsoil, & - & sihnew,sicnew,sihanl,sicanl, & !cwu [+1l] add sihnew,sicnew,sihanl,sicanl + & sihnew,sicnew,sihanl,sicanl, & !cwu [+1l] add sihnew,sicnew,sihanl,sicanl & albanl,snoanl,zoranl,smcanl,stcanl, & & albsea,snosea,zorsea,smcsea,smcice, & & tsfmin,tsfice,albice,zorice,tgice, & @@ -5237,7 +5556,7 @@ subroutine newice(slianl,slifcs,tsfanl,tsffcs,len,lsoil, & real (kind=kind_io8) tgice,albice,zorice,tsfice,albsea,snosea, & & smcice,tsfmin,zorsea,smcsea !cwu [+1l] add sicnew,sihnew - &, sicnew,sihnew + &, sicnew,sihnew integer i,me,kount1,kount2,k,len,lsoil real (kind=kind_io8) slianl(len), slifcs(len), & tsffcs(len),tsfanl(len) @@ -5276,7 +5595,7 @@ subroutine newice(slianl,slifcs,tsfanl,tsffcs,len,lsoil, & do k = 1, lsoil smcanl(i,k) = smcsea !cwu [+1l] set stcanl to tgice (over sea-ice) - stcanl(i,k) = tgice + stcanl(i,k) = tgice enddo !cwu [+2l] set siganl and sicanl sihanl(i) = 0. @@ -6255,6 +6574,7 @@ subroutine setrmsk(kpds5,slmask,igaul,jgaul,wlon,rnlat, & &, gaus,blno, blto, kgds1, kpds4, lbms) use machine , only : kind_io8,kind_io4,kind_dbl_prec use sfccyc_module + use local_sp, only : splat implicit none real (kind=kind_io8) blno,blto,wlon,rnlat,crit,data_max integer i,j,ijmax,jgaul,igaul,kpds5,jmax,imax, kgds1, kspla @@ -6741,6 +7061,7 @@ subroutine ga2la(gauin,imxin,jmxin,regout,imxout,jmxout, & & wlon,rnlat,rlnout,rltout,gaus,blno, blto) use machine , only : kind_io8,kind_io4,kind_dbl_prec use sfccyc_module , only : num_threads + use local_sp, only : splat implicit none integer i1,i2,j2,ishft,i,jj,j1,jtem,jmxout,imxin,jmxin,imxout, & & j,iret @@ -6988,7 +7309,7 @@ subroutine landtyp(vegtype,soiltype,colortype,slptype,slmask,len) implicit none integer i,len real (kind=kind_io8) vegtype(len),soiltype(len),slmask(len) & - &, slptype(len),colortype(len) + &, slptype(len),colortype(len) ! ! make sure that the soil type and veg type are non-zero over land ! @@ -7008,6 +7329,7 @@ end subroutine landtyp subroutine gaulat(gaul,k) ! use machine , only : kind_io8,kind_io4,kind_dbl_prec + use local_sp, only : splat implicit none integer n,k real (kind=kind_io8) radi @@ -7061,7 +7383,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & & vmnclm,vmxclm,slpclm,absclm, & & kpdtsf,kpdwet,kpdsno,kpdzor,kpdalb,kpdais, & & kpdtg3,kpdscv,kpdacn,kpdsmc,kpdstc,kpdveg, & - & kpdvet,kpdsot,kpdsoc,kpdalf,tsfcl0, & + & kpdvet,kpdsot,kpdsoc,kpdalf,tsfcl0, & & kpdvmn,kpdvmx,kpdslp,kpdabs, & & deltsfc, lanom & &, imsk, jmsk, slmskh, outlat, outlon & @@ -7080,7 +7402,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & & jy,mon1,is2,isx,kpd9,is1,l,nn,mon2,mon,is,kpdsno, & & kpdzor,kpdtsf,kpdwet,kpdscv,kpdacn,kpdais,kpdtg3,im,id, & & lugb,iy,len,lsoil,ih,kpdsmc,iprnt,me,m1,m2,k1,k2, & - & kpdvet,kpdsot,kpdsoc,kpdstc,kpdveg,jmsk,imsk,j,ialb & + & kpdvet,kpdsot,kpdsoc,kpdstc,kpdveg,jmsk,imsk,j,ialb & &, kpdvmn,kpdvmx,kpdslp,kpdabs,landice_cat integer kpdalb(4), kpdalf(2) ! @@ -7096,7 +7418,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & & cnpclm(len), & & smcclm(len,lsoil),stcclm(len,lsoil), & & sliclm(len),scvclm(len),vegclm(len), & - & vetclm(len),sotclm(len),socclm(len),alfclm(len,2) & + & vetclm(len),sotclm(len),socclm(len),alfclm(len,2) & &, vmnclm(len),vmxclm(len),slpclm(len),absclm(len) real (kind=kind_io8) slmskh(imsk,jmsk) real (kind=kind_io8) outlat(len), outlon(len) @@ -7136,7 +7458,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & & zor(:,:),wet(:,:), & ais(:,:), acn(:,:), scv(:,:), smc(:,:,:), & tg3(:), alb(:,:,:), alf(:,:), - & vet(:), sot(:), soc(:), tsf2(:), + & vet(:), sot(:), soc(:), tsf2(:), & veg(:,:), stc(:,:,:) &, vmn(:), vmx(:), slp(:), absm(:) ! @@ -7145,7 +7467,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & data mon1s/0/, mon2s/0/, sea1s/0/, sea2s/0/ ! save first, tsf, sno, zor, wet, ais, acn, scv, smc, tg3, - & alb, alf, vet, sot, soc,tsf2, veg, stc, + & alb, alf, vet, sot, soc,tsf2, veg, stc, & vmn, vmx, slp, absm, & mon1s, mon2s, sea1s, sea2s, dayhf, k1, k2, m1, m2, & landice_cat @@ -7200,7 +7522,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & & wet(len,2), ais(len,2), acn(len,2), & scv(len,2), smc(len,lsoil,2), & tg3(len), alb(len,4,2), alf(len,2), - & vet(len), sot(len), soc(len),tsf2(len), + & vet(len), sot(len), soc(len),tsf2(len), !clu [+1l] add vmn, vmx, slp, abs & vmn(len), vmx(len), slp(len), absm(len), & veg(len,2), stc(len,lsoil,2)) @@ -7450,7 +7772,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & &, outlat, outlon, me) landice_cat=13 if (maxval(vet)> 13.0) landice_cat=15 - else + else call fixrdc_tile(fnvetc, tile_num_ch, i_index, j_index, & kpdvet, vet, 1, len, me) landice_cat=15 @@ -7604,7 +7926,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & mon = mon1 if (nn .eq. 2) mon = mon2 !cbosu -!cbosu new snowfree albedo database is monthly. +!cbosu new snowfree albedo database is monthly. if (ialb == 1 .or. ialb == 2) then if ( index(fnalbc, "tileX.nc") == 0) then ! grib file kpd7=-1 @@ -7866,7 +8188,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & ! two for strong zeneith angle dependent (visible and near ir) ! and two for weak zeneith angle dependent (vis ans nir) ! -!cbosu +!cbosu if (ialb == 0) then kpd7=-1 do k = 1, 4 @@ -8077,10 +8399,10 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & z0_season(1) = z0_igbp_min(ivtyp) z0_season(7) = z0_igbp_max(ivtyp) if (outlat(i) < 0.0) then - zorclm(i) = wei1y * z0_season(hyr2) + + zorclm(i) = wei1y * z0_season(hyr2) + & wei2y * z0_season(hyr1) else - zorclm(i) = wei1y * z0_season(hyr1) + + zorclm(i) = wei1y * z0_season(hyr1) + & wei2y * z0_season(hyr2) endif endif @@ -8217,7 +8539,7 @@ subroutine clima(lugb,iy,im,id,ih,fh,len,lsoil,slmskl,slmskw, & !clu ---------------------------------------------------------------------- ! !cbosu diagnostic print - if (me == 0) print*,'monthly albedo weights are ', + if (me == 0) print*,'monthly albedo weights are ', & wei1m,' for k', k1, wei2m, ' for k', k2 if (ialb == 1 .or. ialb == 2) then