SUBROUTINE coefcdrag(klon, knon, nsrf, zxli, & speed, t, q, zgeop, psol, & ts, qsurf, rugos, okri, ri1, & cdram, cdrah, cdran, zri1, pref) USE indice_sol_mod USE lmdz_abort_physic, ONLY: abort_physic USE lmdz_clesphys USE lmdz_yoethf USE lmdz_yomcst IMPLICIT NONE !------------------------------------------------------------------------- ! Objet : calcul des cdrags pour le moment (cdram) et les flux de chaleur ! sensible et latente (cdrah), du cdrag neutre (cdran), ! du nombre de Richardson entre la surface et le niveau de reference ! (zri1) et de la pression au niveau de reference (pref). ! I. Musat, 01.07.2002 !------------------------------------------------------------------------- ! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude) ! knon----input-I- nombre de points pour un type de surface ! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90 ! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li ! speed---input-R- module du vent au 1er niveau du modele ! t-------input-R- temperature de l'air au 1er niveau du modele ! q-------input-R- humidite de l'air au 1er niveau du modele ! zgeop---input-R- geopotentiel au 1er niveau du modele ! psol----input-R- pression au sol ! ts------input-R- temperature de l'air a la surface ! qsurf---input-R- humidite de l'air a la surface ! rugos---input-R- rugosite ! okri----input-L- TRUE si on veut tester le nb. Richardson entre la sfce ! et zref par rapport au Ri entre la sfce et la 1ere couche ! ri1-----input-R- nb. Richardson entre la surface et la 1ere couche ! cdram--output-R- cdrag pour le moment ! cdrah--output-R- cdrag pour les flux de chaleur latente et sensible ! cdran--output-R- cdrag neutre ! zri1---output-R- nb. Richardson entre la surface et la couche zgeop/RG ! pref---output-R- pression au niveau zgeop/RG INTEGER, INTENT(IN) :: klon, knon, nsrf LOGICAL, INTENT(IN) :: zxli REAL, DIMENSION(klon), INTENT(IN) :: speed, t, q, zgeop, psol REAL, DIMENSION(klon), INTENT(IN) :: ts, qsurf, rugos, ri1 LOGICAL, INTENT(IN) :: okri REAL, DIMENSION(klon), INTENT(OUT) :: cdram, cdrah, cdran, zri1, pref !------------------------------------------------------------------------- ! Quelques constantes : REAL, parameter :: RKAR = 0.40, CB = 5.0, CC = 5.0, CD = 5.0, cepdu2 = (0.1)**2 ! Variables locales : INTEGER :: i REAL, DIMENSION(klon) :: zdu2, zdphi, ztsolv, ztvd REAL, DIMENSION(klon) :: zscf, friv, frih, zucf, zcr REAL, DIMENSION(klon) :: zcfm1, zcfh1 REAL, DIMENSION(klon) :: zcfm2, zcfh2 REAL, DIMENSION(klon) :: trm0, trm1 CHARACTER (LEN = 80) :: abort_message CHARACTER (LEN = 20) :: modname = 'coefcdra' !------------------------------------------------------------------------- REAL :: fsta, fins, x fsta(x) = 1.0 / (1.0 + 10.0 * x * (1 + 8.0 * x)) fins(x) = SQRT(1.0 - 18.0 * x) !------------------------------------------------------------------------- abort_message = 'obsolete, remplace par cdrag, use at you own risk' CALL abort_physic(modname, abort_message, 1) DO i = 1, knon zdphi(i) = zgeop(i) zdu2(i) = max(cepdu2, speed(i)**2) pref(i) = exp(log(psol(i)) - zdphi(i) / (RD * t(i) * & (1. + RETV * max(q(i), 0.0)))) ztsolv(i) = ts(i) ! ztvd(i) = t(i) * (psol(i)/pref(i))**RKAPPA ! ztvd(i) = (t(i)+zdphi(i)/RCPD/(1.+RVTMP2*q(i))) & ! *(1.+RETV*q(i)) ztvd(i) = (t(i) + zdphi(i) / RCPD / (1. + RVTMP2 * q(i))) trm0(i) = 1. + RETV * max(qsurf(i), 0.0) trm1(i) = 1. + RETV * max(q(i), 0.0) ztsolv(i) = ztsolv(i) * trm0(i) ztvd(i) = ztvd(i) * trm1(i) zri1(i) = zdphi(i) * (ztvd(i) - ztsolv(i)) / (zdu2(i) * ztvd(i)) ! on teste zri1 par rapport au Richardson de la 1ere couche ri1 !IM +++ IF(1==0) THEN IF (okri) THEN IF (ri1(i)>=0.0.AND.zri1(i)<0.0) THEN zri1(i) = ri1(i) ELSE IF(ri1(i)<0.0.AND.zri1(i)>=0.0) THEN zri1(i) = ri1(i) ENDIF ENDIF ENDIF !IM --- cdran(i) = (RKAR / log(1. + zdphi(i) / (RG * rugos(i))))**2 IF (zri1(i) >= 0.) THEN ! situation stable : pour eviter les inconsistances dans les cas ! tres stables on limite zri1 a 20. cf Hess et al. (1995) zri1(i) = min(20., zri1(i)) IF (.NOT.zxli) THEN zscf(i) = SQRT(1. + CD * ABS(zri1(i))) friv(i) = max(1. / (1. + 2. * CB * zri1(i) / zscf(i)), f_ri_cd_min) zcfm1(i) = cdran(i) * friv(i) frih(i) = max(1. / (1. + 3. * CB * zri1(i) * zscf(i)), f_ri_cd_min) ! zcfh1(i) = cdran(i) * frih(i) zcfh1(i) = f_cdrag_ter * cdran(i) * frih(i) IF(nsrf==is_oce) zcfh1(i) = f_cdrag_oce * cdran(i) * frih(i) cdram(i) = zcfm1(i) cdrah(i) = zcfh1(i) ELSE cdram(i) = cdran(i) * fsta(zri1(i)) cdrah(i) = cdran(i) * fsta(zri1(i)) ENDIF ELSE ! situation instable IF (.NOT.zxli) THEN zucf(i) = 1. / (1. + 3.0 * CB * CC * cdran(i) * SQRT(ABS(zri1(i)) & * (1.0 + zdphi(i) / (RG * rugos(i))))) zcfm2(i) = cdran(i) * max((1. - 2.0 * CB * zri1(i) * zucf(i)), f_ri_cd_min) ! zcfh2(i) = cdran(i)*max((1.-3.0*CB*zri1(i)*zucf(i)),f_ri_cd_min) zcfh2(i) = f_cdrag_ter * cdran(i) * max((1. - 3.0 * CB * zri1(i) * zucf(i)), f_ri_cd_min) cdram(i) = zcfm2(i) cdrah(i) = zcfh2(i) ELSE cdram(i) = cdran(i) * fins(zri1(i)) cdrah(i) = cdran(i) * fins(zri1(i)) ENDIF ! cdrah sur l'ocean cf. Miller et al. (1992) zcr(i) = (0.0016 / (cdran(i) * SQRT(zdu2(i)))) * ABS(ztvd(i) - ztsolv(i)) & **(1. / 3.) ! IF (nsrf.EQ.is_oce) cdrah(i) = cdran(i)*(1.0+zcr(i)**1.25) & ! **(1./1.25) IF (nsrf==is_oce) cdrah(i) = f_cdrag_oce * cdran(i) * (1.0 + zcr(i)**1.25) & **(1. / 1.25) ENDIF END DO END SUBROUTINE coefcdrag