CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
gcycle.F90
1
4
6
7 implicit none
8
9 private
10
11 public gcycle
12
13contains
14
18 subroutine gcycle (me, nthrds, nx, ny, isc, jsc, nsst, tile_num, nlunit, fn_nml, &
19 input_nml_file, lsoil, lsoil_lsm, kice, idate, ialb, isot, ivegsrc, &
20 use_ufo, nst_anl, fhcyc, phour, landfrac, lakefrac, min_seaice, min_lakeice, &
21 frac_grid, smc, slc, stc, smois, sh2o, tslb, tiice, tg3, tref, tsfc, &
22 tsfco, tisfc, hice, fice, facsf, facwf, alvsf, alvwf, alnsf, alnwf, &
23 zorli, zorll, zorlo, weasd, slope, snoalb, canopy, vfrac, vtype, &
24 stype, scolor, shdmin, shdmax, snowd, cv, cvb, cvt, oro, oro_uf, &
25 xlat_d, xlon_d, slmsk, imap, jmap, errmsg, errflg)
26!
27!
28 use machine, only: kind_phys, kind_io8
29 implicit none
30
31 integer, intent(in) :: me, nthrds, nx, ny, isc, jsc, nsst, &
32 tile_num, nlunit, lsoil, lsoil_lsm, kice
33 integer, intent(in) :: idate(:), ialb, isot, ivegsrc
34 character(len = 64), intent(in) :: fn_nml
35 character(len=*), intent(in) :: input_nml_file(:)
36 logical, intent(in) :: use_ufo, nst_anl, frac_grid
37 real(kind=kind_phys), intent(in) :: fhcyc, phour, landfrac(:), lakefrac(:), &
38 min_seaice, min_lakeice, &
39 xlat_d(:), xlon_d(:)
40 real(kind=kind_phys), intent(inout), optional :: &
41 smois(:,:), &
42 sh2o(:,:), &
43 tslb(:,:), &
44 tref(:)
45 real(kind=kind_phys), intent(inout) :: smc(:,:), &
46 slc(:,:), &
47 stc(:,:), &
48 tiice(:,:), &
49 tg3(:), &
50 tsfc(:), &
51 tsfco(:), &
52 tisfc(:), &
53 hice(:), &
54 fice(:), &
55 facsf(:), &
56 facwf(:), &
57 alvsf(:), &
58 alvwf(:), &
59 alnsf(:), &
60 alnwf(:), &
61 zorli(:), &
62 zorll(:), &
63 zorlo(:), &
64 weasd(:), &
65 snoalb(:), &
66 canopy(:), &
67 vfrac(:), &
68 shdmin(:), &
69 shdmax(:), &
70 snowd(:), &
71 cv(:), &
72 cvb(:), &
73 cvt(:), &
74 oro(:), &
75 oro_uf(:), &
76 slmsk(:)
77 integer, intent(inout) :: vtype(:), &
78 stype(:), &
79 scolor(:), &
80 slope(:)
81
82 integer, intent(in) :: imap(:), jmap(:)
83 character(len=*), intent(out) :: errmsg
84 integer, intent(out) :: errflg
85
86!
87! Local variables
88! ---------------
89! real(kind=kind_phys) :: &
90 real(kind=kind_io8) :: &
91 slmskl(nx*ny), &
92 slmskw(nx*ny), &
93 slpfcs(nx*ny), &
94 vegfcs(nx*ny), &
95 sltfcs(nx*ny), &
96 slcfcs(nx*ny), & !soil color
97 tsffcs(nx*ny), &
98 zorfcs(nx*ny), &
99 aisfcs(nx*ny), &
100 alffc1(nx*ny*2), &
101 albfc1(nx*ny*4), &
102 smcfc1(nx*ny*max(lsoil,lsoil_lsm)), &
103 stcfc1(nx*ny*max(lsoil,lsoil_lsm)), &
104 slcfc1(nx*ny*max(lsoil,lsoil_lsm))
105
106
107 real (kind=kind_io8) :: min_ice(nx*ny)
108 integer :: i_indx(nx*ny), j_indx(nx*ny)
109 character(len=6) :: tile_num_ch
110 real(kind=kind_phys) :: sig1t
111 integer :: npts, nb, ix, jx, ls, ios, ll
112 logical :: exists
113
114 ! Initialize CCPP error handling variables
115 errmsg = ''
116 errflg = 0
117
118!
119!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
120!
121! if (Model%me .eq. 0) print *,' nlats=',nlats,' lonsinpe='
122! *,lonsinpe(0,1)
123!
124 tile_num_ch = " "
125 if (tile_num < 10) then
126 write(tile_num_ch, "(a4,i1)") "tile", tile_num
127 else
128 write(tile_num_ch, "(a4,i2)") "tile", tile_num
129 endif
130!
131 sig1t = 0.0_kind_phys
132 npts = nx*ny
133!
134 if ( nsst > 0 ) then
135 tsffcs = tref
136 else
137 tsffcs = tsfco
138 end if
139! integer to real/double precision
140 slpfcs = real(slope)
141 vegfcs = real(vtype)
142 sltfcs = real(stype)
143 slcfcs = real(scolor) !soil color
144!
145 if (frac_grid) then
146 do ix=1,npts
147! if (landfrac(ix) > -1.0e-8_kind_phys) then
148 if (landfrac(ix) > 0.0_kind_phys) then
149 slmskl(ix) = ceiling(landfrac(ix)-1.0e-6_kind_phys)
150 slmskw(ix) = floor(landfrac(ix)+1.0e-6_kind_phys)
151 else
152 if (nint(slmsk(ix)) == 1) then
153 slmskl(ix) = 1.0_kind_phys
154 slmskw(ix) = 1.0_kind_phys
155 else
156 slmskl(ix) = 0.0_kind_phys
157 slmskw(ix) = 0.0_kind_phys
158 endif
159 endif
160 zorfcs(ix) = zorll(ix)
161 if (nint(slmskl(ix)) == 0) then
162 if (slmsk(ix) > 1.99_kind_phys) then
163 zorfcs(ix) = zorli(ix)
164 else
165 zorfcs(ix) = zorlo(ix)
166 endif
167 endif
168 enddo
169 else
170 do ix=1,npts
171 if (nint(slmsk(ix)) == 1) then
172 slmskl(ix) = 1.0_kind_phys
173 slmskw(ix) = 1.0_kind_phys
174 else
175 slmskl(ix) = 0.0_kind_phys
176 slmskw(ix) = 0.0_kind_phys
177 endif
178 zorfcs(ix) = zorll(ix)
179 if (slmsk(ix) > 1.99_kind_phys) then
180 zorfcs(ix) = zorli(ix)
181 elseif (slmsk(ix) < 0.1_kind_phys) then
182 zorfcs(ix) = zorlo(ix)
183 endif
184 enddo
185 endif
186 do ix=1,npts
187 i_indx(ix) = imap(ix) + isc - 1
188 j_indx(ix) = jmap(ix) + jsc - 1
189
190 if (lakefrac(ix) > 0.0_kind_phys) then
191 min_ice(ix) = min_lakeice
192 else
193 min_ice(ix) = min_seaice
194 endif
195
196 IF (slmsk(ix) > 1.99_kind_phys) THEN
197 aisfcs(ix) = 1.0_kind_phys
198 ELSE
199 aisfcs(ix) = 0.0_kind_phys
200 ENDIF
201 !
202 alffc1(ix ) = facsf(ix)
203 alffc1(ix + npts ) = facwf(ix)
204 !
205 albfc1(ix ) = alvsf(ix)
206 albfc1(ix + npts ) = alvwf(ix)
207 albfc1(ix + npts*2) = alnsf(ix)
208 albfc1(ix + npts*3) = alnwf(ix)
209 !
210 do ls = 1,max(lsoil,lsoil_lsm)
211 ll = ix + (ls-1)*npts
212 if (lsoil == lsoil_lsm) then
213 smcfc1(ll) = smc(ix,ls)
214 stcfc1(ll) = stc(ix,ls)
215 slcfc1(ll) = slc(ix,ls)
216 else
217 smcfc1(ll) = smois(ix,ls)
218 stcfc1(ll) = tslb(ix,ls)
219 slcfc1(ll) = sh2o(ix,ls)
220 endif
221 enddo
222!
223 enddo
224!
225#ifndef INTERNAL_FILE_NML
226 inquire (file=trim(fn_nml),exist=exists)
227 if (.not. exists) then
228 write(6,*) 'gcycle:: namelist file: ',trim(fn_nml),' does not exist'
229 errflg = 1
230 errmsg = 'ERROR(gcycle): namelist file: ',trim(fn_nml),' does not exist.'
231 return
232 else
233 open (unit=nlunit, file=trim(fn_nml), action='READ', status='OLD', iostat=ios)
234 rewind(nlunit)
235 endif
236#endif
237 CALL sfccycle (9998, npts, max(lsoil,lsoil_lsm), sig1t, fhcyc, &
238 idate(4), idate(2), idate(3), idate(1), &
239 phour, xlat_d, xlon_d, slmskl, slmskw, &
240 oro, oro_uf, use_ufo, nst_anl, &
241 hice, fice, tisfc, snowd, slcfc1, &
242 shdmin, shdmax, slpfcs, snoalb, tsffcs, &
243 weasd, zorfcs, albfc1, tg3, canopy, &
244 smcfc1, stcfc1, slmsk, aisfcs, &
245 vfrac, vegfcs, sltfcs, slcfcs,alffc1, cv, & !slcfcs: soil color
246 cvb, cvt, me, nthrds, &
247 nlunit, size(input_nml_file), input_nml_file, &
248 min_ice, ialb, isot, ivegsrc, &
249 trim(tile_num_ch), i_indx, j_indx)
250#ifndef INTERNAL_FILE_NML
251 close (nlunit)
252#endif
253!
254 if ( nsst > 0 ) then
255 tref = tsffcs
256 else
257 tsfc = tsffcs
258 tsfco = tsffcs
259 endif
260!
261! real/double precision to integer
262 slope = int(slpfcs)
263 vtype = int(vegfcs)
264 stype = int(sltfcs)
265 scolor = int(slcfcs) !soil color
266!
267 do ix=1,npts
268 zorll(ix) = zorfcs(ix)
269 if (nint(slmskl(ix)) == 0) then
270 if (slmsk(ix) > 1.99_kind_phys) then
271 zorli(ix) = zorfcs(ix)
272 elseif (slmsk(ix) < 0.1_kind_phys) then
273 zorlo(ix) = zorfcs(ix)
274 endif
275 endif
276 !
277 facsf(ix) = alffc1(ix )
278 facwf(ix) = alffc1(ix + npts )
279 !
280 alvsf(ix) = albfc1(ix )
281 alvwf(ix) = albfc1(ix + npts )
282 alnsf(ix) = albfc1(ix + npts*2)
283 alnwf(ix) = albfc1(ix + npts*3)
284 !
285 do ls = 1,max(lsoil,lsoil_lsm)
286 ll = ix + (ls-1)*npts
287 if(lsoil == lsoil_lsm) then
288 smc(ix,ls) = smcfc1(ll)
289 stc(ix,ls) = stcfc1(ll)
290 slc(ix,ls) = slcfc1(ll)
291 else
292 smois(ix,ls) = smcfc1(ll)
293 tslb(ix,ls) = stcfc1(ll)
294 sh2o(ix,ls) = slcfc1(ll)
295 endif
296! if (ls<=kice) tiice(ix,ls) = STCFC1(ll)
297 enddo
298 enddo
299!
300! if (Model%me .eq. 0) print*,'executed gcycle during hour=',fhour
301!
302 RETURN
303 END
304
305end module gcycle_mod
subroutine sfccycle(lugb, len, lsoil, sig1t, deltsfc, iy, im, id, ih, fh, rla, rlo, slmskl, slmskw, orog, orog_uf, use_ufo, nst_anl, sihfcs, sicfcs, sitfcs, swdfcs, slcfcs, vmnfcs, vmxfcs, slpfcs, absfcs, tsffcs, snofcs, zorfcs, albfcs, tg3fcs, cnpfcs, smcfcs, stcfcs, slifcs, aisfcs, vegfcs, vetfcs, sotfcs, socfcs, alffcs, cvfcs, cvbfcs, cvtfcs, me, nthrds, nlunit, sz_nml, input_nml_file, min_ice, ialb, isot, ivegsrc, tile_num_ch, i_index, j_index)
This subroutine reads or interpolates surface climatology data in analysis and forecast mode.
Definition sfcsub.F:90
subroutine, public gcycle(me, nthrds, nx, ny, isc, jsc, nsst, tile_num, nlunit, fn_nml, input_nml_file, lsoil, lsoil_lsm, kice, idate, ialb, isot, ivegsrc, use_ufo, nst_anl, fhcyc, phour, landfrac, lakefrac, min_seaice, min_lakeice, frac_grid, smc, slc, stc, smois, sh2o, tslb, tiice, tg3, tref, tsfc, tsfco, tisfc, hice, fice, facsf, facwf, alvsf, alvwf, alnsf, alnwf, zorli, zorll, zorlo, weasd, slope, snoalb, canopy, vfrac, vtype, stype, scolor, shdmin, shdmax, snowd, cv, cvb, cvt, oro, oro_uf, xlat_d, xlon_d, slmsk, imap, jmap, errmsg, errflg)
This subroutine repopulates specific time-varying surface properties for atmospheric forecast runs.
Definition gcycle.F90:26