- Timestamp:
- 2020-03-28T08:38:26+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r12581_ticket2418/src/OCE/SBC/sbcblk_phy.F90
r12377 r12623 520 520 zCe = zz0*pqst(ji,jj)/zdq 521 521 522 CALL BULK_FORMULA ( pzu, pTs(ji,jj), pqs(ji,jj), pTa(ji,jj), pqa(ji,jj), &523 & zCd, zCh, zCe,&524 & pwnd(ji,jj), pUb(ji,jj), pslp(ji,jj),&525 & pTau(ji,jj), zQsen, zQlat )526 522 CALL BULK_FORMULA_SCLR( pzu, pTs(ji,jj), pqs(ji,jj), pTa(ji,jj), pqa(ji,jj), & 523 & zCd, zCh, zCe, & 524 & pwnd(ji,jj), pUb(ji,jj), pslp(ji,jj), & 525 & pTau(ji,jj), zQsen, zQlat ) 526 527 527 zTs2 = pTs(ji,jj)*pTs(ji,jj) 528 528 zQlw = emiss_w*(prlw(ji,jj) - stefan*zTs2*zTs2) ! Net longwave flux … … 535 535 536 536 537 SUBROUTINE BULK_FORMULA_VCTR( pzu, pTs, pqs, pTa, pqa, & 538 & pCd, pCh, pCe, & 539 & pwnd, pUb, pslp, & 540 & pTau, pQsen, pQlat, pEvap, prhoa ) 537 SUBROUTINE BULK_FORMULA_SCLR( pzu, pTs, pqs, pTa, pqa, & 538 & pCd, pCh, pCe, & 539 & pwnd, pUb, pslp, & 540 & pTau, pQsen, pQlat, & 541 & pEvap, prhoa, pfact_evap ) 542 !!---------------------------------------------------------------------------------- 543 REAL(wp), INTENT(in) :: pzu ! height above the sea-level where all this takes place (normally 10m) 544 REAL(wp), INTENT(in) :: pTs ! water temperature at the air-sea interface [K] 545 REAL(wp), INTENT(in) :: pqs ! satur. spec. hum. at T=pTs [kg/kg] 546 REAL(wp), INTENT(in) :: pTa ! potential air temperature at z=pzu [K] 547 REAL(wp), INTENT(in) :: pqa ! specific humidity at z=pzu [kg/kg] 548 REAL(wp), INTENT(in) :: pCd 549 REAL(wp), INTENT(in) :: pCh 550 REAL(wp), INTENT(in) :: pCe 551 REAL(wp), INTENT(in) :: pwnd ! wind speed module at z=pzu [m/s] 552 REAL(wp), INTENT(in) :: pUb ! bulk wind speed at z=pzu (inc. pot. effect of gustiness etc) [m/s] 553 REAL(wp), INTENT(in) :: pslp ! sea-level atmospheric pressure [Pa] 554 !! 555 REAL(wp), INTENT(out) :: pTau ! module of the wind stress [N/m^2] 556 REAL(wp), INTENT(out) :: pQsen ! [W/m^2] 557 REAL(wp), INTENT(out) :: pQlat ! [W/m^2] 558 !! 559 REAL(wp), INTENT(out), OPTIONAL :: pEvap ! Evaporation [kg/m^2/s] 560 REAL(wp), INTENT(out), OPTIONAL :: prhoa ! Air density at z=pzu [kg/m^3] 561 REAL(wp), INTENT(in) , OPTIONAL :: pfact_evap ! ABOMINATION: corrective factor for evaporation (doing this against my will! /laurent) 562 !! 563 REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap, zfact_evap 564 INTEGER :: jq 565 !!---------------------------------------------------------------------------------- 566 zfact_evap = 1._wp 567 IF( PRESENT(pfact_evap) ) zfact_evap = pfact_evap 568 569 !! Need ztaa, absolute temperature at pzu (formula to estimate rho_air needs absolute temperature, not the potential temperature "pTa") 570 ztaa = pTa ! first guess... 571 DO jq = 1, 4 572 zgamma = gamma_moist( 0.5*(ztaa+pTs) , pqa ) !LOLO: why not "0.5*(pqs+pqa)" rather then "pqa" ??? 573 ztaa = pTa - zgamma*pzu ! Absolute temp. is slightly colder... 574 END DO 575 zrho = rho_air(ztaa, pqa, pslp) 576 zrho = rho_air(ztaa, pqa, pslp-zrho*grav*pzu) ! taking into account that we are pzu m above the sea level where SLP is given! 577 578 zUrho = pUb*MAX(zrho, 1._wp) ! rho*U10 579 580 pTau = zUrho * pCd * pwnd ! Wind stress module 581 582 zevap = zUrho * pCe * (pqa - pqs) 583 pQsen = zUrho * pCh * (pTa - pTs) * cp_air(pqa) 584 pQlat = L_vap(pTs) * zevap 585 586 IF( PRESENT(pEvap) ) pEvap = - zfact_evap * zevap 587 IF( PRESENT(prhoa) ) prhoa = zrho 588 589 END SUBROUTINE BULK_FORMULA_SCLR 590 591 SUBROUTINE BULK_FORMULA_VCTR( pzu, pTs, pqs, pTa, pqa, & 592 & pCd, pCh, pCe, & 593 & pwnd, pUb, pslp, & 594 & pTau, pQsen, pQlat, & 595 & pEvap, prhoa, pfact_evap ) 541 596 !!---------------------------------------------------------------------------------- 542 597 REAL(wp), INTENT(in) :: pzu ! height above the sea-level where all this takes place (normally 10m) … … 558 613 REAL(wp), DIMENSION(jpi,jpj), INTENT(out), OPTIONAL :: pEvap ! Evaporation [kg/m^2/s] 559 614 REAL(wp), DIMENSION(jpi,jpj), INTENT(out), OPTIONAL :: prhoa ! Air density at z=pzu [kg/m^3] 560 !! 561 REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap 562 INTEGER :: ji, jj, jq ! dummy loop indices 563 !!---------------------------------------------------------------------------------- 564 DO_2D_11_11 565 566 !! Need ztaa, absolute temperature at pzu (formula to estimate rho_air needs absolute temperature, not the potential temperature "pTa") 567 ztaa = pTa(ji,jj) ! first guess... 568 DO jq = 1, 4 569 zgamma = gamma_moist( 0.5*(ztaa+pTs(ji,jj)) , pqa(ji,jj) ) 570 ztaa = pTa(ji,jj) - zgamma*pzu ! Absolute temp. is slightly colder... 571 END DO 572 zrho = rho_air(ztaa, pqa(ji,jj), pslp(ji,jj)) 573 zrho = rho_air(ztaa, pqa(ji,jj), pslp(ji,jj)-zrho*grav*pzu) ! taking into account that we are pzu m above the sea level where SLP is given! 574 575 zUrho = pUb(ji,jj)*MAX(zrho, 1._wp) ! rho*U10 576 577 pTau(ji,jj) = zUrho * pCd(ji,jj) * pwnd(ji,jj) ! Wind stress module 578 579 zevap = zUrho * pCe(ji,jj) * (pqa(ji,jj) - pqs(ji,jj)) 580 pQsen(ji,jj) = zUrho * pCh(ji,jj) * (pTa(ji,jj) - pTs(ji,jj)) * cp_air(pqa(ji,jj)) 581 pQlat(ji,jj) = L_vap(pTs(ji,jj)) * zevap 582 583 IF( PRESENT(pEvap) ) pEvap(ji,jj) = - zevap 615 REAL(wp), INTENT(in) , OPTIONAL :: pfact_evap ! ABOMINATION: corrective factor for evaporation (doing this against my will! /laurent) 616 !! 617 REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap, zfact_evap 618 INTEGER :: ji, jj 619 !!---------------------------------------------------------------------------------- 620 zfact_evap = 1._wp 621 IF( PRESENT(pfact_evap) ) zfact_evap = pfact_evap 622 623 DO_2D_11_11 624 625 CALL BULK_FORMULA_SCLR( pzu, pTs(ji,jj), pqs(ji,jj), pTa(ji,jj), pqa(ji,jj), & 626 & pCd(ji,jj), pCh(ji,jj), pCe(ji,jj), & 627 & pwnd(ji,jj), pUb(ji,jj), pslp(ji,jj), & 628 & pTau(ji,jj), pQsen(ji,jj), pQlat(ji,jj), & 629 & pEvap=zevap, prhoa=zrho, pfact_evap=zfact_evap ) 630 631 IF( PRESENT(pEvap) ) pEvap(ji,jj) = zevap 584 632 IF( PRESENT(prhoa) ) prhoa(ji,jj) = zrho 585 633 586 634 END_2D 587 635 END SUBROUTINE BULK_FORMULA_VCTR 588 589 590 SUBROUTINE BULK_FORMULA_SCLR( pzu, pTs, pqs, pTa, pqa, &591 & pCd, pCh, pCe, &592 & pwnd, pUb, pslp, &593 & pTau, pQsen, pQlat, pEvap, prhoa )594 !!----------------------------------------------------------------------------------595 REAL(wp), INTENT(in) :: pzu ! height above the sea-level where all this takes place (normally 10m)596 REAL(wp), INTENT(in) :: pTs ! water temperature at the air-sea interface [K]597 REAL(wp), INTENT(in) :: pqs ! satur. spec. hum. at T=pTs [kg/kg]598 REAL(wp), INTENT(in) :: pTa ! potential air temperature at z=pzu [K]599 REAL(wp), INTENT(in) :: pqa ! specific humidity at z=pzu [kg/kg]600 REAL(wp), INTENT(in) :: pCd601 REAL(wp), INTENT(in) :: pCh602 REAL(wp), INTENT(in) :: pCe603 REAL(wp), INTENT(in) :: pwnd ! wind speed module at z=pzu [m/s]604 REAL(wp), INTENT(in) :: pUb ! bulk wind speed at z=pzu (inc. pot. effect of gustiness etc) [m/s]605 REAL(wp), INTENT(in) :: pslp ! sea-level atmospheric pressure [Pa]606 !!607 REAL(wp), INTENT(out) :: pTau ! module of the wind stress [N/m^2]608 REAL(wp), INTENT(out) :: pQsen ! [W/m^2]609 REAL(wp), INTENT(out) :: pQlat ! [W/m^2]610 !!611 REAL(wp), INTENT(out), OPTIONAL :: pEvap ! Evaporation [kg/m^2/s]612 REAL(wp), INTENT(out), OPTIONAL :: prhoa ! Air density at z=pzu [kg/m^3]613 !!614 REAL(wp) :: ztaa, zgamma, zrho, zUrho, zevap615 INTEGER :: jq616 !!----------------------------------------------------------------------------------617 618 !! Need ztaa, absolute temperature at pzu (formula to estimate rho_air needs absolute temperature, not the potential temperature "pTa")619 ztaa = pTa ! first guess...620 DO jq = 1, 4621 zgamma = gamma_moist( 0.5*(ztaa+pTs) , pqa )622 ztaa = pTa - zgamma*pzu ! Absolute temp. is slightly colder...623 END DO624 zrho = rho_air(ztaa, pqa, pslp)625 zrho = rho_air(ztaa, pqa, pslp-zrho*grav*pzu) ! taking into account that we are pzu m above the sea level where SLP is given!626 627 zUrho = pUb*MAX(zrho, 1._wp) ! rho*U10628 629 pTau = zUrho * pCd * pwnd ! Wind stress module630 631 zevap = zUrho * pCe * (pqa - pqs)632 pQsen = zUrho * pCh * (pTa - pTs) * cp_air(pqa)633 pQlat = L_vap(pTs) * zevap634 635 IF( PRESENT(pEvap) ) pEvap = - zevap636 IF( PRESENT(prhoa) ) prhoa = zrho637 638 END SUBROUTINE BULK_FORMULA_SCLR639 640 641 636 642 637
Note: See TracChangeset
for help on using the changeset viewer.