Skip to content

Commit 6019b52

Browse files
authored
Merge branch 'develop' into feature/zhaobin74/update-cice6-forcing-nml
2 parents 2c2fd04 + c0f5822 commit 6019b52

File tree

3 files changed

+134
-35
lines changed

3 files changed

+134
-35
lines changed

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/ConvPar_GF2020.F90

Lines changed: 9 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -3514,9 +3514,9 @@ SUBROUTINE CUP_gf(its,ite,kts,kte ,itf,ktf, mtp, nmp &
35143514
zdo(i,k )*(hcdo(i,k )-heo_cup(i,k ) ) )*g/dp*edto(i)
35153515

35163516
!---meltglac-------------------------------------------------
3517-
dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))*0.5*(qrco(i,k+1)+qrco(i,k)) &
3518-
- melting(i,k))*g/dp
3519-
3517+
!dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))*0.5*(qrco(i,k+1)+qrco(i,k)) &
3518+
! - melting(i,k))*g/dp
3519+
dellah(i,k) = dellah(i,k) - xlf*melting(i,k)*g/dp
35203520
!-- for output only
35213521
subten_H(i,k) = -(zuo(i,k+1)*(-heo_cup(i,k+1)) - zuo(i,k)*(-heo_cup(i,k)))*g/dp &
35223522
+(zdo(i,k+1)*(-heo_cup(i,k+1)) - zdo(i,k)*(-heo_cup(i,k)))*g/dp*edto(i)
@@ -3678,9 +3678,9 @@ SUBROUTINE CUP_gf(its,ite,kts,kte ,itf,ktf, mtp, nmp &
36783678
+(zdo(i,k+1)*(hcdo(i,k+1)-heo_cup(i,k+1) ) - &
36793679
zdo(i,k )*(hcdo(i,k )-heo_cup(i,k ) ) )*g/dp*edto(i)
36803680

3681-
dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))* &
3682-
0.5*(qrco(i,k+1)+qrco(i,k)) - melting(i,k))*g/dp
3683-
3681+
!dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))* &
3682+
! 0.5*(qrco(i,k+1)+qrco(i,k)) - melting(i,k))*g/dp
3683+
dellah(i,k) = dellah(i,k) - xlf*melting(i,k)*g/dp
36843684
!--- for output only
36853685
subten_H(i,k) = -(zuo(i,k+1)*(-heo_cup(i,k+1)) - zuo(i,k)*(-heo_cup(i,k)))*g/dp &
36863686
+(zdo(i,k+1)*(-heo_cup(i,k+1)) - zdo(i,k)*(-heo_cup(i,k)))*g/dp*edto(i)
@@ -3707,8 +3707,9 @@ SUBROUTINE CUP_gf(its,ite,kts,kte ,itf,ktf, mtp, nmp &
37073707
dellah(i,k) =-( zuo(i,k+1)*hco (i,k+1) - zuo(i,k)*hco (i,k) )*g/dp &
37083708
+( zdo(i,k+1)*hcdo(i,k+1) - zdo(i,k)*hcdo(i,k) )*g/dp*edto(i)
37093709

3710-
dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))* &
3711-
0.5*(qrco(i,k+1)+qrco(i,k)) - melting(i,k))*g/dp
3710+
!dellah(i,k) = dellah(i,k) + xlf*((1.-p_liq_ice(i,k))* &
3711+
! 0.5*(qrco(i,k+1)+qrco(i,k)) - melting(i,k))*g/dp
3712+
dellah(i,k) = dellah(i,k) - xlf*melting(i,k)*g/dp
37123713
!- update with subsidence term from the FCT scheme
37133714
dellah(i,k) = dellah(i,k) + sub_tend(1,k)
37143715
!--- for output only

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/Process_Library.F90

