Changeset 11804 for NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90
- Timestamp:
- 2019-10-25T17:25:19+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90
r11772 r11804 520 520 !! Estimate of potential temperature at z=rn_zqt, based on adiabatic lapse-rate 521 521 !! (see Josey, Gulev & Yu, 2013) / doi=10.1016/B978-0-12-391851-2.00005-2 522 !! (since reanalysis products provide Tat z, not theta !)522 !! (since reanalysis products provide absolute temperature "T" at z, not theta !) 523 523 ztpot = sf(jp_tair)%fnow(:,:,1) + gamma_moist( sf(jp_tair)%fnow(:,:,1), zqair(:,:) ) * rn_zqt 524 524 … … 594 594 END IF ! IF ( ln_skin_cs .OR. ln_skin_wl ) 595 595 596 597 598 599 ! ! Compute true air density :600 IF( ABS(rn_zu - rn_zqt) > 0.01 ) THEN ! At zu: (probably useless to remove rhoa*grav*rn_zu from SLP...)601 rhoa(:,:) = rho_air( t_zu(:,:) , q_zu(:,:) , sf(jp_slp)%fnow(:,:,1) )602 ELSE ! At zt:603 rhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), zqair(:,:), sf(jp_slp)%fnow(:,:,1) )604 END IF605 606 596 !! CALL iom_put( "Cd_oce", Cd_atm) ! output value of pure ocean-atm. transfer coef. 607 597 !! CALL iom_put( "Ch_oce", Ch_atm) ! output value of pure ocean-atm. transfer coef. 608 609 DO jj = 1, jpj ! tau module, i and j component 598 599 IF( ABS(rn_zu - rn_zqt) < 0.1_wp ) THEN 600 !! If zu == zt, then ensuring once for all that: 601 t_zu(:,:) = ztpot(:,:) 602 q_zu(:,:) = zqair(:,:) 603 END IF 604 605 606 ! Turbulent fluxes over ocean => TURB_FLUXES @ sbcblk_phy.F90 607 ! ------------------------------------------------------------- 608 609 CALL TURB_FLUXES( rn_zu, zst(:,:), zsq(:,:), t_zu(:,:), q_zu(:,:), Cd_atm(:,:), Ch_atm(:,:), Ce_atm(:,:), & 610 & wndm(:,:), zU_zu(:,:), sf(jp_slp)%fnow(:,:,1), & 611 & taum(:,:), zqsb(:,:), zqla(:,:), & 612 & pEvap=zevap(:,:), prhoa=rhoa(:,:) ) 613 614 zqla(:,:) = zqla(:,:) * tmask(:,:,1) 615 zqsb(:,:) = zqsb(:,:) * tmask(:,:,1) 616 taum(:,:) = taum(:,:) * tmask(:,:,1) 617 zevap(:,:) = zevap(:,:) * tmask(:,:,1) 618 619 DO jj = 1, jpj ! tau i and j component on T-grid points 610 620 DO ji = 1, jpi 611 zztmp = rhoa(ji,jj) * zU_zu(ji,jj) * Cd_atm(ji,jj) ! using bulk wind speed 612 taum (ji,jj) = zztmp * wndm (ji,jj) 621 zztmp = taum(ji,jj) / wndm(ji,jj) 613 622 zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) 614 623 zwnd_j(ji,jj) = zztmp * zwnd_j(ji,jj) 615 624 END DO 616 625 END DO 617 618 626 ! ! add the HF tau contribution to the wind stress module 619 627 IF( lhftau ) taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1) … … 634 642 CALL lbc_lnk_multi( 'sbcblk', utau, 'U', -1., vtau, 'V', -1. ) 635 643 636 ! Turbulent fluxes over ocean637 ! -----------------------------638 639 ! zqla used as temporary array, for rho*U (common term of bulk formulae):640 zqla(:,:) = rhoa(:,:) * zU_zu(:,:) * tmask(:,:,1)641 642 IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN643 !! q_air and t_air are given at 10m (wind reference height)644 zevap(:,:) = rn_efac*MAX( 0._wp, zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - zqair(:,:)) ) ! Evaporation, using bulk wind speed645 zqsb (:,:) = cp_air(zqair(:,:)) * zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - ztpot(:,:) ) ! Sensible Heat, using bulk wind speed646 ELSE647 !! q_air and t_air are not given at 10m (wind reference height)648 ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!!649 zevap(:,:) = rn_efac*MAX( 0._wp, zqla(:,:)*Ce_atm(:,:)*(zsq(:,:) - q_zu(:,:)) ) ! Evaporation, using bulk wind speed650 zqsb (:,:) = cp_air(zqair(:,:)) * zqla(:,:)*Ch_atm(:,:)*(zst(:,:) - t_zu(:,:) ) ! Sensible Heat, using bulk wind speed651 ENDIF652 653 zqla(:,:) = L_vap(zst(:,:)) * zevap(:,:) ! Latent Heat flux654 655 644 656 645 ! ----------------------------------------------------------------------------- ! … … 658 647 ! ----------------------------------------------------------------------------- ! 659 648 660 !! LB: now after Turbulent fluxes because must use the skin temperature rather that the SST ! (zst is skin temperature if ln_skin_cs==.TRUE. .OR. ln_skin_wl==.TRUE.)649 !! LB: now moved after Turbulent fluxes because must use the skin temperature rather that the SST ! (zst is skin temperature if ln_skin_cs==.TRUE. .OR. ln_skin_wl==.TRUE.) 661 650 zqlw(:,:) = emiss_w * ( sf(jp_qlw)%fnow(:,:,1) - stefan*zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:) ) * tmask(:,:,1) ! Net radiative longwave flux 662 651 … … 681 670 & - sf(jp_prec)%fnow(:,:,1) * rn_pfac ) * tmask(:,:,1) 682 671 ! 683 qns(:,:) = zqlw(:,:) - zqsb(:,:) -zqla(:,:) & ! Downward Non Solar672 qns(:,:) = zqlw(:,:) + zqsb(:,:) + zqla(:,:) & ! Downward Non Solar 684 673 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * rLfus & ! remove latent melting heat for solid precip 685 674 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST !LB??? pst is Celsius !? … … 691 680 ! 692 681 #if defined key_si3 693 qns_oce(:,:) = zqlw(:,:) - zqsb(:,:) -zqla(:,:) ! non solar without emp (only needed by SI3)682 qns_oce(:,:) = zqlw(:,:) + zqsb(:,:) + zqla(:,:) ! non solar without emp (only needed by SI3) 694 683 qsr_oce(:,:) = qsr(:,:) 695 684 #endif … … 698 687 CALL iom_put( "rho_air" , rhoa ) ! output air density (kg/m^3) !#LB 699 688 CALL iom_put( "qlw_oce" , zqlw ) ! output downward longwave heat over the ocean 700 CALL iom_put( "qsb_oce" , -zqsb ) ! output downward sensible heat over the ocean701 CALL iom_put( "qla_oce" , -zqla ) ! output downward latent heat over the ocean689 CALL iom_put( "qsb_oce" , zqsb ) ! output downward sensible heat over the ocean 690 CALL iom_put( "qla_oce" , zqla ) ! output downward latent heat over the ocean 702 691 CALL iom_put( "evap_oce" , zevap ) ! evaporation 703 CALL iom_put( "qemp_oce" , qns-zqlw +zqsb+zqla ) ! output downward heat content of E-P over the ocean692 CALL iom_put( "qemp_oce" , qns-zqlw-zqsb-zqla ) ! output downward heat content of E-P over the ocean 704 693 CALL iom_put( "qns_oce" , qns ) ! output downward non solar heat over the ocean 705 694 CALL iom_put( "qsr_oce" , qsr ) ! output downward solar heat over the ocean
Note: See TracChangeset
for help on using the changeset viewer.