00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 module ice_atmo
00022
00023
00024
00025 use ice_kinds_mod
00026 use ice_constants
00027
00028
00029
00030 implicit none
00031 save
00032
00033 character (len=char_len) ::
00034 atmbndy
00035
00036 logical (kind=log_kind) ::
00037 calc_strair
00038
00039
00040
00041 contains
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051 subroutine atmo_boundary_layer (nx_block, ny_block, &
00052 sfctype, icells, &
00053 indxi, indxj, &
00054 Tsf, potT, &
00055 uatm, vatm, &
00056 uvel, vvel, &
00057 wind, zlvl, &
00058 Qa, rhoa, &
00059 strx, stry, &
00060 Tref, Qref, &
00061 delt, delq, &
00062 lhcoef, shcoef)
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083 integer (kind=int_kind), intent(in) ::
00084 nx_block, ny_block,
00085 icells
00086
00087 integer (kind=int_kind), dimension(nx_block*ny_block),
00088 intent(in) ::
00089 indxi, indxj
00090
00091 character (len=3), intent(in) ::
00092 sfctype
00093
00094 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
00095 Tsf ,
00096 potT ,
00097 uatm ,
00098 vatm ,
00099 uvel ,
00100 vvel ,
00101 wind ,
00102 zlvl ,
00103 Qa ,
00104 rhoa
00105
00106 real (kind=dbl_kind), dimension (nx_block,ny_block),
00107 intent(inout) ::
00108 strx ,
00109 stry
00110
00111 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(out) ::
00112 Tref ,
00113 Qref ,
00114 delt ,
00115 delq ,
00116 shcoef ,
00117 lhcoef
00118
00119
00120
00121 integer (kind=int_kind) ::
00122 k ,
00123 i, j ,
00124 ij
00125
00126 real (kind=dbl_kind) ::
00127 TsfK ,
00128 xqq ,
00129 psimh ,
00130 tau ,
00131 fac ,
00132 al2 ,
00133 psix2 ,
00134 psimhs,
00135 ssq ,
00136 qqq ,
00137 TTT ,
00138 qsat ,
00139 Lheat
00140
00141 real (kind=dbl_kind), dimension (icells) ::
00142 ustar ,
00143 tstar ,
00144 qstar ,
00145 rdn ,
00146 rhn ,
00147 ren ,
00148 rd ,
00149 re ,
00150 rh ,
00151 vmag ,
00152 alz ,
00153 thva ,
00154 cp ,
00155 hol ,
00156 stable,
00157 psixh
00158
00159 real (kind=dbl_kind), parameter ::
00160 cpvir = cp_wv/cp_air-c1,
00161 zTrf = c2 ,
00162 umin = c1
00163
00164
00165 real (kind=dbl_kind) ::
00166 xd ,
00167 psimhu,
00168 psixhu
00169
00170
00171
00172
00173
00174 psimhu(xd) = log((c1+xd*(c2+xd))*(c1+xd*xd)/c8) &
00175 - c2*atan(xd) + pih
00176
00177
00178 psixhu(xd) = c2 * log((c1 + xd*xd)/c2)
00179
00180 al2 = log(zref/zTrf)
00181
00182
00183
00184
00185
00186 do j = 1, ny_block
00187 do i = 1, nx_block
00188 Tref(i,j) = c0
00189 Qref(i,j) = c0
00190 delt(i,j) = c0
00191 delq(i,j) = c0
00192 shcoef(i,j) = c0
00193 lhcoef(i,j) = c0
00194 enddo
00195 enddo
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205
00206 if (sfctype(1:3)=='ice') then
00207
00208 qqq = qqqice
00209 TTT = TTTice
00210 Lheat = Lsub
00211 do ij = 1, icells
00212 i = indxi(ij)
00213 j = indxj(ij)
00214 vmag(ij) = max(umin, wind(i,j))
00215 rdn(ij) = vonkar/log(zref/iceruf)
00216 enddo
00217
00218 elseif (sfctype(1:3)=='ocn') then
00219
00220 qqq = qqqocn
00221 TTT = TTTocn
00222 Lheat = Lvap
00223 do ij = 1, icells
00224 i = indxi(ij)
00225 j = indxj(ij)
00226 vmag(ij) = max(umin, wind(i,j))
00227 rdn(ij) = sqrt(0.0027_dbl_kind/vmag(ij) &
00228 + .000142_dbl_kind + .0000764_dbl_kind*vmag(ij))
00229 enddo
00230
00231 endif
00232
00233 do ij = 1, icells
00234 i = indxi(ij)
00235 j = indxj(ij)
00236
00237
00238
00239
00240
00241 TsfK = Tsf(i,j) + Tffresh
00242 qsat = qqq * exp(-TTT/TsfK)
00243 ssq = qsat / rhoa(i,j)
00244
00245 thva(ij) = potT(i,j) * (c1 + zvir * Qa(i,j))
00246 delt(i,j) = potT(i,j) - TsfK
00247 delq(i,j) = Qa(i,j) - ssq
00248 alz(ij) = log(zlvl(i,j)/zref)
00249 cp(ij) = cp_air*(c1 + cpvir*ssq)
00250
00251
00252
00253
00254
00255
00256 rhn(ij) = rdn(ij)
00257 ren(ij) = rdn(ij)
00258
00259
00260 ustar(ij) = rdn(ij) * vmag(ij)
00261 tstar(ij) = rhn(ij) * delt(i,j)
00262 qstar(ij) = ren(ij) * delq(i,j)
00263
00264 enddo
00265
00266
00267
00268
00269
00270 do k=1,5
00271
00272 do ij = 1, icells
00273 i = indxi(ij)
00274 j = indxj(ij)
00275
00276
00277 hol(ij) = vonkar * gravit * zlvl(i,j) &
00278 * (tstar(ij)/thva(ij) &
00279 + qstar(ij)/(c1/zvir+Qa(i,j))) &
00280 / ustar(ij)**2
00281 hol(ij) = sign( min(abs(hol(ij)),c10), hol(ij) )
00282 stable(ij) = p5 + sign(p5 , hol(ij))
00283 xqq = max(sqrt(abs(c1 - c16*hol(ij))) , c1)
00284 xqq = sqrt(xqq)
00285
00286
00287 psimhs = -(0.7_dbl_kind*hol(ij) &
00288 + 0.75_dbl_kind*(hol(ij)-14.3_dbl_kind) &
00289 * exp(-0.35_dbl_kind*hol(ij)) + 10.7_dbl_kind)
00290 psimh = psimhs*stable(ij) &
00291 + (c1 - stable(ij))*psimhu(xqq)
00292 psixh(ij) = psimhs*stable(ij) &
00293 + (c1 - stable(ij))*psixhu(xqq)
00294
00295
00296 rd(ij) = rdn(ij) / (c1+rdn(ij)/vonkar*(alz(ij)-psimh))
00297 rh(ij) = rhn(ij) / (c1+rhn(ij)/vonkar*(alz(ij)-psixh(ij)))
00298 re(ij) = ren(ij) / (c1+ren(ij)/vonkar*(alz(ij)-psixh(ij)))
00299
00300
00301 ustar(ij) = rd(ij) * vmag(ij)
00302 tstar(ij) = rh(ij) * delt(i,j)
00303 qstar(ij) = re(ij) * delq(i,j)
00304
00305 enddo
00306 enddo
00307
00308 if (calc_strair) then
00309
00310
00311 do j = 1, ny_block
00312 do i = 1, nx_block
00313 strx(i,j) = c0
00314 stry(i,j) = c0
00315 enddo
00316 enddo
00317
00318
00319
00320
00321 do ij = 1, icells
00322 i = indxi(ij)
00323 j = indxj(ij)
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333 tau = rhoa(i,j) * ustar(ij) * rd(ij)
00334 strx(i,j) = tau * (uatm(i,j)-uvel(i,j))
00335 stry(i,j) = tau * (vatm(i,j)-vvel(i,j))
00336
00337 enddo
00338
00339 endif
00340
00341
00342
00343
00344 do ij = 1, icells
00345 i = indxi(ij)
00346 j = indxj(ij)
00347
00348
00349
00350
00351
00352
00353
00354
00355 shcoef(i,j) = rhoa(i,j) * ustar(ij) * cp(ij) * rh(ij) + c1
00356 lhcoef(i,j) = rhoa(i,j) * ustar(ij) * Lheat * re(ij)
00357
00358
00359
00360
00361 hol(ij) = hol(ij)*zTrf/zlvl(i,j)
00362 xqq = max( c1, sqrt(abs(c1-c16*hol(ij))) )
00363 xqq = sqrt(xqq)
00364 psix2 = -c5*hol(ij)*stable(ij) + (c1-stable(ij))*psixhu(xqq)
00365 fac = (rh(ij)/vonkar) &
00366 * (alz(ij) + al2 - psixh(ij) + psix2)
00367 Tref(i,j)= potT(i,j) - delt(i,j)*fac
00368 Tref(i,j)= Tref(i,j) - p01*zTrf
00369 fac = (re(ij)/vonkar) &
00370 * (alz(ij) + al2 - psixh(ij) + psix2)
00371 Qref(i,j)= Qa(i,j) - delq(i,j)*fac
00372
00373 enddo
00374
00375 end subroutine atmo_boundary_layer
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385 subroutine atmo_boundary_const (nx_block, ny_block, &
00386 sfctype, icells, &
00387 indxi, indxj, &
00388 uatm, vatm, &
00389 wind, rhoa, &
00390 strx, stry, &
00391 lhcoef, shcoef)
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405
00406 integer (kind=int_kind), intent(in) ::
00407 nx_block, ny_block,
00408 icells
00409
00410 integer (kind=int_kind), dimension(nx_block*ny_block),
00411 intent(in) ::
00412 indxi, indxj
00413
00414 character (len=3), intent(in) ::
00415 sfctype
00416
00417 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
00418 uatm ,
00419 vatm ,
00420 wind ,
00421 rhoa
00422
00423 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(inout)::
00424 strx ,
00425 stry
00426
00427 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(out)::
00428 shcoef ,
00429 lhcoef
00430
00431
00432
00433 integer (kind=int_kind) ::
00434 i, j,
00435 ij
00436
00437 real (kind=dbl_kind) ::
00438 tau,
00439 Lheat
00440
00441
00442
00443
00444
00445 do j = 1, ny_block
00446 do i = 1, nx_block
00447 shcoef(i,j) = c0
00448 lhcoef(i,j) = c0
00449 enddo
00450 enddo
00451
00452 if (calc_strair) then
00453
00454 do j = 1, ny_block
00455 do i = 1, nx_block
00456 strx(i,j) = c0
00457 stry(i,j) = c0
00458 enddo
00459 enddo
00460
00461
00462
00463
00464 do ij = 1, icells
00465 i = indxi(ij)
00466 j = indxj(ij)
00467
00468
00469
00470
00471
00472 tau = rhoa(i,j) * 0.0012_dbl_kind * wind(i,j)
00473
00474
00475 strx(i,j) = tau * uatm(i,j)
00476 stry(i,j) = tau * vatm(i,j)
00477
00478 enddo
00479
00480 endif
00481
00482
00483
00484
00485
00486 if (sfctype(1:3)=='ice') then
00487 Lheat = Lsub
00488 elseif (sfctype(1:3)=='ocn') then
00489 Lheat = Lvap
00490 endif
00491
00492
00493
00494
00495 do ij = 1, icells
00496 i = indxi(ij)
00497 j = indxj(ij)
00498
00499
00500
00501
00502
00503 shcoef(i,j) = (1.20e-3_dbl_kind)*cp_air*rhoa(i,j)*wind(i,j)
00504 lhcoef(i,j) = (1.50e-3_dbl_kind)*Lheat *rhoa(i,j)*wind(i,j)
00505
00506 enddo
00507
00508 end subroutine atmo_boundary_const
00509
00510
00511
00512 end module ice_atmo
00513
00514