CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_suite_interstitial_2.F90
1
3
5
6 use machine, only: kind_phys
7 real(kind=kind_phys), parameter :: one = 1.0_kind_phys
8 logical :: linit_mod = .false.
9
10 contains
11
15 subroutine gfs_suite_interstitial_2_run (im, levs, lssav, ldiag3d, lsidea, flag_cice, shal_cnv, old_monin, mstrat, &
16 do_shoc, frac_grid, imfshalcnv, dtf, xcosz, adjsfcdsw, adjsfcdlw, cice, pgr, ulwsfc_cice, lwhd, htrsw, htrlw, xmu, ctei_rm, &
17 work1, work2, prsi, tgrs, prsl, qgrs_water_vapor, qgrs_cloud_water, cp, hvap, prslk, suntim, adjsfculw, adjsfculw_lnd, &
18 adjsfculw_ice, adjsfculw_wat, dlwsfc, ulwsfc, psmean, dtend, dtidx, index_of_process_longwave, index_of_process_shortwave, &
19 index_of_process_pbl, index_of_process_dcnv, index_of_process_scnv, index_of_process_mp, index_of_temperature, &
20 ctei_rml, ctei_r, kinver, dry, icy, wet, frland, huge, use_LW_jacobian, htrlwu, errmsg, errflg)
21
22 implicit none
23
24 ! interface variables
25 integer, intent(in ) :: im, levs, imfshalcnv
26 logical, intent(in ) :: lssav, ldiag3d, lsidea, shal_cnv
27 logical, intent(in ) :: old_monin, mstrat, do_shoc, frac_grid, use_LW_jacobian
28 real(kind=kind_phys), intent(in ) :: dtf, cp, hvap
29
30 logical, intent(in ), dimension(:) :: flag_cice
31 real(kind=kind_phys), intent(in ), dimension(:) :: ctei_rm
32 real(kind=kind_phys), intent(in ), dimension(:) :: xcosz, adjsfcdsw, adjsfcdlw, pgr, xmu, work1, work2
33 real(kind=kind_phys), intent(in ), dimension(:), optional :: ulwsfc_cice
34 real(kind=kind_phys), intent(in ), dimension(:) :: cice
35 real(kind=kind_phys), intent(in ), dimension(:,:) :: htrsw, htrlw, tgrs, prsl, qgrs_water_vapor, qgrs_cloud_water, prslk
36 real(kind=kind_phys), intent(in ), dimension(:,:), optional :: htrlwu
37 real(kind=kind_phys), intent(in ), dimension(:,:) :: prsi
38 real(kind=kind_phys), intent(in ), dimension(:,:,:) :: lwhd
39 integer, intent(inout), dimension(:) :: kinver
40 real(kind=kind_phys), intent(inout), dimension(:) :: suntim, dlwsfc, ulwsfc, psmean, ctei_rml, ctei_r
41 real(kind=kind_phys), intent(in ), dimension(:) :: adjsfculw_lnd, adjsfculw_ice, adjsfculw_wat
42 real(kind=kind_phys), intent(inout), dimension(:) :: adjsfculw
43
44 ! dtend is only allocated if ldiag3d is .true.
45 real(kind=kind_phys), optional, intent(inout), dimension(:,:,:) :: dtend
46 integer, intent(in), dimension(:,:) :: dtidx
47 integer, intent(in) :: index_of_process_longwave, index_of_process_shortwave, &
48 index_of_process_pbl, index_of_process_dcnv, index_of_process_scnv, &
49 index_of_process_mp, index_of_temperature
50
51 logical, intent(in ), dimension(:) :: dry, icy, wet
52 real(kind=kind_phys), intent(in ), dimension(:) :: frland
53 real(kind=kind_phys), intent(in ) :: huge
54
55 character(len=*), intent( out) :: errmsg
56 integer, intent( out) :: errflg
57
58 ! local variables
59 real(kind=kind_phys), parameter :: czmin = 0.0001_kind_phys ! cos(89.994)
60 integer :: i, k, idtend
61 real(kind=kind_phys) :: tem1, tem2, tem, hocp
62 logical, dimension(im) :: invrsn
63 real(kind=kind_phys), dimension(im) :: tx1, tx2
64
65 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys, one = 1.0_kind_phys
66 real(kind=kind_phys), parameter :: qmin = 1.0e-10_kind_phys, epsln=1.0e-10_kind_phys
67
68 ! Initialize CCPP error handling variables
69 errmsg = ''
70 errflg = 0
71
72 hocp = hvap/cp
73
74 if (lssav) then ! --- ... accumulate/save output variables
75
76! --- ... sunshine duration time is defined as the length of time (in mdl output
77! interval) that solar radiation falling on a plane perpendicular to the
78! direction of the sun >= 120 w/m2
79
80 do i = 1, im
81 if ( xcosz(i) >= czmin ) then ! zenth angle > 89.994 deg
82 tem1 = adjsfcdsw(i) / xcosz(i)
83 if ( tem1 >= 120.0_kind_phys ) then
84 suntim(i) = suntim(i) + dtf
85 endif
86 endif
87 enddo
88
89! --- ... sfc lw fluxes used by atmospheric model are saved for output
90 if (frac_grid) then
91 do i=1,im
92 tem = (one - frland(i)) * cice(i) ! tem = ice fraction wrt whole cell
93 if (flag_cice(i)) then
94 adjsfculw(i) = adjsfculw_lnd(i) * frland(i) &
95 + ulwsfc_cice(i) * tem &
96 + adjsfculw_wat(i) * (one - frland(i) - tem)
97 else
98 adjsfculw(i) = adjsfculw_lnd(i) * frland(i) &
99 + adjsfculw_ice(i) * tem &
100 + adjsfculw_wat(i) * (one - frland(i) - tem)
101 endif
102 enddo
103 else
104 do i=1,im
105 if (dry(i)) then ! all land
106 adjsfculw(i) = adjsfculw_lnd(i)
107 elseif (icy(i)) then ! ice (and water)
108 tem = one - cice(i)
109 if (flag_cice(i)) then
110 if (wet(i) .and. adjsfculw_wat(i) /= huge) then
111 adjsfculw(i) = ulwsfc_cice(i)*cice(i) + adjsfculw_wat(i)*tem
112 else
113 adjsfculw(i) = ulwsfc_cice(i)
114 endif
115 else
116 if (wet(i) .and. adjsfculw_wat(i) /= huge) then
117 adjsfculw(i) = adjsfculw_ice(i)*cice(i) + adjsfculw_wat(i)*tem
118 else
119 adjsfculw(i) = adjsfculw_ice(i)
120 endif
121 endif
122 else ! all water
123 adjsfculw(i) = adjsfculw_wat(i)
124 endif
125 enddo
126 endif
127
128 do i=1,im
129 dlwsfc(i) = dlwsfc(i) + adjsfcdlw(i)*dtf
130 ulwsfc(i) = ulwsfc(i) + adjsfculw(i)*dtf
131 psmean(i) = psmean(i) + pgr(i)*dtf ! mean surface pressure
132 enddo
133
134 if (ldiag3d) then
135 if (lsidea) then
136 idtend = dtidx(index_of_temperature,index_of_process_longwave)
137 if(idtend>=1) then
138 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,1)*dtf
139 endif
140
141 idtend = dtidx(index_of_temperature,index_of_process_shortwave)
142 if(idtend>=1) then
143 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,2)*dtf
144 endif
145
146 idtend = dtidx(index_of_temperature,index_of_process_pbl)
147 if(idtend>=1) then
148 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,3)*dtf
149 endif
150
151 idtend = dtidx(index_of_temperature,index_of_process_dcnv)
152 if(idtend>=1) then
153 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,4)*dtf
154 endif
155
156 idtend = dtidx(index_of_temperature,index_of_process_scnv)
157 if(idtend>=1) then
158 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,5)*dtf
159 endif
160
161 idtend = dtidx(index_of_temperature,index_of_process_mp)
162 if(idtend>=1) then
163 dtend(:,:,idtend) = dtend(:,:,idtend) + lwhd(:,:,6)*dtf
164 endif
165 else
166 idtend = dtidx(index_of_temperature,index_of_process_longwave)
167 if(idtend>=1) then
168 if (use_lw_jacobian) then
169 dtend(:,:,idtend) = dtend(:,:,idtend) + htrlwu(:,:)*dtf
170 else
171 dtend(:,:,idtend) = dtend(:,:,idtend) + htrlw(:,:)*dtf
172 endif
173 endif
174
175 idtend = dtidx(index_of_temperature,index_of_process_shortwave)
176 if(idtend>=1) then
177 do k=1,levs
178 do i=1,im
179 dtend(i,k,idtend) = dtend(i,k,idtend) + htrsw(i,k)*dtf*xmu(i)
180 enddo
181 enddo
182 endif
183 endif
184 endif
185 endif ! end if_lssav_block
186
187 do i=1, im
188 invrsn(i) = .false.
189 tx1(i) = zero
190 tx2(i) = 10.0_kind_phys
191 ctei_r(i) = 10.0_kind_phys
192 enddo
193
194 if ((((imfshalcnv == 0 .and. shal_cnv) .or. old_monin) .and. mstrat) &
195 .or. do_shoc) then
196 ctei_rml(:) = ctei_rm(1)*work1(:) + ctei_rm(2)*work2(:)
197 do k=1,levs/2
198 do i=1,im
199 if (prsi(i,1)-prsi(i,k+1) < 0.35_kind_phys*prsi(i,1) &
200 .and. (.not. invrsn(i))) then
201 tem = (tgrs(i,k+1) - tgrs(i,k)) &
202 / (prsl(i,k) - prsl(i,k+1))
203
204 if (((tem > 0.0001_kind_phys) .and. (tx1(i) < zero)) .or. &
205 ((tem-abs(tx1(i)) > zero) .and. (tx2(i) < zero))) then
206 invrsn(i) = .true.
207
208 if (qgrs_water_vapor(i,k) > qgrs_water_vapor(i,k+1)) then
209 tem1 = tgrs(i,k+1) + hocp*max(qgrs_water_vapor(i,k+1),qmin)
210 tem2 = tgrs(i,k) + hocp*max(qgrs_water_vapor(i,k),qmin)
211
212 tem1 = tem1 / prslk(i,k+1) - tem2 / prslk(i,k)
213
214! --- ... (cp/l)(deltathetae)/(deltatwater) > ctei_rm -> conditon for CTEI
215 ctei_r(i) = (one/hocp)*tem1/(qgrs_water_vapor(i,k+1)-qgrs_water_vapor(i,k) &
216 + qgrs_cloud_water(i,k+1)-qgrs_cloud_water(i,k))
217 else
218 ctei_r(i) = 10.0_kind_phys
219 endif
220
221 if ( ctei_rml(i) > ctei_r(i) ) then
222 kinver(i) = k
223 else
224 kinver(i) = levs
225 endif
226 endif
227
228 tx2(i) = tx1(i)
229 tx1(i) = tem
230 endif
231 enddo
232 enddo
233 endif
234
235 end subroutine gfs_suite_interstitial_2_run
236
237 end module gfs_suite_interstitial_2