Skip to content
New issue

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

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

Already on GitHub? Sign in to your account

Update gsl/develop from develop 2020/12/08 #61

Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 5 additions & 3 deletions atmos_model.F90
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ module atmos_model_mod
use IPD_driver, only: IPD_initialize, IPD_initialize_rst
use CCPP_driver, only: CCPP_step, non_uniform_blocks

use stochastic_physics_wrapper_mod, only: stochastic_physics_wrapper
use stochastic_physics_wrapper_mod, only: stochastic_physics_wrapper,stochastic_physics_wrapper_end
#else
use IPD_driver, only: IPD_initialize, IPD_initialize_rst, IPD_step
use physics_abstraction_layer, only: time_vary_step, radiation_step1, physics_step1, physics_step2
Expand Down Expand Up @@ -962,6 +962,9 @@ subroutine atmos_model_end (Atmos)
!---- termination routine for atmospheric model ----

call atmosphere_end (Atmos % Time, Atmos%grid, restart_endfcst)

call stochastic_physics_wrapper_end(IPD_Control)

if(restart_endfcst) then
call FV3GFS_restart_write (IPD_Data, IPD_Restart, Atm_block, &
IPD_Control, Atmos%domain)
Expand Down Expand Up @@ -1761,7 +1764,6 @@ subroutine assign_importdata(rc)
nb = Atm_block%blkno(i,j)
ix = Atm_block%ixp(i,j)

IPD_Data(nb)%Sfcprop%fice(ix) = zero
IPD_Data(nb)%Coupling%slimskin_cpl(ix) = IPD_Data(nb)%Sfcprop%slmsk(ix)
ofrac = IPD_Data(nb)%Sfcprop%oceanfrac(ix)
if (ofrac > zero) then
Expand All @@ -1776,7 +1778,7 @@ subroutine assign_importdata(rc)
if (abs(one-ofrac) < epsln) then
IPD_Data(nb)%Sfcprop%slmsk(ix) = zero
IPD_Data(nb)%Coupling%slimskin_cpl(ix) = zero
end if
endif
endif
endif
enddo
Expand Down
2 changes: 1 addition & 1 deletion ccpp/physics
45 changes: 27 additions & 18 deletions io/FV3GFS_io.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1070,14 +1070,7 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain)
endif

if(Model%frac_grid) then ! obtain slmsk from landfrac
!! next 5 lines are temporary till lake model is available
if (Sfcprop(nb)%lakefrac(ix) > zero) then
! Sfcprop(nb)%lakefrac(ix) = nint(Sfcprop(nb)%lakefrac(ix))
Sfcprop(nb)%landfrac(ix) = one - Sfcprop(nb)%lakefrac(ix)
if (Sfcprop(nb)%lakefrac(ix) == zero) Sfcprop(nb)%fice(ix) = zero
endif
Sfcprop(nb)%slmsk(ix) = ceiling(Sfcprop(nb)%landfrac(ix))
if (Sfcprop(nb)%fice(ix) > Model%min_lakeice .and. Sfcprop(nb)%landfrac(ix) == zero) Sfcprop(nb)%slmsk(ix) = 2 ! land dominates ice if co-exist
Sfcprop(nb)%slmsk(ix) = ceiling(Sfcprop(nb)%landfrac(ix)) !nint/floor are options
else ! obtain landfrac from slmsk
if (Sfcprop(nb)%slmsk(ix) > 1.9_r8) then
Sfcprop(nb)%landfrac(ix) = zero
Expand All @@ -1088,16 +1081,32 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain)