Lines changed: 23 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,13 @@ module GEOSmoist_Process_Library
3232
module procedure ICE_FRACTION_1D
3333
module procedure ICE_FRACTION_SC
3434
end interface ICE_FRACTION
35+
36+
! SRF_TYPE constants
37+
integer, parameter :: SRF_TYPE_LAND = 1
38+
integer, parameter :: SRF_TYPE_SNOW = 2
39+
integer, parameter :: SRF_TYPE_ICE = 3
40+
integer, parameter :: SRF_TYPE_OCEAN = 0
41+
3542
! ICE_FRACTION constants
3643
! In anvil/convective clouds
3744
real, parameter :: aT_ICE_ALL = 252.16
@@ -403,29 +410,20 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
403410
ICEFRCT_C = MAX(ICEFRCT_C,0.00)
404411
ICEFRCT_C = ICEFRCT_C**aICEFRPWR
405412
! Sigmoidal functions like figure 6b/6c of Hu et al 2010, doi:10.1029/2009JD012384
406-
if (SRF_TYPE >= 2.0) then
413+
select case (nint(SRF_TYPE))
414+
case (SRF_TYPE_SNOW, SRF_TYPE_ICE)
407415
! Over snow (SRF_TYPE == 2.0) and ice (SRF_TYPE == 3.0)
408-
if (ICE_RADII_PARAM == 1) then
409-
! Jason formula
410-
ICEFRCT_M = 0.00
411-
if ( TEMP <= JiT_ICE_ALL ) then
412-
ICEFRCT_M = 1.000
413-
else if ( (TEMP > JiT_ICE_ALL) .AND. (TEMP <= JiT_ICE_MAX) ) then
414-
ICEFRCT_M = 1.00 - ( TEMP - JiT_ICE_ALL ) / ( JiT_ICE_MAX - JiT_ICE_ALL )
415-
end if
416-
else
417-
ICEFRCT_M = 0.00
418-
if ( TEMP <= iT_ICE_ALL ) then
419-
ICEFRCT_M = 1.000
420-
else if ( (TEMP > iT_ICE_ALL) .AND. (TEMP <= iT_ICE_MAX) ) then
421-
ICEFRCT_M = SIN( 0.5*MAPL_PI*( 1.00 - ( TEMP - iT_ICE_ALL ) / ( iT_ICE_MAX - iT_ICE_ALL ) ) )
422-
end if
416+
ICEFRCT_M = 0.00
417+
if ( TEMP <= iT_ICE_ALL ) then
418+
ICEFRCT_M = 1.000
419+
else if ( (TEMP > iT_ICE_ALL) .AND. (TEMP <= iT_ICE_MAX) ) then
420+
ICEFRCT_M = SIN( 0.5*MAPL_PI*( 1.00 - ( TEMP - iT_ICE_ALL ) / ( iT_ICE_MAX - iT_ICE_ALL ) ) )
423421
end if
424422
ICEFRCT_M = MIN(ICEFRCT_M,1.00)
425423
ICEFRCT_M = MAX(ICEFRCT_M,0.00)
426424
ICEFRCT_M = ICEFRCT_M**iICEFRPWR
427-
else if (SRF_TYPE > 1.0) then
428-
! Over Land
425+
case (SRF_TYPE_LAND)
426+
! Over Land (SRF_TYPE == 1)
429427
ICEFRCT_M = 0.00
430428
if ( TEMP <= lT_ICE_ALL ) then
431429
ICEFRCT_M = 1.000
@@ -435,8 +433,8 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
435433
ICEFRCT_M = MIN(ICEFRCT_M,1.00)
436434
ICEFRCT_M = MAX(ICEFRCT_M,0.00)
437435
ICEFRCT_M = ICEFRCT_M**lICEFRPWR
438-
else
439-
! Over Oceans
436+
case (SRF_TYPE_OCEAN)
437+
! Over Oceans (SRF_TYPE == 0)
440438
ICEFRCT_M = 0.00
441439
if ( TEMP <= oT_ICE_ALL ) then
442440
ICEFRCT_M = 1.000
@@ -446,7 +444,11 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
446444
ICEFRCT_M = MIN(ICEFRCT_M,1.00)
447445
ICEFRCT_M = MAX(ICEFRCT_M,0.00)
448446
ICEFRCT_M = ICEFRCT_M**oICEFRPWR
449-
endif
447+
case default
448+
! You should not be here
449+
print *, 'ICE_FRACTION_SC: Unknown SRF_TYPE = ',SRF_TYPE
450+
error stop
451+
end select
450452
! Combine the Convective and MODIS functions
451453
ICEFRCT = ICEFRCT_M*(1.0-CNV_FRACTION) + ICEFRCT_C*(CNV_FRACTION)
452454
#endif

