From e198c5252ac222be5f5d1a3b1fea138f70fe4421 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Tue, 5 Nov 2019 15:51:23 +0000 Subject: [PATCH 01/18] this PR includes: 1)vlab #69814: post update,2)ufs issue #3, Add 3D reflectivity to restart file and restart reproducibility for regional fv3,3)ufs issue #5, Updates to WW3, 4)vlab #69735, update netcdf time units attribute when iau_offset --- .gitmodules | 4 +-- atmos_cubed_sphere | 2 +- gfsphysics/GFS_layer/GFS_physics_driver.F90 | 2 +- gfsphysics/GFS_layer/GFS_restart.F90 | 11 +++++- gfsphysics/GFS_layer/GFS_typedefs.F90 | 12 ++++++- io/FV3GFS_io.F90 | 9 ++--- io/module_write_netcdf.F90 | 39 ++++++++++++++++++--- io/module_wrt_grid_comp.F90 | 4 +-- io/post_gfs.F90 | 32 +++++++++++++---- 9 files changed, 90 insertions(+), 25 deletions(-) diff --git a/.gitmodules b/.gitmodules index 949d298df..04f5ebec0 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,7 +1,7 @@ [submodule "atmos_cubed_sphere"] path = atmos_cubed_sphere - url = https://github.com/NOAA-EMC/GFDL_atmos_cubed_sphere - branch = dev/emc + url = https://github.com/junwang-noaa/GFDL_atmos_cubed_sphere + branch = regfv3_rst [submodule "ccpp/framework"] path = ccpp/framework url = https://github.com/NCAR/ccpp-framework diff --git a/atmos_cubed_sphere b/atmos_cubed_sphere index 786447c83..cc28aa5aa 160000 --- a/atmos_cubed_sphere +++ b/atmos_cubed_sphere @@ -1 +1 @@ -Subproject commit 786447c8391a6806cd7b869bfa9dca69e3c95a48 +Subproject commit cc28aa5aa94817b24aea033cc19c1792e5927459 diff --git a/gfsphysics/GFS_layer/GFS_physics_driver.F90 b/gfsphysics/GFS_layer/GFS_physics_driver.F90 index 401fbbf86..5b67f7faa 100644 --- a/gfsphysics/GFS_layer/GFS_physics_driver.F90 +++ b/gfsphysics/GFS_layer/GFS_physics_driver.F90 @@ -1149,7 +1149,7 @@ subroutine GFS_physics_driver & if (fice(i) < one) then wet(i) = .true. ! Sfcprop%tsfco(i) = tgice - Sfcprop%tsfco(i) = max(Sfcprop%tisfc(i), tgice) + if (.not. Model%cplflx) Sfcprop%tsfco(i) = max(Sfcprop%tisfc(i), tgice) ! Sfcprop%tsfco(i) = max((Sfcprop%tsfc(i) - fice(i)*sfcprop%tisfc(i)) & ! / (one - fice(i)), tgice) endif diff --git a/gfsphysics/GFS_layer/GFS_restart.F90 b/gfsphysics/GFS_layer/GFS_restart.F90 index eafbcb9ba..a24cc0fc6 100644 --- a/gfsphysics/GFS_layer/GFS_restart.F90 +++ b/gfsphysics/GFS_layer/GFS_restart.F90 @@ -117,6 +117,9 @@ subroutine GFS_restart_populate (Restart, Model, Statein, Stateout, Sfcprop, & #endif Restart%num3d = Model%ntot3d + if(Model%lrefres) then + Restart%num3d = Model%ntot3d+1 + endif #ifdef CCPP ! GF if (Model%imfdeepcnv == 3) then @@ -252,7 +255,13 @@ subroutine GFS_restart_populate (Restart, Model, Statein, Stateout, Sfcprop, & Restart%data(nb,num)%var3p => Tbd(nb)%phy_f3d(:,:,num) enddo enddo - + if (Model%lrefres) then + num = Model%ntot3d+1 + restart%name3d(num) = 'ref_f3d' + do nb = 1,nblks + Restart%data(nb,num)%var3p => IntDiag(nb)%refl_10cm(:,:) + enddo + endif #ifdef CCPP !--- RAP/HRRR-specific variables, 3D num = Model%ntot3d diff --git a/gfsphysics/GFS_layer/GFS_typedefs.F90 b/gfsphysics/GFS_layer/GFS_typedefs.F90 index 7f8239a5a..39520b0d4 100644 --- a/gfsphysics/GFS_layer/GFS_typedefs.F90 +++ b/gfsphysics/GFS_layer/GFS_typedefs.F90 @@ -701,6 +701,9 @@ module GFS_typedefs !--- GFDL microphysical paramters logical :: lgfdlmprad !< flag for GFDL mp scheme and radiation consistency + !--- Thompson,GFDL mp parameter + logical :: lrefres !< flag for radar reflectivity in restart file + !--- land/surface model parameters integer :: lsm !< flag for land surface model lsm=1 for noah lsm integer :: lsm_noah=1 !< flag for NOAH land surface model @@ -2740,6 +2743,9 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & !--- GFDL microphysical parameters logical :: lgfdlmprad = .false. !< flag for GFDLMP radiation interaction + !--- Thompson,GFDL microphysical parameter + logical :: lrefres = .false. !< flag for radar reflectivity in restart file + !--- land/surface model parameters integer :: lsm = 1 !< flag for land surface model to use =0 for osu lsm; =1 for noah lsm; =2 for noah mp lsm; =3 for RUC lsm integer :: lsoil = 4 !< number of soil layers @@ -3023,7 +3029,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & mg_do_graupel, mg_do_hail, mg_nccons, mg_nicons, mg_ngcons, & mg_ncnst, mg_ninst, mg_ngnst, sed_supersat, do_sb_physics, & mg_alf, mg_qcmin, mg_do_ice_gmao, mg_do_liq_liu, & - ltaerosol, lradar, ttendlim, lgfdlmprad, & + ltaerosol, lradar, lrefres, ttendlim, lgfdlmprad, & !--- max hourly avg_max_length, & !--- land/surface model control @@ -3312,6 +3318,8 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, & Model%ttendlim = ttendlim !--- gfdl MP parameters Model%lgfdlmprad = lgfdlmprad +!--- Thompson,GFDL MP parameter + Model%lrefres = lrefres !--- land/surface model parameters Model%lsm = lsm @@ -4310,6 +4318,7 @@ subroutine control_print(Model) print *, ' Thompson microphysical parameters' print *, ' ltaerosol : ', Model%ltaerosol print *, ' lradar : ', Model%lradar + print *, ' lrefres : ', Model%lrefres print *, ' ttendlim : ', Model%ttendlim print *, ' ' endif @@ -4327,6 +4336,7 @@ subroutine control_print(Model) if (Model%imp_physics == Model%imp_physics_gfdl) then print *, ' GFDL microphysical parameters' print *, ' GFDL MP radiation inter: ', Model%lgfdlmprad + print *, ' lrefres : ', Model%lrefres print *, ' ' endif diff --git a/io/FV3GFS_io.F90 b/io/FV3GFS_io.F90 index 3e7b7d2e7..e7f7dbd57 100644 --- a/io/FV3GFS_io.F90 +++ b/io/FV3GFS_io.F90 @@ -1119,15 +1119,10 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain) Sfcprop(nb)%zorll(ix) = Sfcprop(nb)%zorlo(ix) Sfcprop(nb)%zorl(ix) = Sfcprop(nb)%zorlo(ix) Sfcprop(nb)%tsfc(ix) = Sfcprop(nb)%tsfco(ix) - if (Sfcprop(nb)%slmsk(ix) < 0.1 .or. Sfcprop(nb)%slmsk(ix) > 1.9) then + if (Sfcprop(nb)%slmsk(ix) > 1.9) then Sfcprop(nb)%landfrac(ix) = 0.0 - if (Sfcprop(nb)%oro_uf(ix) > 0.01) then - Sfcprop(nb)%lakefrac(ix) = 1.0 ! lake - else - Sfcprop(nb)%lakefrac(ix) = 0.0 ! ocean - endif else - Sfcprop(nb)%landfrac(ix) = 1.0 ! land + Sfcprop(nb)%landfrac(ix) = Sfcprop(nb)%slmsk(ix) endif enddo enddo diff --git a/io/module_write_netcdf.F90 b/io/module_write_netcdf.F90 index 1fce3d8b9..5de1362af 100644 --- a/io/module_write_netcdf.F90 +++ b/io/module_write_netcdf.F90 @@ -22,7 +22,7 @@ module module_write_netcdf contains !---------------------------------------------------------------------------------------- - subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, rc) + subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, idate, rc) ! type(ESMF_FieldBundle), intent(in) :: fieldbundle type(ESMF_FieldBundle), intent(in) :: wrtfb @@ -30,6 +30,7 @@ subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, rc integer, intent(in) :: mpi_comm integer, intent(in) :: mype integer, intent(in) :: im, jm + integer, intent(in) :: idate(7) integer, optional,intent(out) :: rc ! !** local vars @@ -146,7 +147,7 @@ subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, rc call add_dim(ncid, "phalf", phalf_dimid, wrtgrid, rc) end if - call add_dim(ncid, "time", time_dimid, wrtgrid, rc) + call add_dim(ncid, "time", time_dimid, wrtgrid, rc, idate=idate) call get_global_attr(wrtfb, ncid, rc) @@ -503,11 +504,12 @@ subroutine get_grid_attr(grid, prefix, ncid, varid, rc) end subroutine get_grid_attr - subroutine add_dim(ncid, dim_name, dimid, grid, rc) + subroutine add_dim(ncid, dim_name, dimid, grid, idate, rc) integer, intent(in) :: ncid character(len=*), intent(in) :: dim_name - integer, intent(inout) :: dimid + integer, intent(inout) :: dimid type(ESMF_Grid), intent(in) :: grid + integer, intent(in), optional :: idate(7) integer, intent(out) :: rc ! local variable @@ -556,9 +558,38 @@ subroutine add_dim(ncid, dim_name, dimid, grid, rc) call get_grid_attr(grid, dim_name, ncid, dim_varid, rc) + ! if write grid comp changes time units + if ( present (idate) ) then + if ( trim(dim_name) == "time") then + ncerr = nf90_get_att(ncid, dim_varid, 'units', time_units); NC_ERR_STOP(ncerr) + time_units = get_time_units_from_idate(idate) + ncerr = nf90_put_att(ncid, dim_varid, 'units', trim(time_units)); NC_ERR_STOP(ncerr) + endif + endif + end subroutine add_dim ! !---------------------------------------------------------------------------------------- + function get_time_units_from_idate(idate, time_measure) result(time_units) + ! create time units attribute of form 'hours since YYYY-MM-DD HH:MM:SS' + ! from integer array with year,month,day,hour,minute,second + ! optional argument 'time_measure' can be used to change 'hours' to + ! 'days', 'minutes', 'seconds' etc. + character(len=*), intent(in), optional :: time_measure + integer, intent(in) :: idate(7) + character(len=12) :: timechar + character(len=nf90_max_name) :: time_units + if (present(time_measure)) then + timechar = trim(time_measure) + else + timechar = 'hours' + endif + write(time_units,101) idate(1:6) +101 format(' since ',i4.4,'-',i2.2,'-',i2.2,' ',& + i2.2,':',i2.2,':',i2.2) + time_units = trim(adjustl(timechar))//time_units + end function get_time_units_from_idate +! subroutine nccheck(status) use netcdf implicit none diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index d0846be53..f8546614b 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -1386,7 +1386,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) wbeg = MPI_Wtime() call write_netcdf(file_bundle,wrt_int_state%wrtFB(nbdl),trim(filename), & - wrt_mpi_comm,wrt_int_state%mype,imo,jmo,rc) + wrt_mpi_comm,wrt_int_state%mype,imo,jmo,idate,rc) wend = MPI_Wtime() if (lprnt) then write(*,'(A,F10.5,A,I4.2,A,I2.2)')' netcdf Write Time is ',wend-wbeg & @@ -1431,7 +1431,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) wbeg = MPI_Wtime() call write_netcdf(file_bundle,wrt_int_state%wrtFB(nbdl),trim(filename), & - wrt_mpi_comm,wrt_int_state%mype,imo,jmo,rc) + wrt_mpi_comm,wrt_int_state%mype,imo,jmo,idate,rc) wend = MPI_Wtime() if (mype == lead_write_task) then write(*,'(A,F10.5,A,I4.2,A,I2.2)')' netcdf Write Time is ',wend-wbeg & diff --git a/io/post_gfs.F90 b/io/post_gfs.F90 index 7b5a87026..ce6bb281e 100644 --- a/io/post_gfs.F90 +++ b/io/post_gfs.F90 @@ -451,8 +451,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & end do ! ! GFS does not output PD -! pt = 10000. ! this is for 100 hPa added by Moorthi - pt = 0. + pt = ak5(1) ! GFS may not have model derived radar ref. ! TKE @@ -1296,7 +1295,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & endif ! inst incoming sfc longwave - if(trim(fieldname)=='dlwsf') then + if(trim(fieldname)=='dlwrf') then !$omp parallel do private(i,j) do j=jsta,jend do i=ista, iend @@ -1845,6 +1844,16 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & ! time averaged surface clear sky outgoing LW if(trim(fieldname)=='csulf') then + !$omp parallel do private(i,j) + do j=jsta,jend + do i=ista, iend + alwoutc(i,j) = arrayr42d(i,j) + enddo + enddo + endif + + ! time averaged TOA clear sky outgoing LW + if(trim(fieldname)=='csulftoa') then !$omp parallel do private(i,j) do j=jsta,jend do i=ista, iend @@ -1864,7 +1873,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & endif ! time averaged TOA clear sky outgoing SW - if(trim(fieldname)=='csusf') then + if(trim(fieldname)=='csusftoa') then !$omp parallel do private(i,j) do j=jsta,jend do i=ista, iend @@ -2271,7 +2280,6 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & enddo end do -!??? reset pint(lev=1) !$omp parallel do private(i,j) do j=jsta,jend do i=1,im @@ -2282,7 +2290,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & ! print *,'in setvar, pt=',pt,'ak5(lp1)=', ak5(lp1),'ak5(1)=',ak5(1) ! compute alpint - do l=lp1,2,-1 + do l=lp1,1,-1 !$omp parallel do private(i,j) do j=jsta,jend do i=1,im @@ -2321,6 +2329,18 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & enddo enddo +! compute cwm for gfdlmp + if( imp_physics == 11 ) then + do l=1,lm +!$omp parallel do private(i,j) + do j=jsta,jend + do i=ista,iend + cwm(i,j,l)=qqg(i,j,l)+qqs(i,j,l)+qqr(i,j,l)+qqi(i,j,l)+qqw(i,j,l) + enddo + enddo + enddo + endif + ! estimate 2m pres and convert t2m to theta !$omp parallel do private(i,j) do j=jsta,jend From 0e01f028b669d94de540a798f9e019b9eef0a26e Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Wed, 6 Nov 2019 02:12:00 +0000 Subject: [PATCH 02/18] fix the add_dim interface --- io/module_write_netcdf.F90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/io/module_write_netcdf.F90 b/io/module_write_netcdf.F90 index 5de1362af..0a2b63d73 100644 --- a/io/module_write_netcdf.F90 +++ b/io/module_write_netcdf.F90 @@ -504,17 +504,18 @@ subroutine get_grid_attr(grid, prefix, ncid, varid, rc) end subroutine get_grid_attr - subroutine add_dim(ncid, dim_name, dimid, grid, idate, rc) + subroutine add_dim(ncid, dim_name, dimid, grid, rc, idate) integer, intent(in) :: ncid character(len=*), intent(in) :: dim_name integer, intent(inout) :: dimid type(ESMF_Grid), intent(in) :: grid - integer, intent(in), optional :: idate(7) integer, intent(out) :: rc + integer, intent(in), optional :: idate(7) ! local variable integer :: i, attcount, n, dim_varid integer :: ncerr + character(255) :: time_units character(len=ESMF_MAXSTR) :: attName type(ESMF_TypeKind_Flag) :: typekind From d47965cadfa665b0353d675851e3954472fb6857 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Wed, 6 Nov 2019 20:38:00 +0000 Subject: [PATCH 03/18] update gfslmp code changes from Ruiyu --- gfsphysics/physics/gfdl_cloud_microphys.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/gfsphysics/physics/gfdl_cloud_microphys.F90 b/gfsphysics/physics/gfdl_cloud_microphys.F90 index 9a1820465..ba4c814d6 100644 --- a/gfsphysics/physics/gfdl_cloud_microphys.F90 +++ b/gfsphysics/physics/gfdl_cloud_microphys.F90 @@ -4688,7 +4688,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms, real :: n0r = 8.0e6, n0s = 3.0e6, n0g = 4.0e6 real :: alphar = 0.8, alphas = 0.25, alphag = 0.5 real :: gammar = 17.837789, gammas = 8.2850630, gammag = 11.631769 - real :: qmin = 1.0e-12, beta = 1.22 + real :: qmin = 1.0e-12, beta = 1.22, qmin1 = 9.e-6 do k = ks, ke do i = is, ie @@ -4718,7 +4718,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms, ! cloud ice (Heymsfield and Mcfarquhar, 1996) ! ----------------------------------------------------------------------- - if (qmi (i, k) .gt. qmin) then + if (qmi (i, k) .gt. qmin1) then qci (i, k) = dpg * qmi (i, k) * 1.0e3 rei_fac = log (1.0e3 * qmi (i, k) * den (i, k)) if (t (i, k) - tice .lt. - 50) then @@ -4744,7 +4744,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms, ! cloud ice (Wyser, 1998) ! ----------------------------------------------------------------------- - if (qmi (i, k) .gt. qmin) then + if (qmi (i, k) .gt. qmin1) then qci (i, k) = dpg * qmi (i, k) * 1.0e3 bw = - 2. + 1.e-3 * log10 (den (i, k) * qmi (i, k) / rho_0) * max (0.0, tice - t (i, k)) ** 1.5 rei (i, k) = 377.4 + bw * (203.3 + bw * (37.91 + 2.3696 * bw)) @@ -4774,7 +4774,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms, ! snow (Lin et al., 1983) ! ----------------------------------------------------------------------- - if (qms (i, k) .gt. qmin) then + if (qms (i, k) .gt. qmin1) then qcs (i, k) = dpg * qms (i, k) * 1.0e3 lambdas = exp (0.25 * log (pi * rhos * n0s / qms (i, k) / den (i, k))) res (i, k) = 0.5 * exp (log (gammas / 6) / alphas) / lambdas * 1.0e6 From 2d737534a785579b680c5aa773a5bcc0c64c8351 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Thu, 7 Nov 2019 15:10:44 +0000 Subject: [PATCH 04/18] update satmedmfvdifq from Jongil --- gfsphysics/physics/satmedmfvdifq.f | 53 +++++++++++++++--------------- 1 file changed, 27 insertions(+), 26 deletions(-) diff --git a/gfsphysics/physics/satmedmfvdifq.f b/gfsphysics/physics/satmedmfvdifq.f index 1a8be355e..11c047fd0 100644 --- a/gfsphysics/physics/satmedmfvdifq.f +++ b/gfsphysics/physics/satmedmfvdifq.f @@ -151,7 +151,7 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke, & rlmn, rlmn1, rlmx, elmx, & ttend, utend, vtend, qtend, & zfac, zfmin, vk, spdk2, - & tkmin, xkzinv, xkgdx, + & tkmin, tkminx, xkzinv, xkgdx, & zlup, zldn, bsum, & tem, tem1, tem2, & ptem, ptem0, ptem1, ptem2 @@ -176,11 +176,11 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke, parameter(prmin=0.25,prmax=4.0) parameter(pr0=1.0,prtke=1.0,prscu=0.67) parameter(f0=1.e-4,crbmin=0.15,crbmax=0.35) - parameter(tkmin=1.e-9,dspmax=10.0) + parameter(tkmin=1.e-9,tkminx=0.2,dspmax=10.0) parameter(qmin=1.e-8,qlmin=1.e-12,zfmin=1.e-8) parameter(aphi5=5.,aphi16=16.) parameter(elmfac=1.0,elefac=1.0,cql=100.) - parameter(dw2min=1.e-4,dkmax=1000.,xkgdx=25000.) + parameter(dw2min=1.e-4,dkmax=1000.,xkgdx=5000.) parameter(qlcr=3.5e-5,zstblmax=2500.,xkzinv=0.1) parameter(h1=0.33333333) parameter(ck0=0.4,ck1=0.15,ch0=0.4,ch1=0.15) @@ -273,20 +273,20 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke, xkzo(i,k) = 0.0 xkzmo(i,k) = 0.0 if (k < kinver(i)) then -! vertical background diffusivity - ptem = prsi(i,k+1) * tx1(i) - tem1 = 1.0 - ptem - tem2 = tem1 * tem1 * 10.0 - tem2 = min(1.0, exp(-tem2)) - xkzo(i,k) = xkzm_hx(i) * tem2 -! +! minimum turbulent mixing length ptem = prsl(i,k) * tx1(i) tem1 = 1.0 - ptem tem2 = tem1 * tem1 * 2.5 tem2 = min(1.0, exp(-tem2)) rlmnz(i,k)= rlmn * tem2 rlmnz(i,k)= max(rlmnz(i,k), rlmn1) -! vertical background diffusivity for momentum +! vertical background diffusivity + ptem = prsi(i,k+1) * tx1(i) + tem1 = 1.0 - ptem + tem2 = tem1 * tem1 * 10.0 + tem2 = min(1.0, exp(-tem2)) + xkzo(i,k) = xkzm_hx(i) * tem2 +! vertical background diffusivity for momentum if (ptem >= xkzm_s) then xkzmo(i,k) = xkzm_mx(i) kx1(i) = k + 1 @@ -674,20 +674,20 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke, ! ! background diffusivity decreasing with increasing surface layer stability ! - do i = 1, im - if(.not.sfcflg(i)) then - tem = (1. + 5. * rbsoil(i))**2. -! tem = (1. + 5. * zol(i))**2. - frik(i) = 0.1 + 0.9 / tem - endif - enddo -! - do k = 1,km1 - do i=1,im - xkzo(i,k) = frik(i) * xkzo(i,k) - xkzmo(i,k)= frik(i) * xkzmo(i,k) - enddo - enddo +! do i = 1, im +! if(.not.sfcflg(i)) then +! tem = (1. + 5. * rbsoil(i))**2. +!! tem = (1. + 5. * zol(i))**2. +! frik(i) = 0.1 + 0.9 / tem +! endif +! enddo +! +! do k = 1,km1 +! do i=1,im +! xkzo(i,k) = frik(i) * xkzo(i,k) +! xkzmo(i,k)= frik(i) * xkzmo(i,k) +! enddo +! enddo ! ! The background vertical diffusivities in the inversion layers are limited ! to be less than or equal to xkzminv @@ -867,13 +867,14 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke, do i = 1, im if(k == 1) then tem = ckz(i,1) - tem1 = xkzmo(i,1) + tem1 = 0.5 * xkzmo(i,1) else tem = 0.5 * (ckz(i,k-1) + ckz(i,k)) tem1 = 0.5 * (xkzmo(i,k-1) + xkzmo(i,k)) endif ptem = tem1 / (tem * elm(i,k)) tkmnz(i,k) = ptem * ptem + tkmnz(i,k) = min(tkmnz(i,k), tkminx) tkmnz(i,k) = max(tkmnz(i,k), tkmin) enddo enddo From 7e91866d114c17022822946de0bdb6c86dd51de3 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Thu, 7 Nov 2019 20:12:57 +0000 Subject: [PATCH 05/18] update ccpp with gfdlmp change from Ruiyu --- .gitmodules | 3 ++- ccpp/physics | 2 +- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/.gitmodules b/.gitmodules index 04f5ebec0..f64926ad2 100644 --- a/.gitmodules +++ b/.gitmodules @@ -7,4 +7,5 @@ url = https://github.com/NCAR/ccpp-framework [submodule "ccpp/physics"] path = ccpp/physics - url = https://github.com/NCAR/ccpp-physics + url = https://github.com/junwang-noaa/ccpp-physics + branch = gfdlmp_update diff --git a/ccpp/physics b/ccpp/physics index d4b1cd020..594b5db85 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit d4b1cd020f8347147b86d3a18b56c03cb5c57d67 +Subproject commit 594b5db851b98acb68c9b33a9f326b9e203fc1a1 From 3e018222ea70e30b0d5df1234044cfb817ec5e86 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Sun, 10 Nov 2019 00:58:26 +0000 Subject: [PATCH 06/18] update idate and nfhours for iau --- io/module_write_netcdf.F90 | 40 +++------------------ io/module_wrt_grid_comp.F90 | 70 ++++++++++++++++++++++--------------- 2 files changed, 46 insertions(+), 64 deletions(-) diff --git a/io/module_write_netcdf.F90 b/io/module_write_netcdf.F90 index 0a2b63d73..1fce3d8b9 100644 --- a/io/module_write_netcdf.F90 +++ b/io/module_write_netcdf.F90 @@ -22,7 +22,7 @@ module module_write_netcdf contains !---------------------------------------------------------------------------------------- - subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, idate, rc) + subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, rc) ! type(ESMF_FieldBundle), intent(in) :: fieldbundle type(ESMF_FieldBundle), intent(in) :: wrtfb @@ -30,7 +30,6 @@ subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, id integer, intent(in) :: mpi_comm integer, intent(in) :: mype integer, intent(in) :: im, jm - integer, intent(in) :: idate(7) integer, optional,intent(out) :: rc ! !** local vars @@ -147,7 +146,7 @@ subroutine write_netcdf(fieldbundle, wrtfb, filename, mpi_comm, mype, im, jm, id call add_dim(ncid, "phalf", phalf_dimid, wrtgrid, rc) end if - call add_dim(ncid, "time", time_dimid, wrtgrid, rc, idate=idate) + call add_dim(ncid, "time", time_dimid, wrtgrid, rc) call get_global_attr(wrtfb, ncid, rc) @@ -504,18 +503,16 @@ subroutine get_grid_attr(grid, prefix, ncid, varid, rc) end subroutine get_grid_attr - subroutine add_dim(ncid, dim_name, dimid, grid, rc, idate) + subroutine add_dim(ncid, dim_name, dimid, grid, rc) integer, intent(in) :: ncid character(len=*), intent(in) :: dim_name - integer, intent(inout) :: dimid + integer, intent(inout) :: dimid type(ESMF_Grid), intent(in) :: grid integer, intent(out) :: rc - integer, intent(in), optional :: idate(7) ! local variable integer :: i, attcount, n, dim_varid integer :: ncerr - character(255) :: time_units character(len=ESMF_MAXSTR) :: attName type(ESMF_TypeKind_Flag) :: typekind @@ -559,38 +556,9 @@ subroutine add_dim(ncid, dim_name, dimid, grid, rc, idate) call get_grid_attr(grid, dim_name, ncid, dim_varid, rc) - ! if write grid comp changes time units - if ( present (idate) ) then - if ( trim(dim_name) == "time") then - ncerr = nf90_get_att(ncid, dim_varid, 'units', time_units); NC_ERR_STOP(ncerr) - time_units = get_time_units_from_idate(idate) - ncerr = nf90_put_att(ncid, dim_varid, 'units', trim(time_units)); NC_ERR_STOP(ncerr) - endif - endif - end subroutine add_dim ! !---------------------------------------------------------------------------------------- - function get_time_units_from_idate(idate, time_measure) result(time_units) - ! create time units attribute of form 'hours since YYYY-MM-DD HH:MM:SS' - ! from integer array with year,month,day,hour,minute,second - ! optional argument 'time_measure' can be used to change 'hours' to - ! 'days', 'minutes', 'seconds' etc. - character(len=*), intent(in), optional :: time_measure - integer, intent(in) :: idate(7) - character(len=12) :: timechar - character(len=nf90_max_name) :: time_units - if (present(time_measure)) then - timechar = trim(time_measure) - else - timechar = 'hours' - endif - write(time_units,101) idate(1:6) -101 format(' since ',i4.4,'-',i2.2,'-',i2.2,' ',& - i2.2,':',i2.2,':',i2.2) - time_units = trim(adjustl(timechar))//time_units - end function get_time_units_from_idate -! subroutine nccheck(status) use netcdf implicit none diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index f8546614b..49dfdcff3 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -71,6 +71,7 @@ module module_wrt_grid_comp logical,save :: first_init=.false. logical,save :: first_run=.false. logical,save :: first_getlatlon=.true. + logical,save :: change_wrtidate=.false. ! !----------------------------------------------------------------------- ! @@ -173,7 +174,8 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) real(ESMF_KIND_R4) :: valueR4 real(ESMF_KIND_R8) :: valueR8 - integer :: attCount, axeslen, jidx, noutfile + integer :: attCount, axeslen, jidx, idx, noutfile + character(19) :: newdate character(128) :: FBlist_outfilename(100), outfile_name character(128),dimension(:,:), allocatable :: outfilename real(8), dimension(:), allocatable :: slat @@ -183,7 +185,6 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) real(ESMF_KIND_R8) :: geo_lon, geo_lat real(ESMF_KIND_R8) :: lon1_r8, lat1_r8 real(ESMF_KIND_R8) :: x1, y1, x, y - type(ESMF_Time) :: IO_BASETIME_IAU type(ESMF_TimeInterval) :: IAU_offsetTI type(ESMF_DataCopy_Flag) :: copyflag=ESMF_DATACOPY_REFERENCE ! real(8),parameter :: PI=3.14159265358979d0 @@ -469,8 +470,31 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) call ESMF_Finalize(endflag=ESMF_END_ABORT) endif +! +!----------------------------------------------------------------------- +!*** get write grid component initial time from clock +!----------------------------------------------------------------------- +! + call ESMF_ClockGet(clock =CLOCK & !<-- The ESMF Clock + ,startTime=wrt_int_state%IO_BASETIME & !<-- The Clock's starting time + ,rc =RC) - + call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), & + m=idate(5),s=idate(6),rc=rc) +! if (lprnt) write(0,*) 'in wrt initial, io_baseline time=',idate,'rc=',rc + idate(7) = 1 + wrt_int_state%idate = idate + wrt_int_state%fdate = idate + if(iau_offset > 0 ) then + call ESMF_TimeIntervalSet(IAU_offsetTI, h=iau_offset, rc=rc) + wrt_int_state%IO_BASETIME = wrt_int_state%IO_BASETIME + IAU_offsetTI + call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), & + m=idate(5),s=idate(6),rc=rc) + wrt_int_state%idate = idate + change_wrtidate = .true. + if (mype == lead_write_task) print *,'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc + endif +! ! Create field bundle !------------------------------------------------------------------- ! @@ -867,6 +891,17 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return +! add iau here + if ( index(trim(attNameList(i)),'time:units')>0) then + if ( change_wrtidate ) then + idx = index(trim(valueS),' since ') + if(mype == lead_write_task) print *,'in write grid comp, time:unit=',trim(valueS) + write(newdate,'(I4.4,a,I2.2,a,I2.2,a,I2.2,a,I2.2,a,I2.2)') idate(1),'-', & + idate(2),'-',idate(3),' ',idate(4),':',idate(5),':',idate(6) + valueS = valueS(1:idx+6)//newdate + if(mype == lead_write_task) print *,'in write grid comp, new time:unit=',trim(valueS) + endif + endif call ESMF_AttributeSet(wrtgrid, convention="NetCDF", purpose="FV3", & name=trim(attNameList(i)), value=valueS, rc=rc) @@ -1036,28 +1071,6 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return ! !----------------------------------------------------------------------- -!*** SET THE IO_BaseTime TO THE INITIAL CLOCK TIME. -!----------------------------------------------------------------------- -! - call ESMF_ClockGet(clock =CLOCK & !<-- The ESMF Clock - ,startTime=wrt_int_state%IO_BASETIME & !<-- The Clock's starting time - ,rc =RC) - - call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), & - m=idate(5),s=idate(6),rc=rc) -! if (lprnt) write(0,*) 'in wrt initial, io_baseline time=',idate,'rc=',rc - idate(7) = 1 - wrt_int_state%idate = idate - wrt_int_state%fdate = idate - if(iau_offset > 0 ) then - call ESMF_TimeIntervalSet(IAU_offsetTI, h=iau_offset, rc=rc) - IO_BASETIME_IAU = wrt_int_state%IO_BASETIME + IAU_offsetTI - call ESMF_TimeGet(time=IO_BASETIME_IAU,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), & - m=idate(5),s=idate(6),rc=rc) -! if (lprnt) write(0,*) 'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc - endif -! -!----------------------------------------------------------------------- !*** SET THE FIRST HISTORY FILE'S TIME INDEX. !----------------------------------------------------------------------- ! @@ -1251,9 +1264,10 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) ! nf_seconds = nf_hours*3600+nf_minuteS*60+nseconds+real(nseconds_num)/real(nseconds_den) ! shift forecast hour by iau_offset if iau is on. - nf_seconds = nf_seconds - iau_offset*3600 + !nf_seconds = nf_seconds - iau_offset*3600 wrt_int_state%nfhour = nf_seconds/3600. nf_hours = int(nf_seconds/3600.) + if(mype == lead_write_task) print *,'in write grid comp, nf_hours=',nf_hours ! if iau_offset > nf_hours, don't write out anything if (nf_hours < 0) return @@ -1386,7 +1400,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) wbeg = MPI_Wtime() call write_netcdf(file_bundle,wrt_int_state%wrtFB(nbdl),trim(filename), & - wrt_mpi_comm,wrt_int_state%mype,imo,jmo,idate,rc) + wrt_mpi_comm,wrt_int_state%mype,imo,jmo,rc) wend = MPI_Wtime() if (lprnt) then write(*,'(A,F10.5,A,I4.2,A,I2.2)')' netcdf Write Time is ',wend-wbeg & @@ -1431,7 +1445,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) wbeg = MPI_Wtime() call write_netcdf(file_bundle,wrt_int_state%wrtFB(nbdl),trim(filename), & - wrt_mpi_comm,wrt_int_state%mype,imo,jmo,idate,rc) + wrt_mpi_comm,wrt_int_state%mype,imo,jmo,rc) wend = MPI_Wtime() if (mype == lead_write_task) then write(*,'(A,F10.5,A,I4.2,A,I2.2)')' netcdf Write Time is ',wend-wbeg & From 94842bcc0b5798963ab9e0cd74566d930605a287 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Sun, 10 Nov 2019 04:11:51 +0000 Subject: [PATCH 07/18] fix print line --- io/module_wrt_grid_comp.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index 49dfdcff3..3439e89a4 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -492,7 +492,7 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) m=idate(5),s=idate(6),rc=rc) wrt_int_state%idate = idate change_wrtidate = .true. - if (mype == lead_write_task) print *,'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc + if (lprnt) print *,'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc endif ! ! Create field bundle @@ -895,11 +895,11 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) if ( index(trim(attNameList(i)),'time:units')>0) then if ( change_wrtidate ) then idx = index(trim(valueS),' since ') - if(mype == lead_write_task) print *,'in write grid comp, time:unit=',trim(valueS) + if(lprnt) print *,'in write grid comp, time:unit=',trim(valueS) write(newdate,'(I4.4,a,I2.2,a,I2.2,a,I2.2,a,I2.2,a,I2.2)') idate(1),'-', & idate(2),'-',idate(3),' ',idate(4),':',idate(5),':',idate(6) valueS = valueS(1:idx+6)//newdate - if(mype == lead_write_task) print *,'in write grid comp, new time:unit=',trim(valueS) + if(lprnt) print *,'in write grid comp, new time:unit=',trim(valueS) endif endif call ESMF_AttributeSet(wrtgrid, convention="NetCDF", purpose="FV3", & From b71e959dc6e9f96359472ee8be267905b70f0143 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Sun, 10 Nov 2019 20:24:53 +0000 Subject: [PATCH 08/18] initialize ihrst for ihrst --- io/post_gfs.F90 | 28 +++++++++------------------- 1 file changed, 9 insertions(+), 19 deletions(-) diff --git a/io/post_gfs.F90 b/io/post_gfs.F90 index ce6bb281e..d7a2a0d42 100644 --- a/io/post_gfs.F90 +++ b/io/post_gfs.F90 @@ -28,7 +28,7 @@ subroutine post_run_gfs(wrt_int_state,mypei,mpicomp,lead_write, & ! use ctlblk_mod, only : komax,ifhr,ifmin,modelname,datapd,fld_info, & npset,grib,gocart_on,icount_calmict, jsta, & - jend,im, nsoil + jend,im, nsoil, ihrst use gridspec_mod, only : maptype, gridtype use grib2_module, only : gribit2,num_pset,nrecout,first_grbtbl ! @@ -369,7 +369,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & ! integer i, ip1, j, l, k, n, iret, ibdl, rc, kstart, kend integer ista,iend,fieldDimCount,gridDimCount,ncount_field - integer idate(8), jdate(8) + integer jdate(8) logical foundland, foundice, found real(4) rinc(5) real tlmh,RADI,TMP,ES,TV,RHOAIR,tem,tstart,dtp @@ -641,26 +641,16 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & enddo ! ! get inital date - idate = 0 - idate(1) = wrt_int_state%idate(1) - idate(2) = wrt_int_state%idate(2) - idate(3) = wrt_int_state%idate(3) - idate(5) = wrt_int_state%idate(4) - idate(6) = wrt_int_state%idate(5) sdat(1) = wrt_int_state%idate(2) !month sdat(2) = wrt_int_state%idate(3) !day sdat(3) = wrt_int_state%idate(1) !year - jdate = 0 - jdate(1) = wrt_int_state%fdate(1) - jdate(2) = wrt_int_state%fdate(2) - jdate(3) = wrt_int_state%fdate(3) !jdate(4): time zone - jdate(5) = wrt_int_state%fdate(4) - jdate(6) = wrt_int_state%fdate(5) - idat(1) = wrt_int_state%idate(2) - idat(2) = wrt_int_state%idate(3) - idat(3) = wrt_int_state%idate(1) - idat(4) = IFHR - idat(5) = IFMIN + ihrst = wrt_int_state%idate(4) !hour + + idat(1) = wrt_int_state%fdate(2) + idat(2) = wrt_int_state%fdate(3) + idat(3) = wrt_int_state%fdate(1) + idat(4) = wrt_int_state%fdate(4) + idat(5) = wrt_int_state%fdate(5) ! ! if(mype==0) print *,'jdate=',jdate,'idate=',idate,'sdat=',sdat ! CALL W3DIFDAT(JDATE,IDATE,0,RINC) From 2eb8b2cd7c63bcaae1157ccf9ef418b67933d9c3 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Mon, 11 Nov 2019 02:56:16 +0000 Subject: [PATCH 09/18] initialize ihrst for inline post --- io/post_gfs.F90 | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/io/post_gfs.F90 b/io/post_gfs.F90 index d7a2a0d42..3ea19c555 100644 --- a/io/post_gfs.F90 +++ b/io/post_gfs.F90 @@ -28,7 +28,7 @@ subroutine post_run_gfs(wrt_int_state,mypei,mpicomp,lead_write, & ! use ctlblk_mod, only : komax,ifhr,ifmin,modelname,datapd,fld_info, & npset,grib,gocart_on,icount_calmict, jsta, & - jend,im, nsoil, ihrst + jend,im, nsoil use gridspec_mod, only : maptype, gridtype use grib2_module, only : gribit2,num_pset,nrecout,first_grbtbl ! @@ -342,7 +342,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & tprec, tclod, trdlw, trdsw, tsrfc, tmaxmin, theat, & ardlw, ardsw, asrfc, avrain, avcnvc, iSF_SURFACE_PHYSICS,& td3d, idat, sdat, ifhr, ifmin, dt, nphs, dtq2, pt_tbl, & - alsl, spl + alsl, spl, ihrst use params_mod, only: erad, dtr, capa, p1000 use gridspec_mod,only: latstart, latlast, lonstart, lonlast, cenlon, cenlat use lookup_mod, only: thl, plq, ptbl, ttbl, rdq, rdth, rdp, rdthe, pl, & @@ -652,7 +652,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & idat(4) = wrt_int_state%fdate(4) idat(5) = wrt_int_state%fdate(5) ! -! if(mype==0) print *,'jdate=',jdate,'idate=',idate,'sdat=',sdat + if(mype==0) print *,'idat=',idat,'sdat=',sdat,'ihrst=',ihrst ! CALL W3DIFDAT(JDATE,IDATE,0,RINC) ! ! if(mype==0)print *,' rinc=',rinc @@ -690,7 +690,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & call ESMF_FieldBundleGet(wrt_int_state%wrtFB(ibdl),fieldName='land',isPresent=found, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) return ! bail out - if(mype==0) print *,'ibdl=',ibdl,'land, found=',found +! if(mype==0) print *,'ibdl=',ibdl,'land, found=',found if (found) then call ESMF_FieldBundleGet(wrt_int_state%wrtFB(ibdl),'land',field=theField, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & @@ -714,7 +714,7 @@ subroutine set_postvars_gfs(wrt_int_state,mpicomp,setvar_atmfile, & call ESMF_FieldBundleGet(wrt_int_state%wrtFB(ibdl),'icec',isPresent=found, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) return ! bail out - if(mype==0) print *,'ibdl=',ibdl,'ice, found=',found +! if(mype==0) print *,'ibdl=',ibdl,'ice, found=',found if (found) then call ESMF_FieldBundleGet(wrt_int_state%wrtFB(ibdl),'icec',field=theField, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & From 248196cd055c67b0b16c1cbc2ea946b8f385d027 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Tue, 12 Nov 2019 19:41:30 +0000 Subject: [PATCH 10/18] add comments for write grid comp idate change when iau is enabled --- io/module_wrt_grid_comp.F90 | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index 3439e89a4..405af2841 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -485,6 +485,7 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) idate(7) = 1 wrt_int_state%idate = idate wrt_int_state%fdate = idate +! update IO-BASETIME and idate on write grid comp when IAU is enabled if(iau_offset > 0 ) then call ESMF_TimeIntervalSet(IAU_offsetTI, h=iau_offset, rc=rc) wrt_int_state%IO_BASETIME = wrt_int_state%IO_BASETIME + IAU_offsetTI @@ -891,7 +892,7 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return -! add iau here +! update the time:units when idate on write grid component is changed if ( index(trim(attNameList(i)),'time:units')>0) then if ( change_wrtidate ) then idx = index(trim(valueS),' since ') @@ -1263,8 +1264,6 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) ! 'nseconds_num=',nseconds_num,nseconds_den,'mype=',mype ! nf_seconds = nf_hours*3600+nf_minuteS*60+nseconds+real(nseconds_num)/real(nseconds_den) - ! shift forecast hour by iau_offset if iau is on. - !nf_seconds = nf_seconds - iau_offset*3600 wrt_int_state%nfhour = nf_seconds/3600. nf_hours = int(nf_seconds/3600.) if(mype == lead_write_task) print *,'in write grid comp, nf_hours=',nf_hours From bddc70cf49acd4e39213120be38f116750fe5804 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Fri, 15 Nov 2019 15:15:45 +0000 Subject: [PATCH 11/18] update fv3 dycore and ccpp physics repo --- .gitmodules | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/.gitmodules b/.gitmodules index f64926ad2..81687f8ab 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,11 +1,9 @@ [submodule "atmos_cubed_sphere"] path = atmos_cubed_sphere - url = https://github.com/junwang-noaa/GFDL_atmos_cubed_sphere - branch = regfv3_rst + url = https://github.com/NOAA-EMC/GFDL_atmos_cubed_sphere [submodule "ccpp/framework"] path = ccpp/framework url = https://github.com/NCAR/ccpp-framework [submodule "ccpp/physics"] path = ccpp/physics - url = https://github.com/junwang-noaa/ccpp-physics - branch = gfdlmp_update + url = https://github.com/NCAR/ccpp-physics From 97d12bc6594854674c34b837bfbbf7ae408785a9 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Fri, 15 Nov 2019 15:35:19 +0000 Subject: [PATCH 12/18] checkout out dev/emc from EMC fv3 dycore repo and master branch in NCAR ccpp physics repo --- .gitmodules | 1 + atmos_cubed_sphere | 2 +- ccpp/physics | 2 +- 3 files changed, 3 insertions(+), 2 deletions(-) diff --git a/.gitmodules b/.gitmodules index 81687f8ab..949d298df 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,6 +1,7 @@ [submodule "atmos_cubed_sphere"] path = atmos_cubed_sphere url = https://github.com/NOAA-EMC/GFDL_atmos_cubed_sphere + branch = dev/emc [submodule "ccpp/framework"] path = ccpp/framework url = https://github.com/NCAR/ccpp-framework diff --git a/atmos_cubed_sphere b/atmos_cubed_sphere index cc28aa5aa..f967c71b8 160000 --- a/atmos_cubed_sphere +++ b/atmos_cubed_sphere @@ -1 +1 @@ -Subproject commit cc28aa5aa94817b24aea033cc19c1792e5927459 +Subproject commit f967c71b8f59c22c86de2ead1074e85e3ccf97b4 diff --git a/ccpp/physics b/ccpp/physics index 594b5db85..c0bec7173 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit 594b5db851b98acb68c9b33a9f326b9e203fc1a1 +Subproject commit c0bec7173d10ca1460986a8dea81681736a26e56 From 5bfd3c607241b71dd139a5367ce70ede6607af38 Mon Sep 17 00:00:00 2001 From: "Jun.Wang" Date: Fri, 15 Nov 2019 17:57:09 +0000 Subject: [PATCH 13/18] update submodules atmos_cubed_sphere and ccpp/physics for fv3atm #8 --- atmos_cubed_sphere | 2 +- ccpp/physics | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/atmos_cubed_sphere b/atmos_cubed_sphere index f967c71b8..452333aff 160000 --- a/atmos_cubed_sphere +++ b/atmos_cubed_sphere @@ -1 +1 @@ -Subproject commit f967c71b8f59c22c86de2ead1074e85e3ccf97b4 +Subproject commit 452333aff37f4c4348504bb9d485d98072e984e5 diff --git a/ccpp/physics b/ccpp/physics index c0bec7173..904a43354 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit c0bec7173d10ca1460986a8dea81681736a26e56 +Subproject commit 904a43354bc5922558dffc9b88acbc854c20e605 From ea91f650177a326846a0e257772d7249804f0df7 Mon Sep 17 00:00:00 2001 From: Bin Liu Date: Mon, 25 Nov 2019 14:56:09 -0500 Subject: [PATCH 14/18] The HAFS related developments for the write_grid_component (#10) * output nested domain on cubed sphere grid * output native cubed sphere grid for regional fv3 and nested fv3 * add HWRF PBL and surface drag for HAFS (sfc_diff.f, moninedmf.f, GFS_typedefs.F90, GFS_physics_driver.F90 modified) * add back HWRF PBL and sfc drag after master revert fractional land/sea mask * Enable outputing a grid from the write grid component, which is larger than the native computation grid of the nested/regional domain. a. Remove the limitation of the output grid being completely inside the model native compuation grid. b. Mask out the output grid points outside the computation grid with missing values. * update HAFS with HWRF PBL/sfc drag * add option to switch GWD on/off * support/HAFS: Clean up before merging features back the develop branch. * support/HAFS: Change the missing value in FV3GFS_io to default real type. --- fv3_cap.F90 | 7 +- io/FV3GFS_io.F90 | 3 +- io/module_fv3_io_def.F90 | 1 + io/module_wrt_grid_comp.F90 | 336 ++++++++++++++++++++++++++++++++++-- 4 files changed, 328 insertions(+), 19 deletions(-) diff --git a/fv3_cap.F90 b/fv3_cap.F90 index e8e482099..d04abdce2 100644 --- a/fv3_cap.F90 +++ b/fv3_cap.F90 @@ -32,7 +32,7 @@ module fv3gfs_cap_mod cplprint_flag,output_1st_tstep_rst, & first_kdt - use module_fv3_io_def, only: num_pes_fcst,write_groups, & + use module_fv3_io_def, only: num_pes_fcst,write_groups,app_domain, & num_files, filename_base, & wrttasks_per_group, n_group, & lead_wrttask, last_wrttask, & @@ -322,6 +322,10 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) label ='write_tasks_per_group:',rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return + CALL ESMF_ConfigGetAttribute(config=CF,value=app_domain, default="global", & + label ='app_domain:',rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return + if(mype == 0) print *,'af nems config,restart_interval=',restart_interval, & 'quilting=',quilting,'write_groups=',write_groups,wrttasks_per_group, & 'calendar=',trim(calendar),'calendar_type=',calendar_type @@ -692,6 +696,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) isrctermprocessing = 1 call ESMF_FieldBundleRegridStore(fcstFB(j), wrtFB(j,i), & regridMethod=regridmethod, routehandle=routehandle(j,i), & + unmappedaction=ESMF_UNMAPPEDACTION_IGNORE, & srcTermProcessing=isrctermprocessing, rc=rc) ! if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return diff --git a/io/FV3GFS_io.F90 b/io/FV3GFS_io.F90 index e7f7dbd57..4b2d1e426 100644 --- a/io/FV3GFS_io.F90 +++ b/io/FV3GFS_io.F90 @@ -101,8 +101,7 @@ module FV3GFS_io_mod logical :: uwork_set = .false. character(128) :: uwindname integer, parameter, public :: DIAG_SIZE = 500 -! real(kind=kind_phys), parameter :: missing_value = 1.d30 - real(kind=kind_phys), parameter :: missing_value = 9.99e20 + real, parameter :: missing_value = 9.99e20 real, parameter:: stndrd_atmos_ps = 101325. real, parameter:: stndrd_atmos_lapse = 0.0065 diff --git a/io/module_fv3_io_def.F90 b/io/module_fv3_io_def.F90 index d039ccc73..17e9f308d 100644 --- a/io/module_fv3_io_def.F90 +++ b/io/module_fv3_io_def.F90 @@ -15,6 +15,7 @@ module module_fv3_io_def logical :: write_nemsioflip logical :: write_fsyncflag integer :: num_files + character(255) :: app_domain character(255) :: output_grid character(255) :: output_file integer :: imo,jmo diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index 405af2841..aad3991d0 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -32,7 +32,7 @@ module module_wrt_grid_comp use esmf use write_internal_state use module_fv3_io_def, only : num_pes_fcst,lead_wrttask, last_wrttask, & - n_group, num_files, & + n_group, num_files, app_domain, & filename_base, output_grid, output_file, & imo, jmo, write_nemsioflip, & nsout => nsout_io, & @@ -71,6 +71,7 @@ module module_wrt_grid_comp logical,save :: first_init=.false. logical,save :: first_run=.false. logical,save :: first_getlatlon=.true. + logical,save :: first_getmaskwrt=.true. !<-- for mask the output grid of the write comp logical,save :: change_wrtidate=.false. ! !----------------------------------------------------------------------- @@ -147,10 +148,12 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) integer :: ISTAT, tl, i, j, n, k integer,dimension(2,6) :: decomptile + integer,dimension(2) :: regDecomp !define delayout for the nest grid integer :: fieldCount integer :: vm_mpi_comm character(40) :: fieldName, axesname,longname type(ESMF_Config) :: cf + type(ESMF_DELayout) :: delayout type(ESMF_Grid) :: wrtGrid, fcstGrid type(ESMF_Array) :: array_work, array type(ESMF_FieldBundle) :: fieldbdl_work @@ -163,6 +166,7 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) type(ESMF_TypeKind_Flag) :: typekind character(len=80), allocatable :: fieldnamelist(:) integer :: fieldDimCount, gridDimCount + integer, allocatable :: petMap(:) integer, allocatable :: gridToFieldMap(:) integer, allocatable :: ungriddedLBound(:) integer, allocatable :: ungriddedUBound(:) @@ -196,6 +200,7 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) logical :: lprnt !test integer myattCount + real(ESMF_KIND_R8),dimension(:,:), pointer :: glatPtr, glonPtr ! !----------------------------------------------------------------------- !*********************************************************************** @@ -279,24 +284,47 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc) if ( trim(output_grid) == 'cubed_sphere_grid' ) then mytile = mod(wrt_int_state%mype,ntasks)+1 - do tl=1,6 - decomptile(1,tl) = 1 - decomptile(2,tl) = jidx - enddo - - call ESMF_AttributeGet(imp_state_write, convention="NetCDF", purpose="FV3", & - name="gridfile", value=gridfile, rc=rc) + if ( trim(app_domain) == 'global' ) then + do tl=1,6 + decomptile(1,tl) = 1 + decomptile(2,tl) = jidx + enddo + call ESMF_AttributeGet(imp_state_write, convention="NetCDF", purpose="FV3", & + name="gridfile", value=gridfile, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return - if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return - - CALL ESMF_LogWrite("wrtComp: gridfile:"//trim(gridfile),ESMF_LOGMSG_INFO,rc=rc) + CALL ESMF_LogWrite("wrtComp: gridfile:"//trim(gridfile),ESMF_LOGMSG_INFO,rc=rc) + wrtgrid = ESMF_GridCreateMosaic(filename="INPUT/"//trim(gridfile), & + regDecompPTile=decomptile,tileFilePath="INPUT/", & + staggerlocList=(/ESMF_STAGGERLOC_CENTER, ESMF_STAGGERLOC_CORNER/), & + name='wrt_grid', rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return + else + if(trim(app_domain) == 'nested') then + gridfile='grid.nest02.tile7.nc' + else if(trim(app_domain) == 'regional') then + gridfile='grid.tile7.halo0.nc' + endif + regDecomp(1) = 1 + regDecomp(2) = ntasks + allocate(petMap(ntasks)) + do i=1, ntasks + petMap(i) = i-1 + enddo + delayout = ESMF_DELayoutCreate(petMap=petMap, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return - wrtgrid = ESMF_GridCreateMosaic(filename="INPUT/"//trim(gridfile), & - regDecompPTile=decomptile,tileFilePath="INPUT/", & - staggerlocList=(/ESMF_STAGGERLOC_CENTER, ESMF_STAGGERLOC_CORNER/), & - name='wrt_grid', rc=rc) + ! create the nest Grid by reading it from file but use DELayout + wrtGrid = ESMF_GridCreate(filename="INPUT/"//trim(gridfile), & + fileformat=ESMF_FILEFORMAT_GRIDSPEC, regDecomp=regDecomp, & + delayout=delayout, isSphere=.false., indexflag=ESMF_INDEX_DELOCAL, & + rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return - if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return + print *,'in nested/regional cubed_sphere grid, regDecomp=',regDecomp,' PetMap=',petMap(1),petMap(ntasks), & + 'gridfile=',trim(gridfile) + deallocate(petMap) + endif else if ( trim(output_grid) == 'gaussian_grid') then wrtgrid = ESMF_GridCreate1PeriDim(minIndex=(/1,1/), & @@ -1440,6 +1468,19 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc) trim(output_grid) == 'rotated_latlon' .or. & trim(output_grid) == 'lambert_conformal') then + !mask fields according to sfc pressure + !if (mype == lead_write_task) print *,'before mask_fields' + wbeg = MPI_Wtime() + call ESMF_LogWrite("before mask_fields for wrt field bundle", ESMF_LOGMSG_INFO, rc=rc) + !call mask_fields(wrt_int_state%wrtFB(nbdl),rc) + call mask_fields(file_bundle,rc) + !if (mype == lead_write_task) print *,'after mask_fields' + call ESMF_LogWrite("after mask_fields for wrt field bundle", ESMF_LOGMSG_INFO, rc=rc) + wend = MPI_Wtime() + if (mype == lead_write_task) then + write(*,'(A,F10.5,A,I4.2,A,I2.2)')' mask_fields time is ',wend-wbeg + endif + if (trim(output_file) == 'netcdf') then wbeg = MPI_Wtime() @@ -1749,6 +1790,269 @@ subroutine recover_fields(file_bundle,rc) end subroutine recover_fields ! !----------------------------------------------------------------------- +! + subroutine mask_fields(file_bundle,rc) + + type(ESMF_FieldBundle), intent(in) :: file_bundle + integer, intent(out), optional :: rc +! + integer i,j,k,ifld,fieldCount,nstt,nend,fieldDimCount,gridDimCount + integer istart,iend,jstart,jend,kstart,kend,km + type(ESMF_Grid) fieldGrid + type(ESMF_TypeKind_Flag) typekind + type(ESMF_TypeKind_Flag) attTypeKind + character(len=ESMF_MAXSTR) fieldName + type(ESMF_Field), allocatable :: fcstField(:) + real(ESMF_KIND_R4), dimension(:,:), pointer :: var2dPtr2dr4 + real(ESMF_KIND_R4), dimension(:,:,:), pointer :: var3dPtr3dr4 + real(ESMF_KIND_R4), dimension(:,:,:), pointer :: vect3dPtr2dr4 + real(ESMF_KIND_R4), dimension(:,:,:,:), pointer :: vect4dPtr3dr4 + real(ESMF_KIND_R4), dimension(:,:), allocatable :: maskwrt + + logical :: mvispresent=.false. + real(ESMF_KIND_R4) :: missing_value_r4=-1.e+10 + real(ESMF_KIND_R8) :: missing_value_r8=9.99e20 + character(len=ESMF_MAXSTR) :: msg + + save maskwrt + + call ESMF_LogWrite("call mask field on wrt comp",ESMF_LOGMSG_INFO,rc=RC) + +! get fieldCount + call ESMF_FieldBundleGet(file_bundle, fieldCount=fieldCount, & + grid=fieldGrid, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, file=__FILE__)) return ! bail out +! get gridDimCount + call ESMF_GridGet(fieldgrid, dimCount=gridDimCount, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, file=__FILE__)) return ! bail out + + allocate(fcstField(fieldCount)) + call ESMF_LogWrite("call mask field get fcstField",ESMF_LOGMSG_INFO,rc=RC) + call ESMF_FieldBundleGet(file_bundle, fieldList=fcstField, itemorderflag=ESMF_ITEMORDER_ADDORDER, rc=rc) + +! generate the maskwrt according to surface pressure + if( first_getmaskwrt ) then + + do ifld=1,fieldCount + !call ESMF_LogWrite("call mask field get fieldname, type dimcount",ESMF_LOGMSG_INFO,rc=RC) + call ESMF_FieldGet(fcstField(ifld),name=fieldName,typekind=typekind,dimCount=fieldDimCount, rc=rc) + !write(msg,*) 'fieldName,typekind,fieldDimCount=',trim(fieldName),typekind,fieldDimCount + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + if (.not. allocated(maskwrt)) then + if ( typekind == ESMF_TYPEKIND_R4 .and. fieldDimCount == gridDimCount) then + call ESMF_FieldGet(fcstField(ifld),localDe=0, farrayPtr=var2dPtr2dr4, rc=rc) + istart = lbound(var2dPtr2dr4,1) + iend = ubound(var2dPtr2dr4,1) + jstart = lbound(var2dPtr2dr4,2) + jend = ubound(var2dPtr2dr4,2) + allocate(maskwrt(istart:iend,jstart:jend)) + maskwrt(istart:iend,jstart:jend)=1.0 + endif + endif + if(index(trim(fieldName),"pressfc")>0) then + call ESMF_FieldGet(fcstField(ifld),localDe=0, farrayPtr=var2dPtr2dr4, rc=rc) + istart = lbound(var2dPtr2dr4,1) + iend = ubound(var2dPtr2dr4,1) + jstart = lbound(var2dPtr2dr4,2) + jend = ubound(var2dPtr2dr4,2) + if (.not. allocated(maskwrt)) then + allocate(maskwrt(istart:iend,jstart:jend)) + maskwrt(istart:iend,jstart:jend)=1.0 + endif +!$omp parallel do default(shared) private(i,j) + do j=jstart, jend + do i=istart, iend + if(abs(var2dPtr2dr4(i,j)-0.) < 1.0e-6) maskwrt(i,j)=0. + enddo + enddo + call ESMF_LogWrite("call mask field pressfc found, maskwrt generated",ESMF_LOGMSG_INFO,rc=RC) + exit + endif + enddo + first_getmaskwrt = .false. + + endif !first_getmaskwrt + +! loop to mask all fields according to maskwrt + do ifld=1,fieldCount + !call ESMF_LogWrite("call mask field get fieldname, type dimcount",ESMF_LOGMSG_INFO,rc=RC) + call ESMF_FieldGet(fcstField(ifld),name=fieldName,typekind=typekind,dimCount=fieldDimCount, rc=rc) + !write(msg,*) 'fieldName,typekind,fieldDimCount=',trim(fieldName),typekind,fieldDimCount + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + ! For vector fields + if(index(trim(fieldName),"vector")>0) then + ! Only work on ESMF_TYPEKIND_R4 fields for now + if ( typekind == ESMF_TYPEKIND_R4 ) then + ! 3-d vector fields with 4-d arrays + if( fieldDimCount > gridDimCount+1 ) then + !call ESMF_LogWrite("call mask field get vector 3d farray",ESMF_LOGMSG_INFO,rc=RC) + call ESMF_FieldGet(fcstField(ifld), localDe=0, farrayPtr=vect4dPtr3dr4, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, file=__FILE__)) return ! bail out + if( ubound(vect4dPtr3dr4,1)-lbound(vect4dPtr3dr4,1)+1/=3 ) then + rc=991 + print *,'ERROR, 3D the vector dimension /= 3, rc=',rc + exit + endif + ! Get the _FillValue from the field attribute if exists + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", typekind=attTypeKind, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,attTypeKind,isPresent=',trim(fieldName),attTypeKind,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + if ( mvispresent ) then + if (attTypeKind==ESMF_TYPEKIND_R4) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r4, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r4,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + else if (attTypeKind==ESMF_TYPEKIND_R8) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r8, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r8,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + endif + istart = lbound(vect4dPtr3dr4,2) + iend = ubound(vect4dPtr3dr4,2) + jstart = lbound(vect4dPtr3dr4,3) + jend = ubound(vect4dPtr3dr4,3) + kstart = lbound(vect4dPtr3dr4,4) + kend = ubound(vect4dPtr3dr4,4) +!$omp parallel do default(shared) private(i,j,k) + do k=kstart,kend + do j=jstart, jend + do i=istart, iend + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R4) vect4dPtr3dr4(:,i,j,k)=missing_value_r4 + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R8) vect4dPtr3dr4(:,i,j,k)=missing_value_r8 + enddo + enddo + enddo + endif !mvispresent + ! 2-d vector fields with 3-d arrays + else + call ESMF_FieldGet(fcstField(ifld), localDe=0, farrayPtr=vect3dPtr2dr4, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, file=__FILE__)) return ! bail out + if( ubound(vect3dPtr2dr4,1)-lbound(vect3dPtr2dr4,1)+1 /= 3 ) then + rc=991 + print *,'ERROR, 2D the vector dimension /= 3, rc=',rc + exit + endif + ! Get the _FillValue from the field attribute if exists + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", typekind=attTypeKind, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,attTypeKind,isPresent=',trim(fieldName),attTypeKind,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + if ( mvispresent ) then + if (attTypeKind==ESMF_TYPEKIND_R4) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r4, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r4,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + else if (attTypeKind==ESMF_TYPEKIND_R8) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r8, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r8,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + endif + istart = lbound(vect3dPtr2dr4,2) + iend = ubound(vect3dPtr2dr4,2) + jstart = lbound(vect3dPtr2dr4,3) + jend = ubound(vect3dPtr2dr4,3) +!$omp parallel do default(shared) private(i,j) + do j=jstart, jend + do i=istart, iend + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R4) vect3dPtr2dr4(:,i,j)=missing_value_r4 + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R8) vect3dPtr2dr4(:,i,j)=missing_value_r8 + enddo + enddo + endif !mvispresent + endif + endif +! For non-vector fields + else + ! Only work on ESMF_TYPEKIND_R4 fields for now + if ( typekind == ESMF_TYPEKIND_R4 ) then + ! 2-d fields + if(fieldDimCount == gridDimCount) then + call ESMF_FieldGet(fcstField(ifld),localDe=0, farrayPtr=var2dPtr2dr4, rc=rc) + ! Get the _FillValue from the field attribute if exists + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", typekind=attTypeKind, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,attTypeKind,isPresent=',trim(fieldName),attTypeKind,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + if ( mvispresent ) then + if (attTypeKind==ESMF_TYPEKIND_R4) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r4, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r4,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + else if (attTypeKind==ESMF_TYPEKIND_R8) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r8, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r8,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + endif + istart = lbound(var2dPtr2dr4,1) + iend = ubound(var2dPtr2dr4,1) + jstart = lbound(var2dPtr2dr4,2) + jend = ubound(var2dPtr2dr4,2) +!$omp parallel do default(shared) private(i,j) + do j=jstart, jend + do i=istart, iend + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R4) var2dPtr2dr4(i,j)=missing_value_r4 + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R8) var2dPtr2dr4(i,j)=missing_value_r8 + enddo + enddo + endif !mvispresent + ! 3-d fields + else if(fieldDimCount == gridDimCount+1) then + call ESMF_FieldGet(fcstField(ifld),localDe=0, farrayPtr=var3dPtr3dr4, rc=rc) + ! Get the _FillValue from the field attribute if exists + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", typekind=attTypeKind, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,attTypeKind,isPresent=',trim(fieldName),attTypeKind,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + if ( mvispresent ) then + if (attTypeKind==ESMF_TYPEKIND_R4) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r4, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r4,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + else if (attTypeKind==ESMF_TYPEKIND_R8) then + call ESMF_AttributeGet(fcstField(ifld), convention="NetCDF", purpose="FV3", & + name="_FillValue", value=missing_value_r8, isPresent=mvispresent, rc=rc) + !write(msg,*) 'fieldName,_FillValue,isPresent=',trim(fieldName),missing_value_r8,mvispresent + !call ESMF_LogWrite("call mask field: "//trim(msg),ESMF_LOGMSG_INFO,rc=RC) + endif + istart = lbound(var3dPtr3dr4,1) + iend = ubound(var3dPtr3dr4,1) + jstart = lbound(var3dPtr3dr4,2) + jend = ubound(var3dPtr3dr4,2) + kstart = lbound(var3dPtr3dr4,3) + kend = ubound(var3dPtr3dr4,3) +!$omp parallel do default(shared) private(i,j,k) + do k=kstart,kend + do j=jstart, jend + do i=istart, iend + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R4) var3dPtr3dr4(i,j,k)=missing_value_r4 + if (maskwrt(i,j)<1.0 .and. attTypeKind==ESMF_TYPEKIND_R8) var3dPtr3dr4(i,j,k)=missing_value_r8 + enddo + enddo + enddo + endif !mvispresent + endif + endif + endif + enddo +! + deallocate(fcstField) + rc = 0 + + end subroutine mask_fields +! +!----------------------------------------------------------------------- ! subroutine ESMFproto_FieldBundleWrite(fieldbundle, fileName, & From 9caa172aaf8757dadc2590d2599706e00225bf51 Mon Sep 17 00:00:00 2001 From: Dom Heinzeller Date: Wed, 27 Nov 2019 14:33:32 -0700 Subject: [PATCH 15/18] Update submodule pointer for ccpp/physics --- ccpp/physics | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ccpp/physics b/ccpp/physics index 2e0d55754..c11b6e8bf 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit 2e0d55754ed8763a95f2aa035b2cb755db710350 +Subproject commit c11b6e8bf60f3d23be761b9aa4665d6611b9d7e0 From ac71e234cbecf8e7b516e6fafc7872d1bf45b009 Mon Sep 17 00:00:00 2001 From: Dom Heinzeller Date: Fri, 29 Nov 2019 17:03:49 -0700 Subject: [PATCH 16/18] ccpp/CMakeLists.txt: add flag '-Wall' to DEBUG flags for GNU compiler --- ccpp/CMakeLists.txt | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/ccpp/CMakeLists.txt b/ccpp/CMakeLists.txt index 1ee3fce3d..b95b2bc4f 100644 --- a/ccpp/CMakeLists.txt +++ b/ccpp/CMakeLists.txt @@ -89,8 +89,10 @@ if (${CMAKE_Fortran_COMPILER_ID} MATCHES "GNU") set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} -fcray-pointer -ffree-line-length-none -fno-range-check") set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} -fbacktrace -cpp") if (${CMAKE_BUILD_TYPE} MATCHES "Debug") + set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wall") + set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall") set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} -fno-unsafe-math-optimizations -frounding-math -fsignaling-nans") - set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} -ffpe-trap=invalid,zero,overflow -fcheck=bounds -fbacktrace -fno-range-check") + set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} -ffpe-trap=invalid,zero,overflow -fcheck=bounds -fbacktrace -fno-range-check -Wall") elseif (${CMAKE_BUILD_TYPE} MATCHES "Bitforbit") set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}") endif (${CMAKE_BUILD_TYPE} MATCHES "Debug") From 2de4377fddeaf2c8d58202e8fef54274b313092c Mon Sep 17 00:00:00 2001 From: Dom Heinzeller Date: Mon, 2 Dec 2019 11:02:38 -0700 Subject: [PATCH 17/18] Update submodule pointers for ccpp-framework and ccpp-physics --- ccpp/framework | 2 +- ccpp/physics | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/ccpp/framework b/ccpp/framework index f7ae56dba..1d1027267 160000 --- a/ccpp/framework +++ b/ccpp/framework @@ -1 +1 @@ -Subproject commit f7ae56dbaaff833222c77a313d22f894bd4167ac +Subproject commit 1d102726766a550b2c2292f0e964280c952ceb2c diff --git a/ccpp/physics b/ccpp/physics index c11b6e8bf..308a19747 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit c11b6e8bf60f3d23be761b9aa4665d6611b9d7e0 +Subproject commit 308a1974745b673ff6095b528bd3e35aebe5448c From eae67862e08a1285d4e2029aa8312f6d6d3ec7d8 Mon Sep 17 00:00:00 2001 From: Dom Heinzeller Date: Tue, 3 Dec 2019 15:22:19 -0700 Subject: [PATCH 18/18] Update submodule pointers for atmos_cubed_sphere, ccpp/framework, ccpp/physics --- atmos_cubed_sphere | 2 +- ccpp/framework | 2 +- ccpp/physics | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/atmos_cubed_sphere b/atmos_cubed_sphere index 68576a61f..3a4dfd8c6 160000 --- a/atmos_cubed_sphere +++ b/atmos_cubed_sphere @@ -1 +1 @@ -Subproject commit 68576a61f2a57236931e5ef1b8bd73a4256b0a5f +Subproject commit 3a4dfd8c6c4ceb8cec06397f25cb229ecd98065b diff --git a/ccpp/framework b/ccpp/framework index 1d1027267..42596226f 160000 --- a/ccpp/framework +++ b/ccpp/framework @@ -1 +1 @@ -Subproject commit 1d102726766a550b2c2292f0e964280c952ceb2c +Subproject commit 42596226fcbbf4c151829ec99ff0649adc09a5df diff --git a/ccpp/physics b/ccpp/physics index 308a19747..f895fc0aa 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit 308a1974745b673ff6095b528bd3e35aebe5448c +Subproject commit f895fc0aa011808cece39164443627a8126cc318