if (Sfcprop(nb)%lakefrac(ix) > zero) then
Sfcprop(nb)%oceanfrac(ix) = zero ! lake & ocean don't coexist in a cell
! if (Sfcprop(nb)%fice(ix) < Model%min_lakeice) then
! Sfcprop(nb)%fice(ix) = zero
! if (Sfcprop(nb)%slmsk(ix) == 2) Sfcprop(nb)%slmsk(ix) = 0
! endif
if (Sfcprop(nb)%slmsk(ix) /= one) then
if (Sfcprop(nb)%fice(ix) >= Model%min_lakeice) then
if (Sfcprop(nb)%slmsk(ix) < 1.9_r8) &
write(*,'(a,2i3,3f6.2)') 'reset lake slmsk=2 at nb,ix=' &
,nb,ix,Sfcprop(nb)%fice(ix),Sfcprop(nb)%slmsk(ix),Sfcprop(nb)%lakefrac(ix)
Sfcprop(nb)%slmsk(ix) = 2.
else if (Sfcprop(nb)%slmsk(ix) > 1.e-7) then
write(*,'(a,2i3,3f6.2)') 'reset lake slmsk=0 at nb,ix=' &
,nb,ix,Sfcprop(nb)%fice(ix),Sfcprop(nb)%slmsk(ix),Sfcprop(nb)%lakefrac(ix)
Sfcprop(nb)%slmsk(ix) = zero
end if
end if
else
Sfcprop(nb)%oceanfrac(ix) = one - Sfcprop(nb)%landfrac(ix)
! if (Sfcprop(nb)%fice(ix) < Model%min_seaice) then
! Sfcprop(nb)%fice(ix) = zero
! if (Sfcprop(nb)%slmsk(ix) == 2) Sfcprop(nb)%slmsk(ix) = 0
! endif
if (Sfcprop(nb)%slmsk(ix) /= one) then
if (Sfcprop(nb)%fice(ix) >= Model%min_seaice) then
if (Sfcprop(nb)%slmsk(ix) < 1.9_r8) &
write(*,'(a,2i3,3f6.2)') 'reset sea slmsk=2 at nb,ix=' &
,nb,ix,Sfcprop(nb)%fice(ix),Sfcprop(nb)%slmsk(ix),Sfcprop(nb)%landfrac(ix)
Sfcprop(nb)%slmsk(ix) = 2.
else if (Sfcprop(nb)%slmsk(ix) > 1.e-7) then
write(*,'(a,2i3,4f6.2)') 'reset sea slmsk=0 at nb,ix=' &
,nb,ix,Sfcprop(nb)%fice(ix),Sfcprop(nb)%slmsk(ix),Sfcprop(nb)%landfrac(ix)
Sfcprop(nb)%slmsk(ix) = zero
end if
end if
endif
!
!--- NSSTM variables
Expand Down Expand Up @@ -1351,7 +1360,7 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain)
endif