GEOSmkiau_GridComp/GEOS_mkiauGridComp.F90

Lines changed: 102 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -96,6 +96,7 @@ subroutine SetServices ( GC, RC )
9696
type (ESMF_Config) :: CF
9797

9898
logical :: BLEND_AT_PBL
99+
logical :: BLEND_QV_AT_TP
99100
#ifdef PYMKIAU_INTEGRATION
100101
! IEEE trapping see below
101102
logical :: halting_mode(5)
@@ -125,6 +126,11 @@ subroutine SetServices ( GC, RC )
125126
call MAPL_GetResource(MAPL, BLEND_AT_PBL, LABEL="REPLAY_BLEND_AT_PBL:", default=.FALSE., RC=status)
126127
VERIFY_(STATUS)
127128

129+
call MAPL_GetResource(MAPL, BLEND_QV_AT_TP, LABEL="REPLAY_BLEND_QV_AT_TP:", default=.FALSE., RC=status)
130+
VERIFY_(STATUS)
131+
132+
if ( BLEND_AT_PBL ) BLEND_QV_AT_TP = .FALSE.
133+
128134
! Set the Run entry points (phase 1 for regular IAU and phase 2 for clearing
129135
! --------------------------------------------------------------------------
130136

@@ -251,6 +257,17 @@ subroutine SetServices ( GC, RC )
251257
VERIFY_(STATUS)
252258
endif
253259

260+
if( BLEND_QV_AT_TP ) then
261+
call MAPL_AddImportSpec(GC, &
262+
SHORT_NAME = 'TROPP_BLENDED', &
263+
LONG_NAME = 'tropopause_pressure_based_on_blended_estimate', &
264+
UNITS = 'Pa', &
265+
DIMS = MAPL_DimsHorzOnly, &
266+
VLOCATION = MAPL_VLocationNone, &
267+
RC=STATUS )
268+
VERIFY_(STATUS)
269+
endif
270+
254271

255272
! !EXPORT STATE:
256273
! --------------
@@ -597,6 +614,7 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
597614
real, pointer, dimension(:,:,:) :: pdum1 => null()
598615
real, pointer, dimension(:,:,:) :: pdum2 => null()
599616
real, pointer, dimension(:,:) :: blnpp => null()
617+
real, pointer, dimension(:,:) :: tropp => null()
600618

601619
real, allocatable, dimension(:,:,:) :: du_fix
602620
real, allocatable, dimension(:,:,:) :: dv_fix
@@ -662,6 +680,7 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
662680
real :: FACP1, FACP0, FACM1, FACM2
663681
real :: DAMPBEG, DAMPEND
664682
logical :: BLEND_AT_PBL
683+
logical :: BLEND_QV_AT_TP
665684
integer :: i,j,L,n
666685
integer :: nt,nvars,natts
667686
integer :: nymd, nhms
@@ -894,6 +913,9 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC )
894913
call MAPL_GetResource(MAPL, BLEND_AT_PBL, LABEL="REPLAY_BLEND_AT_PBL:", default=.FALSE., RC=status)
895914
VERIFY_(STATUS)
896915

916+
call MAPL_GetResource(MAPL, BLEND_QV_AT_TP, LABEL="REPLAY_BLEND_QV_AT_TP:", default=.FALSE., RC=status)
917+
VERIFY_(STATUS)
918+
897919
CREMAP = ESMF_UtilStringUpperCase(CREMAP)
898920
FIXWIND = ESMF_UtilStringUpperCase(FIXWIND)
899921
DOWINDFIX = trim(FIXWIND)=="YES"
@@ -2145,7 +2167,7 @@ subroutine handleANA_
21452167
! **** with option to blend QV specially, starting at tropopause. ****
21462168
! **********************************************************************
21472169

