Changeset 827 for branches/dev_002_LIM/NEMO/OPA_SRC/SBC/flxblk.F90
- Timestamp:
- 2008-03-06T15:39:49+01:00 (16 years ago)
- File:
-
- 1 edited
-
branches/dev_002_LIM/NEMO/OPA_SRC/SBC/flxblk.F90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_002_LIM/NEMO/OPA_SRC/SBC/flxblk.F90
r826 r827 25 25 USE albedo 26 26 USE prtctl ! Print control 27 27 #if defined key_lim3 28 USE ice 29 #elif defined key_lim2 30 USE ice_2 31 #endif 28 32 IMPLICIT NONE 29 33 PRIVATE … … 79 83 80 84 CONTAINS 85 #if ! defined key_lim2 86 87 SUBROUTINE flx_blk( psst ) 88 !!--------------------------------------------------------------------------- 89 !! *** ROUTINE flx_blk *** 90 !! 91 !! ** Purpose : Computation of the heat fluxes at ocean and snow/ice 92 !! surface the solar heat at ocean and snow/ice surfaces and the 93 !! sensitivity of total heat fluxes to the SST variations 94 !! 95 !! ** Method : The flux of heat at the ice and ocean surfaces are derived 96 !! from semi-empirical ( or bulk ) formulae which relate the flux to 97 !! the properties of the surface and of the lower atmosphere. Here, we 98 !! follow the work of Oberhuber, 1988 99 !! 100 !! ** Action : call flx_blk_albedo to compute ocean and ice albedo 101 !! computation of snow precipitation 102 !! computation of solar flux at the ocean and ice surfaces 103 !! computation of the long-wave radiation for the ocean and sea/ice 104 !! computation of turbulent heat fluxes over water and ice 105 !! computation of evaporation over water 106 !! computation of total heat fluxes sensitivity over ice (dQ/dT) 107 !! computation of latent heat flux sensitivity over ice (dQla/dT) 108 !! 109 !! History : 110 !! 8.0 ! 97-06 (Louvain-La-Neuve) Original code 111 !! 8.5 ! 02-09 (C. Ethe , G. Madec ) F90: Free form and module 112 !!---------------------------------------------------------------------- 113 !! * Arguments 114 REAL(wp), INTENT( in ), DIMENSION(jpi,jpj) :: & 115 & psst ! Sea Surface Temperature 116 117 !! * Local variables 118 INTEGER :: & 119 ji, jj, jl, jt , & ! dummy loop indices 120 indaet , & ! = -1, 0, 1 for odd, normal and leap years resp. 121 iday , & ! integer part of day 122 indxb , & ! index for budyko coefficient 123 indxc ! index for cloud depth coefficient 124 125 REAL(wp) :: & 126 zalat , zclat , & ! latitude in degrees 127 zmt1, zmt2, zmt3 , & ! tempory air temperatures variables 128 ztatm3, ztatm4 , & ! power 3 and 4 of air temperature 129 z4tatm3 , & ! 4 * ztatm3 130 zcmue , & ! cosine of local solar altitude 131 zcmue2 , & ! root of zcmue1 132 zscmue , & ! square-root of zcmue1 133 zpcmue , & ! zcmue1**1.4 134 zdecl , & ! solar declination 135 zsdecl , zcdecl , & ! sine and cosine of solar declination 136 zalbo , & ! albedo of sea-water 137 zalbi , & ! albedo of ice 138 ztamr , & ! air temperature minus triple point of water (rtt) 139 ztaevbk , & ! part of net longwave radiation 140 zevi , zevo , & ! vapour pressure of ice and ocean 141 zind1,zind2,zind3 , & ! switch for testing the values of air temperature 142 zinda , & ! switch for testing the values of sea ice cover 143 zpis2 , & ! pi / 2 144 z2pi ! 2 * pi 145 146 REAL(wp) :: & 147 zxday , & ! day of year 148 zdist , & ! distance between the sun and the earth during the year 149 zdaycor , & ! corr. factor to take into account the variation of 150 ! ! zday when calc. the solar rad. 151 zesi, zeso , & ! vapour pressure of ice and ocean at saturation 152 zesi2 , & ! root of zesi 153 zqsato , & ! humidity close to the ocean surface (at saturation) 154 zqsati , & ! humidity close to the ice surface (at saturation) 155 zqsati2 , & ! root of zqsati 156 zdesidt , & ! derivative of zesi, function of ice temperature 157 zdteta , & ! diff. betw. sst and air temperature 158 zdeltaq , & ! diff. betw. spec. hum. and hum. close to the surface 159 ztvmoy, zobouks , & ! tempory scalars 160 zpsims, zpsihs, zpsils, zobouku, zxins, zpsimu , & 161 zpsihu, zpsilu, zstab,zpsim, zpsih, zpsil , & 162 zvatmg, zcmn, zchn, zcln, zcmcmn, zdenum , & 163 zdtetar, ztvmoyr, zlxins, zcmn2, zchcm, zclcm , zcoef 164 165 REAL(wp) :: & 166 zrhova , & ! air density per wind speed 167 zcsho , zcleo , & ! transfer coefficient over ocean 168 zcshi , zclei , & ! transfer coefficient over ice-free 169 zrhovacleo , & ! air density per wind speed per transfer coef. 170 zrhovacsho, zrhovaclei, zrhovacshi, & 171 ztice3 , & ! power 3 of ice temperature 172 zticemb, zticemb2 , & ! tempory air temperatures variables 173 zdqlw_ice , & ! sensitivity of long-wave flux over ice 174 zdqsb_ice , & ! sensitivity of sensible heat flux over ice 175 zdqla_ice , & ! sensitivity of latent heat flux over ice 176 zdl, zdr ! fractionnal part of latitude 177 REAL(wp), DIMENSION(jpi,jpj) :: & 178 zpatm , & ! atmospheric pressure 179 zqatm , & ! specific humidity 180 zes , & ! vapour pressure at saturation 181 zev, zevsqr , & ! vapour pressure and his square-root 182 zrhoa , & ! air density 183 ztatm , & ! air temperature in Kelvins 184 zfrld , & ! fraction of sea ice cover 185 zcatm1 , & ! fraction of cloud 186 zcldeff ! correction factor to account cloud effect 187 REAL(wp), DIMENSION(jpi,jpj) :: & 188 zalbocsd , & ! albedo of ocean 189 zalboos , & ! albedo of ocean under overcast sky 190 zalbomu , & ! albedo of ocean when zcmue is 0.4 191 zqsro , & ! solar radiation over ocean 192 zqsrics , & ! solar radiation over ice under clear sky 193 zqsrios , & ! solar radiation over ice under overcast sky 194 zcldcor , & ! cloud correction 195 zlsrise, zlsset , & ! sunrise and sunset 196 zlmunoon , & ! local noon solar altitude 197 zdlha , & ! length of the ninstr segments of the solar day 198 zps , & ! sine of latitude per sine of solar decli. 199 zpc ! cosine of latitude per cosine of solar decli. 200 REAL(wp), DIMENSION(jpi,jpj,jpl) :: & 201 zalbics , & ! albedo of ice under clear sky 202 zalbios ! albedo of ice under overcast sky 203 204 REAL(wp), DIMENSION(jpi,jpj) :: & 205 zqlw_oce , & ! long-wave heat flux over ocean 206 zqla_oce , & ! latent heat flux over ocean 207 zqsb_oce ! sensible heat flux over ocean 208 209 REAL(wp), DIMENSION(jpi,jpj,jpl) :: & 210 zqlw_ice , & ! long-wave heat flux over ice 211 zqla_ice , & ! latent heat flux over ice 212 zqsb_ice ! sensible heat flux over ice 213 214 REAL(wp), DIMENSION(jpi,jpj,jpintsr) :: & 215 zlha , & ! local hour angle 216 zalbocs , & ! tempory var. of ocean albedo under clear sky 217 zsqsro , & ! tempory var. of solar rad. over ocean 218 zsqsrics , & ! temp. var. of solar rad. over ice under clear sky 219 zsqsrios ! temp. var. of solar rad. over ice under overcast sky 220 !!--------------------------------------------------------------------- 221 222 !--------------------- 223 ! Initilization ! 224 !--------------------- 225 #if ! defined key_lim3 226 tn_ice(:,:) = psst(:,:) 227 #endif 228 229 ! Determine cloud optical depths as a function of latitude (Chou et al., 1981). 230 ! and the correction factor for taking into account the effect of clouds 231 !------------------------------------------------------ 232 IF( lbulk_init ) THEN 233 DO jj = 1, jpj 234 DO ji = 1 , jpi 235 zalat = ( 90.e0 - ABS( gphit(ji,jj) ) ) / 5.e0 236 zclat = ( 95.e0 - gphit(ji,jj) ) / 10.e0 237 indxb = 1 + INT( zalat ) 238 ! correction factor to account for the effect of clouds 239 sbudyko(ji,jj) = budyko(indxb) 240 indxc = 1 + INT( zclat ) 241 zdl = zclat - INT( zclat ) 242 zdr = 1.0 - zdl 243 stauc(ji,jj) = zdr * tauco( indxc ) + zdl * tauco( indxc + 1 ) 244 END DO 245 END DO 246 IF( nleapy == 1 ) THEN 247 yearday = 366.e0 248 ELSE IF( nleapy == 0 ) THEN 249 yearday = 365.e0 250 ELSEIF( nleapy == 30) THEN 251 yearday = 360.e0 252 ENDIF 253 lbulk_init = .FALSE. 254 ENDIF 255 256 zqlw_oce(:,:) = 0.e0 257 zqla_oce(:,:) = 0.e0 258 zqsb_oce(:,:) = 0.e0 259 zqlw_ice(:,:,:) = 0.e0 260 zqla_ice(:,:,:) = 0.e0 261 zqsb_ice(:,:,:) = 0.e0 262 263 zpis2 = rpi / 2. 264 z2pi = 2. * rpi 265 266 !CDIR NOVERRCHK 267 DO jj = 1, jpj 268 !CDIR NOVERRCHK 269 DO ji = 1, jpi 270 271 ztatm (ji,jj) = 273.15 + tatm (ji,jj) ! air temperature in Kelvins 272 zcatm1(ji,jj) = 1.0 - catm (ji,jj) ! fractional cloud cover 273 zfrld (ji,jj) = 1.0 - freeze(ji,jj) ! fractional sea ice cover 274 zpatm(ji,jj) = 101000. ! pressure 275 276 ! Computation of air density, obtained from the equation of state for dry air. 277 zrhoa(ji,jj) = zpatm(ji,jj) / ( 287.04 * ztatm(ji,jj) ) 278 279 ! zes : Saturation water vapour 280 ztamr = ztatm(ji,jj) - rtt 281 zmt1 = SIGN( 17.269, ztamr ) 282 zmt2 = SIGN( 21.875, ztamr ) 283 zmt3 = SIGN( 28.200, -ztamr ) 284 zes(ji,jj) = 611.0 * EXP ( ABS( ztamr ) * MIN ( zmt1, zmt2 ) & 285 & / ( ztatm(ji,jj) - 35.86 + MAX( zzero, zmt3 ) ) ) 286 287 ! zev : vapour pressure (hatm is relative humidity) 288 zev(ji,jj) = hatm(ji,jj) * zes(ji,jj) 289 ! square-root of vapour pressure 290 !CDIR NOVERRCHK 291 zevsqr(ji,jj) = SQRT( zev(ji,jj) * 0.01 ) 292 ! zqapb : specific humidity 293 zqatm(ji,jj) = 0.622 * zev(ji,jj) / ( zpatm(ji,jj) - 0.378 * zev(ji,jj) ) 294 295 296 !---------------------------------------------------- 297 ! Computation of snow precipitation (Ledley, 1985) | 298 !---------------------------------------------------- 299 300 zmt1 = 253.0 - ztatm(ji,jj) 301 zmt2 = ( 272.0 - ztatm(ji,jj) ) / 38.0 302 zmt3 = ( 281.0 - ztatm(ji,jj) ) / 18.0 303 zind1 = MAX( zzero, SIGN( zone, zmt1 ) ) 304 zind2 = MAX( zzero, SIGN( zone, zmt2 ) ) 305 zind3 = MAX( zzero, SIGN( zone, zmt3 ) ) 306 ! total precipitation 307 tprecip(ji,jj) = watm(ji,jj) 308 ! solid (snow) precipitation 309 sprecip(ji,jj) = tprecip(ji,jj) * & 310 & ( zind1 & 311 & + ( 1.0 - zind1 ) * ( zind2 * ( 0.5 + zmt2 ) + ( 1.0 - zind2 ) * zind3 * zmt3 ) ) 312 END DO 313 END DO 314 315 !---------------------------------------------------------- 316 ! Computation of albedo (need to calculates heat fluxes)| 317 !----------------------------------------------------------- 318 319 CALL flx_blk_albedo( zalbios, zalboos, zalbics, zalbomu ) 320 321 !------------------------------------- 322 ! Computation of solar irradiance. | 323 !---------------------------------------- 324 indaet = 1 325 ! compution of the day of the year at which the fluxes have to be calculate 326 !--The date corresponds to the middle of the time step. 327 zxday=nday_year + rdtbs2/rday 328 329 iday = INT( zxday ) 330 331 IF(ln_ctl) CALL prt_ctl_info('declin : iday ', ivar1=iday, clinfo2=' nfbulk= ', ivar2=nfbulk) 332 333 ! computation of the solar declination, his sine and his cosine 334 CALL flx_blk_declin( indaet, iday, zdecl ) 335 336 zdecl = zdecl * rad 337 zsdecl = SIN( zdecl ) 338 zcdecl = COS( zdecl ) 339 340 ! correction factor added for computation of shortwave flux to take into account the variation of 341 ! the distance between the sun and the earth during the year (Oberhuber 1988) 342 zdist = zxday * z2pi / yearday 343 zdaycor = 1.0 + 0.0013 * SIN( zdist ) + 0.0342 * COS( zdist ) 344 345 !CDIR NOVERRCHK 346 DO jj = 1, jpj 347 !CDIR NOVERRCHK 348 DO ji = 1, jpi 349 ! product of sine of latitude and sine of solar declination 350 zps (ji,jj) = SIN( gphit(ji,jj) * rad ) * zsdecl 351 ! product of cosine of latitude and cosine of solar declination 352 zpc (ji,jj) = COS( gphit(ji,jj) * rad ) * zcdecl 353 ! computation of the both local time of sunrise and sunset 354 zlsrise (ji,jj) = ACOS( - SIGN( zone, zps(ji,jj) ) * MIN( zone, SIGN( zone, zps(ji,jj) ) & 355 & * ( zps(ji,jj) / zpc(ji,jj) ) ) ) 356 zlsset (ji,jj) = - zlsrise(ji,jj) 357 ! dividing the solar day into jpintsr segments of length zdlha 358 zdlha (ji,jj) = ( zlsrise(ji,jj) - zlsset(ji,jj) ) / REAL( jpintsr ) 359 ! computation of the local noon solar altitude 360 zlmunoon(ji,jj) = ASIN ( ( zps(ji,jj) + zpc(ji,jj) ) ) / rad 361 362 ! cloud correction taken from Reed (1977) (imposed lower than 1) 363 zcldcor (ji,jj) = MIN( zone, ( 1.e0 - 0.62 * catm(ji,jj) + 0.0019 * zlmunoon(ji,jj) ) ) 364 END DO 365 END DO 366 367 ! Computation of solar heat flux at each time of the day between sunrise and sunset. 368 ! We do this to a better optimisation of the code 369 !------------------------------------------------------ 370 DO jl = 1, jpl 371 372 !CDIR NOVERRCHK 373 DO jt = 1, jpintsr 374 zcoef = FLOAT( jt ) - 0.5 375 !CDIR NOVERRCHK 376 DO jj = 1, jpj 377 !CDIR NOVERRCHK 378 DO ji = 1, jpi 379 ! local hour angle 380 zlha (ji,jj,jt) = COS ( zlsrise(ji,jj) - zcoef * zdlha(ji,jj) ) 381 382 ! cosine of local solar altitude 383 zcmue = MAX ( zzero , zps(ji,jj) + zpc(ji,jj) * zlha (ji,jj,jt) ) 384 zcmue2 = 1368.0 * zcmue * zcmue 385 zscmue = SQRT ( zcmue ) 386 zpcmue = zcmue**1.4 387 ! computation of sea-water albedo (Payne, 1972) 388 zalbocs(ji,jj,jt) = 0.05 / ( 1.1 * zpcmue + 0.15 ) 389 zalbo = zcatm1(ji,jj) * zalbocs(ji,jj,jt) + catm(ji,jj) * zalboos(ji,jj) 390 ! solar heat flux absorbed at ocean surfaces (Zillman, 1972) 391 zevo = zev(ji,jj) * 1.0e-05 392 zsqsro(ji,jj,jt) = ( 1.0 - zalbo ) * zdlha(ji,jj) * zcmue2 & 393 / ( ( zcmue + 2.7 ) * zevo + 1.085 * zcmue + 0.10 ) 394 ! solar heat flux absorbed at sea/ice surfaces 395 ! Formulation of Shine and Crane, 1984 adapted for high albedo surfaces 396 397 ! For clear sky 398 zevi = zevo 399 zalbi = zalbics(ji,jj,jl) 400 zsqsrics(ji,jj,jt) = ( 1.0 - zalbi ) * zdlha(ji,jj) * zcmue2 & 401 & / ( ( 1.0 + zcmue ) * zevi + 1.2 * zcmue + 0.0455 ) 402 403 ! For overcast sky 404 zalbi = zalbios(ji,jj,jl) 405 zsqsrios(ji,jj,jt) = zdlha(ji,jj) * & 406 & ( ( 53.5 + 1274.5 * zcmue ) * zscmue * ( 1.0 - 0.996 * zalbi ) ) & 407 & / ( 1.0 + 0.139 * stauc(ji,jj) * ( 1.0 - 0.9435 * zalbi ) ) 408 END DO 409 END DO 410 END DO 411 412 413 ! Computation of daily (between sunrise and sunset) solar heat flux absorbed 414 ! at the ocean and snow/ice surfaces. 415 !-------------------------------------------------------------------- 416 417 zalbocsd(:,:) = 0.e0 418 zqsro (:,:) = 0.e0 419 zqsrics (:,:) = 0.e0 420 zqsrios (:,:) = 0.e0 421 422 DO jt = 1, jpintsr 423 # if defined key_vectopt_loop 424 DO ji = 1, jpij 425 zalbocsd(ji,1) = zalbocsd(ji,1) + zdlha (ji,1) * zalbocs(ji,1,jt) & 426 & / MAX( 2.0 * zlsrise(ji,1) , zeps0 ) 427 zqsro (ji,1) = zqsro (ji,1) + zsqsro (ji,1,jt) 428 zqsrics (ji,1) = zqsrics (ji,1) + zsqsrics(ji,1,jt) 429 zqsrios (ji,1) = zqsrios (ji,1) + zsqsrios(ji,1,jt) 430 END DO 431 # else 432 DO jj = 1, jpj 433 DO ji = 1, jpi 434 zalbocsd(ji,jj) = zalbocsd(ji,jj) + zdlha(ji,jj) * zalbocs(ji,jj,jt) & 435 & / MAX( 2.0 * zlsrise(ji,jj) , zeps0 ) 436 zqsro (ji,jj) = zqsro (ji,jj) + zsqsro (ji,jj,jt) 437 zqsrics(ji,jj) = zqsrics (ji,jj) + zsqsrics(ji,jj,jt) 438 zqsrios(ji,jj) = zqsrios (ji,jj) + zsqsrios(ji,jj,jt) 439 END DO 440 END DO 441 # endif 442 END DO 443 444 DO jj = 1, jpj 445 DO ji = 1, jpi 446 447 !------------------------------------------- 448 ! Computation of shortwave radiation. 449 !------------------------------------------- 450 451 ! the solar heat flux absorbed at ocean and snow/ice surfaces 452 !------------------------------------------------------------ 453 454 ! For snow/ice 455 qsr_ice(ji,jj,jl) = ( zcatm1(ji,jj) * zqsrics(ji,jj) + catm(ji,jj) * zqsrios(ji,jj) ) / z2pi 456 457 ! Taking into account the ellipsity of the earth orbit 458 !----------------------------------------------------- 459 460 qsr_ice(ji,jj,jl) = qsr_ice(ji,jj,jl) * zdaycor 461 !--------------------------------------------------------------------------- 462 ! Computation of long-wave radiation ( Berliand 1952 ; all latitudes ) 463 !--------------------------------------------------------------------------- 464 465 ! tempory variables 466 ztatm3 = ztatm(ji,jj) * ztatm(ji,jj) * ztatm(ji,jj) 467 ztatm4 = ztatm3 * ztatm(ji,jj) 468 z4tatm3 = 4. * ztatm3 469 zcldeff(ji,jj) = 1.0 - sbudyko(ji,jj) * catm(ji,jj) * catm(ji,jj) 470 ztaevbk = ztatm4 * zcldeff(ji,jj) * ( 0.39 - 0.05 * zevsqr(ji,jj) ) 471 472 ! Long-Wave for Ice 473 !---------------------- 474 zqlw_ice(ji,jj,jl) = - emic * stefan * ( ztaevbk + z4tatm3 * ( t_su(ji,jj,jl) - ztatm(ji,jj) ) ) 475 476 END DO !ji 477 END DO !jj 478 479 END DO !jl 480 481 DO jj = 1, jpj 482 DO ji = 1, jpi 483 484 ! fraction of net shortwave radiation which is not absorbed in the 485 ! thin surface layer and penetrates inside the ice cover 486 ! ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 487 !------------------------------------------------------------------ 488 489 fr1_i0(ji,jj) = 0.18 * zcatm1(ji,jj) + 0.35 * catm(ji,jj) 490 fr2_i0(ji,jj) = 0.82 * zcatm1(ji,jj) + 0.65 * catm(ji,jj) 491 492 ! the solar heat flux absorbed at ocean and snow/ice surfaces 493 !------------------------------------------------------------ 494 ! For ocean 495 qsr_oce(ji,jj) = srgamma * zcldcor(ji,jj) * zqsro(ji,jj) / z2pi 496 zinda = SIGN( zone , -( -0.5 - zfrld(ji,jj) ) ) 497 zinda = 1.0 - MAX( zzero , zinda ) 498 qsr_oce(ji,jj) = ( 1.- zinda ) * qsr_oce(ji,jj) 499 500 ! Taking into account the ellipsity of the earth orbit 501 !----------------------------------------------------- 502 qsr_oce(ji,jj) = qsr_oce(ji,jj) * zdaycor 503 504 !--------------------------------------------------------------------------- 505 ! Computation of long-wave radiation ( Berliand 1952 ; all latitudes ) 506 !--------------------------------------------------------------------------- 507 508 ! tempory variables 509 ztatm3 = ztatm(ji,jj) * ztatm(ji,jj) * ztatm(ji,jj) 510 ztatm4 = ztatm3 * ztatm(ji,jj) 511 z4tatm3 = 4. * ztatm3 512 zcldeff(ji,jj) = 1.0 - sbudyko(ji,jj) * catm(ji,jj) * catm(ji,jj) 513 ztaevbk = ztatm4 * zcldeff(ji,jj) * ( 0.39 - 0.05 * zevsqr(ji,jj) ) 514 515 ! Long-Wave for Ocean 516 !----------------------- 517 zqlw_oce(ji,jj) = - emic * stefan * ( ztaevbk + z4tatm3 * ( psst (ji,jj) - ztatm(ji,jj) ) ) 518 519 END DO 520 END DO 521 522 !---------------------------------------- 523 ! Computation of turbulent heat fluxes ( Latent and sensible ) 524 !---------------------------------------- 525 !CDIR NOVERRCHK 526 DO jj = 2 , jpjm1 527 !CDIR NOVERRCHK 528 DO ji = 1, jpi 529 530 ! Turbulent heat fluxes over water 531 !---------------------------------- 532 533 ! zeso : vapour pressure at saturation of ocean 534 ! zqsato : humidity close to the ocean surface (at saturation) 535 zeso = 611.0 * EXP ( 17.2693884 * ( psst(ji,jj) - rtt ) * tmask(ji,jj,1) / ( psst(ji,jj) - 35.86 ) ) 536 zqsato = ( 0.622 * zeso ) / ( zpatm(ji,jj) - 0.378 * zeso ) 537 538 ! Drag coefficients from Large and Pond (1981,1982) 539 !--------------------------------------------------- 540 541 ! Stability parameters 542 zdteta = psst(ji,jj) - ztatm(ji,jj) 543 zdeltaq = zqatm(ji,jj) - zqsato 544 ztvmoy = ztatm(ji,jj) * ( 1. + 2.2e-3 * ztatm(ji,jj) * zqatm(ji,jj) ) 545 zdenum = MAX( vatm(ji,jj) * vatm(ji,jj) * ztvmoy, zeps ) 546 zdtetar = zdteta / zdenum 547 ztvmoyr = ztvmoy * ztvmoy * zdeltaq / zdenum 548 549 ! For stable atmospheric conditions 550 zobouks = -70.0 * 10. * ( zdtetar + 3.2e-3 * ztvmoyr ) 551 zobouks = MAX( zzero , zobouks ) 552 zpsims = -7.0 * zobouks 553 zpsihs = zpsims 554 zpsils = zpsims 555 556 ! For unstable atmospheric conditions 557 zobouku = -100.0 * 10.0 * ( zdtetar + 2.2e-3 * ztvmoyr ) 558 zobouku = MIN( zzero , zobouku ) 559 zxins = ( 1. - 16. * zobouku )**0.25 560 zlxins = LOG( ( 1. + zxins * zxins ) / 2. ) 561 zpsimu = 2. * LOG( ( 1 + zxins ) / 2. ) + zlxins - 2. * ATAN( zxins ) + zpis2 562 zpsihu = 2. * zlxins 563 zpsilu = zpsihu 564 565 ! computation of intermediate values 566 zstab = MAX( zzero , SIGN( zone , zdteta ) ) 567 zpsim = zstab * zpsimu + (1.0 - zstab ) * zpsims 568 zpsih = zstab * zpsihu + (1.0 - zstab ) * zpsihs 569 zpsil = zpsih 570 571 zvatmg = MAX( 0.032 * 1.5e-3 * vatm(ji,jj) * vatm(ji,jj) / grav, zeps ) 572 573 zcmn = vkarmn / LOG ( 10. / zvatmg ) 574 zcmn2 = zcmn * zcmn 575 zchn = 0.0327 * zcmn 576 zcln = 0.0346 * zcmn 577 zcmcmn = 1 / ( 1 - zcmn * zpsim / vkarmn ) 578 zchcm = zcmcmn / ( 1 - zchn * zpsih / ( vkarmn * zcmn ) ) 579 zclcm = zchcm 580 581 582 ! Transfer cofficient zcsho and zcleo over ocean according to Large and Pond (1981,1982) 583 !-------------------------------------------------------------- 584 zcsho = zchn * zchcm 585 zcleo = zcln * zclcm 586 587 588 ! Computation of sensible and latent fluxes over Ocean 589 !---------------------------------------------------------------- 590 591 ! computation of intermediate values 592 zrhova = zrhoa(ji,jj) * vatm(ji,jj) 593 zrhovacsho = zrhova * zcsho 594 zrhovacleo = zrhova * zcleo 595 596 ! sensible heat flux 597 zqsb_oce(ji,jj) = zrhovacsho * 1004.0 * ( psst(ji,jj) - ztatm(ji,jj) ) 598 599 ! latent heat flux 600 zqla_oce(ji,jj) = MAX(0.e0, zrhovacleo * 2.5e+06 * ( zqsato - zqatm(ji,jj) ) ) 601 602 ! Calculate evaporation over water. (kg/m2/s) 603 !------------------------------------------------- 604 evap(ji,jj) = zqla_oce(ji,jj) / cevap 605 606 END DO !ji 607 END DO !jj 608 609 DO jl = 1, jpl 610 !CDIR NOVERRCHK 611 DO jj = 2 , jpjm1 612 !CDIR NOVERRCHK 613 DO ji = 1, jpi 614 615 ! Turbulent heat fluxes over snow/ice. 616 !-------------------------------------------------- 617 618 ! zesi : vapour pressure at saturation of ice 619 ! zqsati : humidity close to the ice surface (at saturation) 620 zesi = 611.0 * EXP ( 21.8745587 * tmask(ji,jj,1) & ! tmask needed to avoid overflow in the exponential 621 & * ( t_su(ji,jj,jl) - rtt )/ ( t_su(ji,jj,jl)- 7.66 ) ) 622 zqsati = ( 0.622 * zesi ) / ( zpatm(ji,jj) - 0.378 * zesi ) 623 624 ! computation of intermediate values 625 zticemb = t_su(ji,jj,jl) - 7.66 626 zticemb2 = zticemb * zticemb 627 ztice3 = t_su(ji,jj,jl) * t_su(ji,jj,jl) * t_su(ji,jj,jl) 628 zqsati2 = zqsati * zqsati 629 zesi2 = zesi * zesi 630 zdesidt = zesi * ( 9.5 * LOG( 10.0 ) * ( rtt - 7.66 ) / zticemb2 ) 631 632 ! Transfer cofficient zcshi and zclei over ice. Assumed to be constant Parkinson 1979 ; Maykut 1982 633 !-------------------------------------------------------------------- 634 zcshi = 1.75e-03 635 zclei = zcshi 636 637 ! Computation of sensible and latent fluxes over ice 638 !---------------------------------------------------------------- 639 640 ! computation of intermediate values 641 zrhova = zrhoa(ji,jj) * vatm(ji,jj) 642 zrhovaclei = zrhova * zcshi * 2.834e+06 643 zrhovacshi = zrhova * zclei * 1004.0 644 645 ! sensible heat flux 646 zqsb_ice(ji,jj,jl) = zrhovacshi * ( t_su(ji,jj,jl) - ztatm(ji,jj) ) 647 648 ! latent heat flux 649 zqla_ice(ji,jj,jl) = zrhovaclei * ( zqsati - zqatm(ji,jj) ) 650 qla_ice (ji,jj,jl) = MAX(0.e0, zqla_ice(ji,jj,jl) ) 651 652 ! Computation of sensitivity of non solar fluxes (dQ/dT) 653 !--------------------------------------------------------------- 654 655 ! computation of long-wave, sensible and latent flux sensitivity 656 zdqlw_ice = 4.0 * emic * stefan * ztice3 657 zdqsb_ice = zrhovacshi 658 zdqla_ice = zrhovaclei * ( zdesidt * ( zqsati2 / zesi2 ) * ( zpatm(ji,jj) / 0.622 ) ) 659 660 ! total non solar sensitivity 661 dqns_ice(ji,jj,jl) = -( zdqlw_ice + zdqsb_ice + zdqla_ice ) 662 663 ! latent flux sensitivity 664 dqla_ice(ji,jj,jl) = zdqla_ice 665 666 END DO 667 END DO 668 END DO !jl 669 670 ! total non solar heat flux over ice 671 qnsr_ice(:,:,:) = zqlw_ice(:,:,:) - zqsb_ice(:,:,:) - zqla_ice(:,:,:) 672 ! total non solar heat flux over water 673 qnsr_oce(:,:) = zqlw_oce(:,:) - zqsb_oce(:,:) - zqla_oce(:,:) 674 675 ! solid precipitations ( kg/m2/day -> kg/m2/s) 676 tprecip(:,:) = tprecip (:,:) / rday 677 ! snow precipitations ( kg/m2/day -> kg/m2/s) 678 sprecip(:,:) = sprecip (:,:) / rday 679 680 CALL lbc_lnk( qsr_oce (:,:) , 'T', 1. ) 681 CALL lbc_lnk( qnsr_oce(:,:) , 'T', 1. ) 682 CALL lbc_lnk( fr1_i0 (:,:) , 'T', 1. ) 683 CALL lbc_lnk( fr2_i0 (:,:) , 'T', 1. ) 684 CALL lbc_lnk( tprecip (:,:) , 'T', 1. ) 685 CALL lbc_lnk( sprecip (:,:) , 'T', 1. ) 686 CALL lbc_lnk( evap (:,:) , 'T', 1. ) 687 DO jl = 1, jpl 688 CALL lbc_lnk( qsr_ice (:,:,jl) , 'T', 1. ) 689 CALL lbc_lnk( qnsr_ice(:,:,jl) , 'T', 1. ) 690 CALL lbc_lnk( dqns_ice(:,:,jl) , 'T', 1. ) 691 CALL lbc_lnk( qla_ice (:,:,jl) , 'T', 1. ) 692 CALL lbc_lnk( dqla_ice(:,:,jl) , 'T', 1. ) 693 END DO 694 695 qsr_oce (:,:) = qsr_oce (:,:)*tmask(:,:,1) 696 qnsr_oce(:,:) = qnsr_oce(:,:)*tmask(:,:,1) 697 DO jl = 1, jpl 698 qsr_ice (:,:,jl) = qsr_ice (:,:,jl)*tmask(:,:,1) 699 qnsr_ice(:,:,jl) = qnsr_ice(:,:,jl)*tmask(:,:,1) 700 qla_ice (:,:,jl) = qla_ice (:,:,jl)*tmask(:,:,1) 701 dqns_ice(:,:,jl) = dqns_ice(:,:,jl)*tmask(:,:,1) 702 dqla_ice(:,:,jl) = dqla_ice(:,:,jl)*tmask(:,:,1) 703 END DO 704 fr1_i0 (:,:) = fr1_i0 (:,:)*tmask(:,:,1) 705 fr2_i0 (:,:) = fr2_i0 (:,:)*tmask(:,:,1) 706 tprecip (:,:) = tprecip (:,:)*tmask(:,:,1) 707 sprecip (:,:) = sprecip (:,:)*tmask(:,:,1) 708 evap (:,:) = evap (:,:)*tmask(:,:,1) 709 710 711 END SUBROUTINE flx_blk 712 713 714 #elif defined key_lim2 81 715 82 716 SUBROUTINE flx_blk( psst ) … … 212 846 !!--------------------------------------------------------------------- 213 847 214 !---------------------215 ! Initilization !216 !---------------------217 #if ! defined key_lim3218 tn_ice(:,:) = psst(:,:)219 #endif220 221 848 ! Determine cloud optical depths as a function of latitude (Chou et al., 1981). 222 849 ! and the correction factor for taking into account the effect of clouds … … 680 1307 681 1308 END SUBROUTINE flx_blk 1309 #endif 682 1310 683 1311
Note: See TracChangeset
for help on using the changeset viewer.