if (sfc_var2(i,j,nvar_s2m) < -9990.0_r8) then
if (Model%me == Model%master ) call mpp_error(NOTE, 'gfs_driver::surface_props_input - computing zorli')
if (Model%me == Model%master ) call mpp_error(NOTE, 'gfs_driver::surface_props_input - computing zorlw')
!$omp parallel do default(shared) private(nb, ix)
do nb = 1, Atm_block%nblks
do ix = 1, Atm_block%blksz(nb)
Expand All @@ -1366,7 +1375,7 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain)
!$omp parallel do default(shared) private(nb, ix, tem, tem1)
do nb = 1, Atm_block%nblks
do ix = 1, Atm_block%blksz(nb)
Sfcprop(nb)%tsfco(ix) = max(con_tice, Sfcprop(nb)%tsfco(ix))
if( Model%phour < 1.e-7) Sfcprop(nb)%tsfco(ix) = max(con_tice, Sfcprop(nb)%tsfco(ix)) ! this may break restart reproducibility
tem1 = one - Sfcprop(nb)%landfrac(ix)
tem = tem1 * Sfcprop(nb)%fice(ix) ! tem = ice fraction wrt whole cell
Sfcprop(nb)%zorl(ix) = Sfcprop(nb)%zorll(ix) * Sfcprop(nb)%landfrac(ix) &
Expand Down
111 changes: 64 additions & 47 deletions stochastic_physics/stochastic_physics_wrapper.F90
Original file line number Diff line number Diff line change
Expand Up @@ -92,25 +92,43 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr)
return
endif
end if
allocate(xlat(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(xlon(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
if (GFS_Control%do_sppt) then
allocate(sppt_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if (GFS_Control%do_shum) then
allocate(shum_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if (GFS_Control%do_skeb) then
allocate(skebu_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
allocate(skebv_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if ( GFS_Control%lndp_type .EQ. 2 ) then ! this scheme updates through forecast
allocate(sfc_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%n_var_lndp))
end if
if (GFS_Control%lndp_type .EQ. 2) then ! save wts, and apply lndp scheme
allocate(smc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(slc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(stc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(stype(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(vfrac(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
endif

do nb=1,Atm_block%nblks
xlat(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlat(:)
xlon(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlon(:)
end do

if ( GFS_Control%lndp_type .EQ. 1 ) then ! this scheme sets perts once
! Copy blocked data into contiguous arrays; no need to copy sfc_wts in (intent out)
allocate(xlat(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(xlon(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(sfc_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%n_var_lndp))
do nb=1,Atm_block%nblks
xlat(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlat(:)
xlon(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlon(:)
end do
call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%phour, GFS_Control%blksz, xlat=xlat, xlon=xlon, &
sppt_wts=sppt_wts, shum_wts=shum_wts, skebu_wts=skebu_wts, skebv_wts=skebv_wts, sfc_wts=sfc_wts, &
nthreads=nthreads)
! Copy contiguous data back; no need to copy xlat/xlon, these are intent(in) - just deallocate
do nb=1,Atm_block%nblks
GFS_Data(nb)%Coupling%sfc_wts(:,:) = sfc_wts(nb,1:GFS_Control%blksz(nb),:)
end do
deallocate(xlat)
deallocate(xlon)
deallocate(sfc_wts)
end if
! Consistency check for cellular automata
Expand All @@ -126,27 +144,6 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr)
else initalize_stochastic_physics

if (GFS_Control%do_sppt .OR. GFS_Control%do_shum .OR. GFS_Control%do_skeb .OR. (GFS_Control%lndp_type .EQ. 2) ) then
! Copy blocked data into contiguous arrays; no need to copy weights in (intent(out))
allocate(xlat(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(xlon(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
do nb=1,Atm_block%nblks
xlat(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlat(:)
xlon(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Grid%xlon(:)
end do
if (GFS_Control%do_sppt) then
allocate(sppt_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if (GFS_Control%do_shum) then
allocate(shum_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if (GFS_Control%do_skeb) then
allocate(skebu_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
allocate(skebv_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%levs))
end if
if ( GFS_Control%lndp_type .EQ. 2 ) then ! this scheme updates through forecast
allocate(sfc_wts(1:Atm_block%nblks,maxval(GFS_Control%blksz),1:GFS_Control%n_var_lndp))
end if

call run_stochastic_physics(GFS_Control%levs, GFS_Control%kdt, GFS_Control%phour, GFS_Control%blksz, xlat=xlat, xlon=xlon, &
sppt_wts=sppt_wts, shum_wts=shum_wts, skebu_wts=skebu_wts, skebv_wts=skebv_wts, sfc_wts=sfc_wts, &
nthreads=nthreads)
Expand All @@ -155,32 +152,23 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr)
do nb=1,Atm_block%nblks
GFS_Data(nb)%Coupling%sppt_wts(:,:) = sppt_wts(nb,1:GFS_Control%blksz(nb),:)
end do
deallocate(sppt_wts)
end if
if (GFS_Control%do_shum) then
do nb=1,Atm_block%nblks
GFS_Data(nb)%Coupling%shum_wts(:,:) = shum_wts(nb,1:GFS_Control%blksz(nb),:)
end do
deallocate(shum_wts)
end if
if (GFS_Control%do_skeb) then
do nb=1,Atm_block%nblks
GFS_Data(nb)%Coupling%skebu_wts(:,:) = skebu_wts(nb,1:GFS_Control%blksz(nb),:)
GFS_Data(nb)%Coupling%skebv_wts(:,:) = skebv_wts(nb,1:GFS_Control%blksz(nb),:)
end do
deallocate(skebu_wts)
deallocate(skebv_wts)
end if
if (GFS_Control%lndp_type .EQ. 2) then ! save wts, and apply lndp scheme
do nb=1,Atm_block%nblks
GFS_Data(nb)%Coupling%sfc_wts(:,:) = sfc_wts(nb,1:GFS_Control%blksz(nb),:)
end do

allocate(smc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(slc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(stc(1:Atm_block%nblks,maxval(GFS_Control%blksz),GFS_Control%lsoil))
allocate(stype(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
allocate(vfrac(1:Atm_block%nblks,maxval(GFS_Control%blksz)))
do nb=1,Atm_block%nblks
stype(nb,1:GFS_Control%blksz(nb)) = GFS_Data(nb)%Sfcprop%stype(:)
smc(nb,1:GFS_Control%blksz(nb),:) = GFS_Data(nb)%Sfcprop%smc(:,:)
Expand All @@ -202,21 +190,13 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr)
write(6,*) 'call to GFS_apply_lndp failed'
return
endif
deallocate(stype)
deallocate(sfc_wts)
do nb=1,Atm_block%nblks
GFS_Data(nb)%Sfcprop%smc(:,:) = smc(nb,1:GFS_Control%blksz(nb),:)
GFS_Data(nb)%Sfcprop%slc(:,:) = slc(nb,1:GFS_Control%blksz(nb),:)
GFS_Data(nb)%Sfcprop%stc(:,:) = stc(nb,1:GFS_Control%blksz(nb),:)
GFS_Data(nb)%Sfcprop%vfrac(:) = vfrac(nb,1:GFS_Control%blksz(nb))
enddo
deallocate(smc)
deallocate(slc)
deallocate(stc)
deallocate(vfrac)
endif ! lndp block
deallocate(xlat)
deallocate(xlon)
end if

endif initalize_stochastic_physics
Expand Down Expand Up @@ -309,4 +289,41 @@ subroutine stochastic_physics_wrapper (GFS_Control, GFS_Data, Atm_block, ierr)

end subroutine stochastic_physics_wrapper


subroutine stochastic_physics_wrapper_end (GFS_Control)

use GFS_typedefs, only: GFS_control_type, GFS_data_type
use stochastic_physics, only: finalize_stochastic_physics

implicit none

type(GFS_control_type), intent(inout) :: GFS_Control

if (GFS_Control%do_sppt .OR. GFS_Control%do_shum .OR. GFS_Control%do_skeb .OR. (GFS_Control%lndp_type .GT. 0) ) then
if (allocated(xlat)) deallocate(xlat)
if (allocated(xlon)) deallocate(xlon)
if (GFS_Control%do_sppt) then
if (allocated(sppt_wts)) deallocate(sppt_wts)
end if
if (GFS_Control%do_shum) then
if (allocated(shum_wts)) deallocate(shum_wts)
end if
if (GFS_Control%do_skeb) then
if (allocated(skebu_wts)) deallocate(skebu_wts)
if (allocated(skebv_wts)) deallocate(skebv_wts)
end if
if ( GFS_Control%lndp_type .EQ. 2 ) then ! this scheme updates through forecast
if (allocated(sfc_wts)) deallocate(sfc_wts)
end if
if (GFS_Control%lndp_type .EQ. 2) then ! save wts, and apply lndp scheme
if (allocated(smc)) deallocate(smc)
if (allocated(slc)) deallocate(slc)
if (allocated(stc)) deallocate(stc)
if (allocated(stype)) deallocate(stype)
if (allocated(vfrac)) deallocate(vfrac)
endif
call finalize_stochastic_physics()
endif
end subroutine stochastic_physics_wrapper_end

end module stochastic_physics_wrapper_mod