2148-
if( DAMPBEG.ne.DAMPEND .or. BLEND_AT_PBL ) then
2170+
if( DAMPBEG.ne.DAMPEND .or. BLEND_AT_PBL .or. BLEND_QV_AT_TP ) then
21492171

21502172
if(first .and. MAPL_AM_I_ROOT()) then
21512173
if(DAMPBEG.ne.DAMPEND) then
@@ -2158,9 +2180,21 @@ subroutine handleANA_
21582180
else
21592181
print *, 'No blending based on PBL'
21602182
endif
2183+
if(BLEND_QV_AT_TP) then
2184+
if(BLEND_AT_PBL) then
2185+
print *, 'Blending at PBL supercedes QV blending at TROPP'
2186+
else
2187+
print *, 'Blending ANA and BKG QV based on TROPP'
2188+
endif
2189+
else
2190+
print *, 'No blending of QV based on TROPP'
2191+
endif
21612192
print *
21622193
endif
21632194

2195+
! Enforce priority now that warning has been logged
2196+
if ( BLEND_AT_PBL ) BLEND_QV_AT_TP = .FALSE.
2197+
21642198
if( BLEND_AT_PBL ) then
21652199
allocate ( pdum1(IMbkg,JMbkg,1) )
21662200
allocate ( pdum2(IM, JM, 1) )
@@ -2178,13 +2212,30 @@ subroutine handleANA_
21782212
endif
21792213
blnpp => pdum2(:,:,1)
21802214
endif
2215+
if( BLEND_QV_AT_TP ) then
2216+
allocate ( pdum1(IMbkg,JMbkg,1) )
2217+
allocate ( pdum2(IM, JM, 1) )
2218+
pdum1=0.0
2219+
2220+
call MAPL_GetPointer(import, ptr2d, 'TROPP_BLENDED', RC=STATUS)
2221+
VERIFY_(STATUS)
2222+
pdum1(:,:,1) = ptr2d
2223+
2224+
if (trim(GRIDINC)=="ANA" .and. do_transforms) then
2225+
call mkiau_internal_state%bkg2ana_regridder%regrid(pdum1, pdum2, RC=STATUS)
2226+
VERIFY_(STATUS)
2227+
else
2228+
pdum2=pdum1
2229+
endif
2230+
tropp => pdum2(:,:,1)
2231+
endif
21812232

21822233
call blend ( ple_ana,u_ana,v_ana,t_ana,q_ana,o3_ana, &
21832234
ple_bkg,u_bkg,v_bkg,t_bkg,q_bkg,o3_bkg, &
21842235
im,jm,LMbkg, DAMPBEG,DAMPEND, BLEND_AT_PBL, &
2185-
blnpp=blnpp )
2236+
BLEND_QV_AT_TP, blnpp=blnpp, tropp=tropp )
21862237

2187-
if( BLEND_AT_PBL ) then
2238+
if( BLEND_AT_PBL .or. BLEND_QV_AT_TP ) then
21882239
deallocate ( pdum1 )
21892240
deallocate ( pdum2 )
21902241
endif
@@ -2628,24 +2679,30 @@ end subroutine CLEAR
26282679
subroutine blend ( plea,ua,va,ta,qa,oa, &
26292680
pleb,ub,vb,tb,qb,ob, &
26302681
im,jm,lm, pabove,pbelow, &
2631-
BLEND_AT_PBL, blnpp )
2682+
BLEND_AT_PBL, &
2683+
BLEND_QV_AT_TP, &
2684+
blnpp, tropp )
26322685

26332686
! Blends Anaylsis and Background values.
2634-
! This routine is called if pabove /= pbelow or BLEND_AT_PBL
2687+
! This routine is called if pabove /= pbelow or BLEND_AT_PBL or BLEND_QV_AT_TP
26352688
! ***************************************************************************
26362689

