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

First reconciliation PR from production/RRFS.v1 #226

Merged
merged 4 commits into from
Oct 18, 2024
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
14 changes: 8 additions & 6 deletions physics/CONV/Grell_Freitas/cu_gf_deep.F90
Original file line number Diff line number Diff line change
Expand Up @@ -425,9 +425,9 @@ subroutine cu_gf_deep_run( &
integer :: turn,pmin_lev(its:ite),start_level(its:ite),ktopkeep(its:ite)
real(kind=kind_phys), dimension (its:ite,kts:kte) :: dtempdz
integer, dimension (its:ite,kts:kte) :: k_inv_layers
real(kind=kind_phys), dimension (its:ite) :: c0 ! HCB
real(kind=kind_phys), dimension (its:ite) :: c0, rrfs_factor ! HCB
real(kind=kind_phys), dimension (its:ite,kts:kte) :: c0t3d ! hli for smoke/dust wet scavenging
!$acc declare create(pmin_lev,start_level,ktopkeep,dtempdz,k_inv_layers,c0,c0t3d)
!$acc declare create(pmin_lev,start_level,ktopkeep,dtempdz,k_inv_layers,c0,rrfs_factor,c0t3d)

! rainevap from sas
real(kind=kind_phys) zuh2(40)
Expand Down Expand Up @@ -486,6 +486,7 @@ subroutine cu_gf_deep_run( &
! Set cloud water to rain water conversion rate (c0)
!$acc kernels
c0(:)=0.004
rrfs_factor(:)=1.
do i=its,itf
xland1(i)=int(xland(i)+.0001) ! 1.
if(xland(i).gt.1.5 .or. xland(i).lt.0.5)then
Expand All @@ -495,6 +496,7 @@ subroutine cu_gf_deep_run( &
if(imid.eq.1)then
c0(i)=0.002
endif
if(kdt.le.(4500./dtime))rrfs_factor(i)=1.-(float(kdt)/(4500./dtime)-1.)**2
enddo
!$acc end kernels

Expand Down Expand Up @@ -591,7 +593,6 @@ subroutine cu_gf_deep_run( &
sig(i)=(1.-frh)**2
!frh_out(i) = frh
if(forcing(i,7).eq.0.)sig(i)=1.
if(kdt.le.(3600./dtime))sig(i)=1.
frh_out(i) = frh*sig(i)
enddo
!$acc end kernels
Expand Down Expand Up @@ -2029,7 +2030,7 @@ subroutine cu_gf_deep_run( &
zuo,pre,pwo_ens,xmb,ktop, &
edto,pwdo,'deep',ierr2,ierr3, &
po_cup,pr_ens,maxens3, &
sig,closure_n,xland1,xmbm_in,xmbs_in, &
sig,closure_n,xland1,xmbm_in,xmbs_in,rrfs_factor, &
ichoice,imid,ipr,itf,ktf, &
its,ite, kts,kte, &
dicycle,xf_dicycle )
Expand Down Expand Up @@ -4056,7 +4057,7 @@ subroutine cup_output_ens_3d(xff_mid,xf_ens,ierr,dellat,dellaq,dellaqc, &
zu,pre,pw,xmb,ktop, &
edt,pwd,name,ierr2,ierr3,p_cup,pr_ens, &
maxens3, &
sig,closure_n,xland1,xmbm_in,xmbs_in, &
sig,closure_n,xland1,xmbm_in,xmbs_in,rrfs_factor, &
ichoice,imid,ipr,itf,ktf, &
its,ite, kts,kte, &
dicycle,xf_dicycle )
Expand Down Expand Up @@ -4118,7 +4119,7 @@ subroutine cup_output_ens_3d(xff_mid,xf_ens,ierr,dellat,dellaq,dellaqc, &
,intent (inout) :: &
ierr,ierr2,ierr3
integer, intent(in) :: dicycle
real(kind=kind_phys), intent(in), dimension (its:ite) :: xf_dicycle
real(kind=kind_phys), intent(in), dimension (its:ite) :: xf_dicycle, rrfs_factor
!$acc declare copyin(zu,pwd,p_cup,sig,xmbm_in,xmbs_in,edt,xff_mid,dellat,dellaqc,dellaq,pw,ktop,xland1,xf_dicycle)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe irrelevant, but if xf_dicycle is part to the $acc directive, should rrfs_factor be added?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@dustinswales If this caused a problem, do you think that this type of thing would be caught in the SCM nvhpc CI test? Just curious. I.e. is this something that the compiler would catch or would it just be a run-time issue (failure or performance hit)?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have no clue if this would get caught at buildtime/runtime.

!$acc declare copy(xf_ens,pr_ens,outtem,outq,outqc,pre,xmb,closure_n,ierr,ierr2,ierr3)
!
Expand Down Expand Up @@ -4198,6 +4199,7 @@ subroutine cup_output_ens_3d(xff_mid,xf_ens,ierr,dellat,dellaq,dellaqc, &
clos_wei=16./max(1.,closure_n(i))
xmb_ave(i)=min(xmb_ave(i),100.)
xmb(i)=clos_wei*sig(i)*xmb_ave(i)
if(dx(i)<dx_thresh) xmb(i)=rrfs_factor(i)*xmb(i)

if(xmb(i) < 1.e-16)then
ierr(i)=19
Expand Down
11 changes: 7 additions & 4 deletions physics/CONV/Grell_Freitas/cu_gf_driver.F90
Original file line number Diff line number Diff line change
Expand Up @@ -883,6 +883,13 @@ subroutine cu_gf_driver_run(ntracer,garea,im,km,dt,flag_init,flag_restart,&
cutenm(i)=0.
endif ! pret > 0

maxupmf(i)=0.
if(forcing2(i,6).gt.0.)then
maxupmf(i)=maxval(xmb(i)*zu(i,kts:ktf)/forcing2(i,6))
endif
if (xland(i)==0)then ! cu precip rate (mm/h)
if((maxupmf(i).lt.0.1) .or. (pret(i)*3600.lt.0.05)) pret(i)=0.
endif
if(pret(i).gt.0.)then
cuten(i)=1.
cutenm(i)=0.
Expand Down Expand Up @@ -999,10 +1006,6 @@ subroutine cu_gf_driver_run(ntracer,garea,im,km,dt,flag_init,flag_restart,&
gdc(i,15,10)=qfx(i)
gdc(i,16,10)=pret(i)*3600.

maxupmf(i)=0.
if(forcing2(i,6).gt.0.)then
maxupmf(i)=maxval(xmb(i)*zu(i,kts:ktf)/forcing2(i,6))
endif

if(ktop(i).gt.2 .and.pret(i).gt.0.)dt_mf(i,ktop(i)-1)=ud_mf(i,ktop(i))
endif
Expand Down
8 changes: 4 additions & 4 deletions physics/PBL/MYNN_EDMF/module_bl_mynn.F90
Original file line number Diff line number Diff line change
Expand Up @@ -2005,7 +2005,7 @@ SUBROUTINE mym_length ( &
ugrid = sqrt(u1(kts)**2 + v1(kts)**2)
uonset= 15.
wt_u = (1.0 - min(max(ugrid - uonset, 0.0)/30.0, 0.5))
cns = 2.7 !was 3.5
cns = 3.5
alp1 = 0.23
alp2 = 0.3
alp3 = 2.5 * wt_u !taper off bouyancy enhancement in shear-driven pbls
Expand Down Expand Up @@ -2039,7 +2039,7 @@ SUBROUTINE mym_length ( &
zwk = zw(k)
DO WHILE (zwk .LE. zi2+h1)
dzk = 0.5*( dz(k)+dz(k-1) )
qdz = min(max( qkw(k)-qmin, 0.03 ), 30.0)*dzk
qdz = min(max( qkw(k)-qmin, 0.02 ), 30.0)*dzk
elt = elt +qdz*zwk
vsc = vsc +qdz
k = k+1
Expand Down Expand Up @@ -5036,7 +5036,7 @@ SUBROUTINE mynn_tendencies(kts,kte,i, &
IF (FLAG_QI) THEN
DO k=kts,kte
Dth(k)=(thl(k) + xlvcp/exner(k)*sqc2(k) &
& + xlscp/exner(k)*(sqi2(k)+sqs(k)) &
& + xlscp/exner(k)*(sqi2(k)) & !+sqs(k)) &
& - th(k))/delt
!Use form from Tripoli and Cotton (1981) with their
!suggested min temperature to improve accuracy:
Expand Down Expand Up @@ -6057,7 +6057,7 @@ SUBROUTINE DMP_mf( &
if ((landsea-1.5).LT.0) then !land
acfac = .5*tanh((fltv2 - 0.02)/0.05) + .5
else !water
acfac = .5*tanh((fltv2 - 0.01)/0.03) + .5
acfac = .5*tanh((fltv2 - 0.015)/0.04) + .5
endif
!add a windspeed-dependent adjustment to acfac that tapers off
!the mass-flux scheme linearly above sfc wind speeds of 10 m/s.
Expand Down
3 changes: 1 addition & 2 deletions physics/SFC_Models/Land/RUC/lsm_ruc.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1312,8 +1312,7 @@ subroutine lsm_ruc_run & ! inputs

! --- ... accumulated total runoff and surface runoff
runoff(i) = runoff(i) + (drain(i)+runof(i)) * delt ! accum total kg m-2
!srunoff(i) = srunoff(i) + runof(i) * delt ! accum surface kg m-2
srunoff(i) = acrunoff(i,j) ! accum surface kg m-2
srunoff(i) = srunoff(i) + runof(i) * delt ! accum surface kg m-2

! --- ... accumulated frozen precipitation (accumulation in lsmruc)
snowfallac_lnd(i) = snfallac_lnd(i,j) ! accum kg m-2
Expand Down
4 changes: 2 additions & 2 deletions physics/SFC_Models/Land/RUC/module_sf_ruclsm.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1740,7 +1740,7 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
!-- will reduce warm bias in western Canada
!-- and US West coast, where max snow albedo is low (0.3-0.5).
!print *,'ALB increase to 0.7',alb_snow,snhei,snhei_crit,albsn,i,j
!ALBsn = 0.7_kind_phys
ALBsn = 0.7_kind_phys
endif

Emiss= emissn
Expand All @@ -1753,7 +1753,7 @@ SUBROUTINE SFCTMP (debug_print, delt,ktau,conflx,i,j, & !--- input varia
!-- will reduce warm bias in western Canada
!-- and US West coast, where max snow albedo is low (0.3-0.5).
!print *,'ALB increase to 0.7',alb_snow,snhei,snhei_crit,albsn,i,j
!ALBsn = 0.7_kind_phys
ALBsn = 0.7_kind_phys
!print *,'NO mosaic ALB increase to 0.7',alb_snow,snhei,snhei_crit,alb,i,j
endif

Expand Down
42 changes: 25 additions & 17 deletions physics/smoke_dust/module_smoke_plumerise.F90
Original file line number Diff line number Diff line change
Expand Up @@ -109,12 +109,6 @@ subroutine plumerise(m1,m2,m3,ia,iz,ja,jz, &
! print *,' Plumerise_scalar 1',ncall
coms => get_thread_coms()

IF (frp_inst<frp_threshold) THEN
k1=1
k2=2
!return
END IF

! print *,' Plumerise_scalar 2',m1
j=1
i=1
Expand Down Expand Up @@ -169,12 +163,20 @@ subroutine plumerise(m1,m2,m3,ia,iz,ja,jz, &
WRITE(1000+mpiid,*) 'inside plumerise: xlat,xlong,curr_secs,imm,FRP,burnt_area ', lat, long, int(curr_secs), imm, FRP,burnt_area
END IF

IF (frp_inst<frp_threshold) THEN
k1=1
k2=2
!exit
return
END IF

!- get fire properties (burned area, plume radius, heating rates ...)
call get_fire_properties(coms,imm,iveg_ag,burnt_area,FRP,errmsg,errflg)
if(errflg/=0) return


!------ generates the plume rise ------
call makeplume (coms,kmt,ztopmax(imm),ixx,imm)
call makeplume (coms,kmt,ztopmax(imm),ixx,imm,mpiid)

IF ( dbg_opt .and. (icall .le. n_dbg_lines) .and. (frp_inst .ge. frp_threshold) ) then
WRITE(1000+mpiid,*) 'inside plumerise after makeplume:xlat,xlong,curr_secs,imm,kmt,ztopmax(imm) ', lat, long, int(curr_secs), imm,kmt, ztopmax(imm)
Expand Down Expand Up @@ -562,7 +564,7 @@ subroutine get_fire_properties(coms,imm,iveg_ag,burnt_area,FRP,errmsg,errflg)
end subroutine get_fire_properties
!-------------------------------------------------------------------------------
!
SUBROUTINE MAKEPLUME (coms,kmt,ztopmax,ixx,imm)
SUBROUTINE MAKEPLUME (coms,kmt,ztopmax,ixx,imm,mpiid)
!
! *********************************************************************
!
Expand Down Expand Up @@ -621,22 +623,23 @@ SUBROUTINE MAKEPLUME (coms,kmt,ztopmax,ixx,imm)
!
!
!**********************************************************************
!**********************************************************************
!use module_zero_plumegen_coms
implicit none
!logical :: endspace
!**********************************************************************
!use module_zero_plumegen_coms
implicit none
!logical :: endspace
type(plumegen_coms), pointer :: coms
character (len=10) :: varn
integer :: izprint, iconv, itime, k, kk, kkmax, deltak,ilastprint,kmt &
,ixx,nrectotal,i_micro,n_sub_step
real(kind=kind_phys) :: vc, g, r, cp, eps, &
tmelt, heatsubl, heatfus, heatcond, tfreeze, &
ztopmax, wmax, rmaxtime, es, esat, heat,dt_save !ESAT_PR,
character (len=2) :: cixx
character (len=2) :: cixx
integer, intent(in) :: mpiid
! Set threshold to be the same as dz=100., the constant grid spacing of plume grid model(meters) found in set_grid()
REAL(kind=kind_phys) :: DELZ_THRESOLD = 100.

INTEGER :: imm
INTEGER :: imm, dtknt

! real(kind=kind_phys), external:: esat_pr!
!
Expand All @@ -654,6 +657,7 @@ SUBROUTINE MAKEPLUME (coms,kmt,ztopmax,ixx,imm)
coms%viscosity = 500.!- coms%viscosity constant (original value: 0.001)

nrectotal=150
dtknt = 0
!
!*************** PROBLEM SETUP AND INITIAL CONDITIONS *****************
coms%mintime = 1
Expand Down Expand Up @@ -697,9 +701,13 @@ SUBROUTINE MAKEPLUME (coms,kmt,ztopmax,ixx,imm)
!sam 81 format('nm1=',I0,' from kmt=',I0,' kkmax=',I0,' deltak=',I0)
!sam write(0,81) coms%nm1,kmt,kkmax,deltak
!-- set timestep
!coms%dt = (coms%zm(2)-coms%zm(1)) / (coms%tstpf * wmax)
coms%dt = min(5.,(coms%zm(2)-coms%zm(1)) / (coms%tstpf * wmax))

!coms%dt = (coms%zm(2)-coms%zm(1)) / (coms%tstpf * wmax) i
coms%dt = max(0.01,min(5.,(coms%zm(2)-coms%zm(1)) / (coms%tstpf * wmax)))
dtknt = dtknt + 1
! if (coms%dt .ne. 5.)then
! WRITE(1000+mpiid,*) 'dtknt,zm(2),zm(1) ', dtknt,coms%zm(2),coms%zm(1)
! WRITE(1000+mpiid,*) 'coms%tstpf,wmax,dt =', coms%tstpf,wmax,coms%dt
! endif
!-- elapsed time, sec
coms%time = coms%time+coms%dt
!-- elapsed time, minutes
Expand Down
7 changes: 5 additions & 2 deletions physics/smoke_dust/rrfs_smoke_wrapper.F90
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,7 @@ subroutine rrfs_smoke_wrapper_run(im, flag_init, kte, kme, ktau, dt, garea, land
real(kind_phys), dimension(:,:), intent(in), optional :: emi_ant_in
real(kind_phys), dimension(:), intent(in) :: u10m, v10m, ustar, dswsfc, &
recmol, garea, rlat,rlon, tskin, pb2d, zorl, snow, &
rain_cpl, rainc_cpl, hf2d, t2m, dpt2m
rain_cpl, rainc_cpl, hf2d, t2m, dpt2m
real(kind_phys), dimension(:,:), intent(in) :: vegtype_frac
real(kind_phys), dimension(:,:), intent(in) :: ph3d, pr3d
real(kind_phys), dimension(:,:), intent(in) :: phl3d, prl3d, tk3d, us3d, vs3d, spechum, w
Expand Down Expand Up @@ -347,7 +347,7 @@ subroutine rrfs_smoke_wrapper_run(im, flag_init, kte, kme, ktau, dt, garea, land
do k=kts,kte
do i=its,ite
! ebu is divided by coef_bb_dc since it is applied in the output
ebu(i,k,1)=ebu_smoke(i,k) / coef_bb_dc(i,1)
ebu(i,k,1)=ebu_smoke(i,k) / MAX(1.E-4,coef_bb_dc(i,1))
enddo
enddo
ENDIF
Expand Down Expand Up @@ -752,6 +752,9 @@ subroutine rrfs_smoke_prep( &
moist = 0._kind_phys
chem = 0._kind_phys
z_at_w = 0._kind_phys
if ( ebb_dcycle == 1 ) then
coef_bb_dc = 1._kind_phys
endif

do i=its,ite
u10 (i,1)=u10m (i)
Expand Down
Loading