26372690
implicit none
26382691
integer, intent(IN) :: im,jm,lm
26392692
real, intent(IN) :: pabove,pbelow
26402693
logical, intent(IN) :: BLEND_AT_PBL
2694+
logical, intent(IN) :: BLEND_QV_AT_TP
26412695

2696+
! Background values
26422697
real, intent(IN) :: pleb(im,jm,lm+1)
26432698
real, intent(IN) :: ub(im,jm,lm)
26442699
real, intent(IN) :: vb(im,jm,lm)
26452700
real, intent(IN) :: tb(im,jm,lm)
26462701
real, intent(IN) :: qb(im,jm,lm)
26472702
real, intent(IN) :: ob(im,jm,lm)
26482703

2704+
! IN: Anaylsis values
2705+
! OUT: Blended values
26492706
real, intent(INOUT) :: plea(im,jm,lm+1)
26502707
real, intent(INOUT) :: ua(im,jm,lm)
26512708
real, intent(INOUT) :: va(im,jm,lm)
@@ -2654,6 +2711,7 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
26542711
real, intent(INOUT) :: oa(im,jm,lm)
26552712

26562713
real, intent(IN), optional, pointer :: blnpp(:,:) ! blending pressure when BLEND_AT_PBL is TRUE
2714+
real, intent(IN), optional, pointer :: tropp(:,:) ! Tropopause Pressure used when BLEND_QV_AT_TP is TRUE
26572715

26582716
! Locals
26592717
! ------
@@ -2670,6 +2728,9 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
26702728
real pabove_BL,pbelow_BL
26712729
real bl_press
26722730

2731+
real pabove_QV,pbelow_QV ! compute from tropp
2732+
real tp_press
2733+
26732734
real alf,eps,p
26742735
integer i,j,L
26752736

@@ -2709,7 +2770,7 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
27092770
ua(i,j,L) = ub(i,j,L) + alf*( ua(i,j,L)- ub(i,j,L) )
27102771
va(i,j,L) = vb(i,j,L) + alf*( va(i,j,L)- vb(i,j,L) )
27112772
oa(i,j,L) = ob(i,j,L) + alf*( oa(i,j,L)- ob(i,j,L) )
2712-
qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )
2773+
IF (.NOT. BLEND_QV_AT_TP) qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )
27132774
enddo
27142775
enddo
27152776
enddo
@@ -2779,6 +2840,41 @@ subroutine blend ( plea,ua,va,ta,qa,oa, &
27792840
enddo
27802841
endif
27812842

2843+
! Use Analysis values of water vapor in the troposphere
2844+
! Relax to Background values in the stratosphere
2845+
! -----------------------------------------------------
2846+
if ( BLEND_QV_AT_TP ) then
2847+
do j=1,jm
2848+
do i=1,im
2849+
2850+
IF ( tropp(i,j) == MAPL_UNDEF ) THEN
2851+
tp_press = 100.0 * 100.0 ! 100 hPa
2852+
ELSE
2853+
tp_press = tropp(i,j)
2854+
ENDIF
2855+
2856+
pabove_QV = tp_press * 0.5
2857+
pbelow_QV = tp_press * 1.0
2858+
2859+
do L=1,lm
2860+
p = 0.5*( pleb(i,j,L)+pleb(i,j,L+1) )
2861+
if( p.le.pabove_QV ) then
2862+
alf = 0.0 ! use the background value
2863+
else if( p.gt.pabove_QV .and. p.le.pbelow_QV ) then
2864+
alf = (LOG(p) -LOG(pabove_QV))/ &
2865+
(LOG(pbelow_QV)-LOG(pabove_QV))
2866+
else
2867+
alf = 1.0 ! use the analysis value
2868+
endif
2869+
2870+
qa(i,j,L) = qb(i,j,L) + alf*( qa(i,j,L)- qb(i,j,L) )
2871+
2872+
enddo
2873+
2874+
enddo
2875+
enddo
2876+
endif
2877+
27822878
return
27832879
end subroutine blend
27842880

0 commit comments

Comments
 (0)