Changeset 407 for LMDZ.3.3/branches/rel-LF/libf/phylmd
- Timestamp:
- Oct 15, 2002, 4:59:04 PM (22 years ago)
- File:
-
- 1 edited
-
LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F (modified) (14 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F
r400 r407 16 16 . flux_t,flux_q,flux_u,flux_v,cdragh,cdragm, 17 17 . dflux_t,dflux_q, 18 . zcoefh,zu1,zv1 )18 . zcoefh,zu1,zv1, t2m, q2m, u10m, v10m) 19 19 cAA . itr, tr, flux_surf, d_tr) 20 20 cAA REM: … … 35 35 USE ioipsl 36 36 USE interface_surf 37 USE stdlevvar_int 37 38 IMPLICIT none 38 39 c====================================================================== … … 169 170 c 170 171 #include "YOMCST.h" 172 #include "YOETHF.h" 173 #include "FCTTRE.h" 171 174 REAL u1lay(klon), v1lay(klon) 172 175 REAL delp(klon,klev) … … 203 206 integer idayref 204 207 #include "temps.h" 205 208 REAL t2m(klon,nbsrf), q2m(klon,nbsrf) 209 REAL u10m(klon,nbsrf), v10m(klon,nbsrf) 210 c 211 REAL yt2m(klon), yq2m(klon), yu10m(klon) 212 c 213 REAL uzon(klon), vmer(klon) 214 REAL tair1(klon), qair1(klon), tairsol(klon) 215 REAL psfce(klon), patm(klon) 216 c 217 REAL qairsol(klon), zgeo1(klon) 218 REAL rugo1(klon) 219 c 220 LOGICAL zxli ! utiliser un jeu de fonctions simples 221 PARAMETER (zxli=.FALSE.) 222 c 223 REAL zt, zqs, zdelta, zcor 224 REAL t_coup 225 PARAMETER(t_coup=273.15) 226 C 227 PRINT*,'IMclmain klon=',klon 206 228 IF (first_appel) THEN 207 229 first_appel=.false. … … 276 298 yv1 = 0.0 277 299 yrads = 0.0 278 ypaprs = 0.0279 300 ypaprs = 0.0 280 301 ypplay = 0.0 … … 573 594 ENDDO 574 595 c 596 c 597 #undef T2m 598 #ifdef T2m 599 ccc diagnostic t,q a 2m et u, v a 10m 600 c 601 DO j=1, knon 602 i = ni(j) 603 uzon(j) = yu(j,1) + y_d_u(j,1) 604 vmer(j) = yv(j,1) + y_d_v(j,1) 605 tair1(j) = yt(j,1) + y_d_t(j,1) 606 qair1(j) = yq(j,1) + y_d_q(j,1) 607 zgeo1(j) = RD * tair1(j) / (0.5*(ypaprs(j,1)+ypplay(j,1))) 608 & * (ypaprs(j,1)-ypplay(j,1)) 609 tairsol(j) = yts(j) + y_d_ts(j) 610 rugo1(j) = yrugos(j) 611 IF(nsrf.EQ.is_oce) THEN 612 rugo1(j) = rugos(i,nsrf) 613 ENDIF 614 psfce(j)=ypaprs(j,1) 615 patm(j)=ypplay(j,1) 616 c 617 IF (nsrf.EQ.1) THEN 618 qairsol(j) = yqsol(j) 619 ELSE IF(nsrf.GT.1) THEN 620 zt = ts(i,nsrf) 621 IF (thermcep) THEN 622 zdelta = MAX(0.,SIGN(1.,RTT-zt)) 623 zqs = R2ES * FOEEW(zt,zdelta) / ypplay(j,1) 624 zqs = MIN(0.5,zqs) 625 zcor = 1./(1.-RETV*zqs) 626 zqs = zqs*zcor 627 ELSE 628 IF (zt .LT. t_coup) THEN 629 zqs = qsats(zt) / ypplay(j,1) 630 ELSE 631 zqs = qsatl(zt) / ypplay(j,1) 632 ENDIF 633 ENDIF 634 qairsol(j) = zqs 635 ENDIF 636 ENDDO 637 c 638 IF(nsrf.EQ.3) THEN 639 j=1465 640 WRITE(*,*)' INstO',klon,knon,nsrf,zxli,uzon(j),vmer(j), 641 & tair1(j),qair1(j),zgeo1(j),tairsol(j),qairsol(j),rugo1(j), 642 & psfce(j),patm(j) 643 ENDIF 644 c 645 CALL stdlevvar(klon, knon, nsrf, zxli, 646 & uzon, vmer, tair1, qair1, zgeo1, 647 & tairsol, qairsol, rugo1, psfce, patm, 648 & yt2m, yq2m, yu10m) 649 650 c 651 IF(nsrf.EQ.3) THEN 652 j=1465 653 WRITE(*,*)' OUstd',klon,knon,nsrf,zxli,uzon(j),vmer(j), 654 & tair1(j),qair1(j),zgeo1(j),tairsol(j),qairsol(j),rugo1(j), 655 & psfce(j),patm(j) 656 WRITE(*,*)' tqu',yt2m(j),yq2m(j),yu10m(j) 657 ENDIF 658 c 659 DO j=1, knon 660 i = ni(j) 661 t2m(i,nsrf)=yt2m(j) 662 663 IF(nsrf.EQ.3) THEN 664 IF(j.EQ.1465) THEN 665 WRITE(*,*) 't2m APRES stdlev',j,i,tair1(j),t2m(i,nsrf), 666 $ tairsol(j),rlon(i),rlat(i) 667 ENDIF 668 ENDIF 669 c 670 q2m(i,nsrf)=yq2m(j) 671 c 672 c u10m, v10m : composantes du vent a 10m sans spirale de Ekman 673 u10m(i,nsrf)=(yu10m(j) * uzon(j))/sqrt(uzon(j)**2+vmer(j)**2) 674 v10m(i,nsrf)=(yu10m(j) * vmer(j))/sqrt(uzon(j)**2+vmer(j)**2) 675 c 676 ENDDO 677 #else 678 DO j=1, knon 679 t2m(i,nsrf)=0. 680 q2m(i,nsrf)=0. 681 u10m(i,nsrf)=0. 682 v10m(i,nsrf)=0. 683 ENDDO 684 #endif 575 685 99999 CONTINUE 576 686 c … … 1028 1138 . u,v,t,q, 1029 1139 . pcfm, pcfh) 1140 USE clcdrag_int 1030 1141 IMPLICIT none 1031 1142 c====================================================================== … … 1070 1181 REAL cepdu2, ckap, cb, cc, cd, clam 1071 1182 PARAMETER (cepdu2 =(0.1)**2) 1072 PARAMETER ( ckap=0.35)1183 PARAMETER (CKAP=0.4) 1073 1184 PARAMETER (cb=5.0) 1074 1185 PARAMETER (cc=5.0) … … 1103 1214 REAL zri(klon) 1104 1215 REAL zl2(klon) 1105 REAL zcfm1(klon), zcfm2(klon) 1106 REAL zcfh1(klon), zcfh2(klon) 1107 REAL zdphi, zdu2, ztvd, ztvu, ztsolv, zcdn 1108 REAL zscf, zucf, zcr 1216 1217 REAL u1(klon), v1(klon), t1(klon), q1(klon), z1(klon) 1218 REAL pcfm1(klon), pcfh1(klon) 1219 c 1220 REAL zdphi, zdu2, ztvd, ztvu, zcdn 1221 REAL zscf 1109 1222 REAL zt, zq, zdelta, zcvm5, zcor, zqs, zfr, zdqs 1110 1223 REAL z2geomf, zalh2, zalm2, zscfh, zscfm … … 1116 1229 c essai qsurf 1117 1230 real qsurf(klon) 1118 real friv, frih1119 1231 c 1120 1232 LOGICAL appel1er … … 1188 1300 c Calculer le frottement au sol (Cdrag) 1189 1301 c 1190 CALL clcdrag(knon, nsrf, zxli, u, v, t, q, zgeop, 1191 . ts, qsurf, rugos, pcfm, pcfh, zcdn, zri) 1302 c DO i = 1, knon 1303 DO i = 1, klon 1304 u1(i) = u(i,1) 1305 v1(i) = v(i,1) 1306 t1(i) = t(i,1) 1307 q1(i) = q(i,1) 1308 z1(i) = zgeop(i,1) 1309 ENDDO 1310 c 1311 CALL clcdrag(klon, knon, nsrf, zxli, 1312 $ u1, v1, t1, q1, z1, 1313 $ ts, qsurf, rugos, 1314 $ pcfm1, pcfh1) 1315 C 1316 DO i = 1, knon 1317 pcfm(i,1)=pcfm1(i) 1318 pcfh(i,1)=pcfh1(i) 1319 ENDDO 1192 1320 c 1193 1321 c Calculer les coefficients turbulents dans l'atmosphere … … 1196 1324 itop(i) = isommet 1197 1325 ENDDO 1326 1327 PRINT*,' isommet=',isommet,' knon=',knon 1198 1328 1199 1329 DO k = 2, isommet … … 1310 1440 END 1311 1441 1312 SUBROUTINE clcdrag(knon, nsrf, zxli,1313 . u, v, t, q, zgeop,1314 . ts, qsurf, rugos,1315 . pcfm, pcfh, zcdn, zri)1316 c ================================================================= c1317 c Objet : calcul cdrags pour le moment et les flux chaleur sensible, latente (pcfm,pcfh)1318 c et du nombre de Richardson zri1319 c ================================================================= c1320 IMPLICIT NONE1321 #include "dimensions.h"1322 #include "dimphy.h"1323 #include "YOMCST.h"1324 #include "YOETHF.h"1325 #include "indicesol.h"1326 c1327 INTEGER knon, nsrf1328 REAL ts(klon), qsurf(klon)1329 REAL u(klon,klev), v(klon,klev), t(klon,klev), q(klon,klev)1330 REAL zgeop(klon,klev)1331 REAL rugos(klon), zri(klon)1332 c1333 REAL zcdn1334 REAL pcfm(klon,klev), pcfh(klon,klev)1335 c1336 c Quelques constantes et options:1337 c1338 REAL ckap, cb, cc, cd, cepdu21339 PARAMETER (ckap=0.35)1340 PARAMETER (cb=5.0)1341 PARAMETER (cc=5.0)1342 PARAMETER (cd=5.0)1343 PARAMETER (cepdu2 =(0.1)**2)1344 c1345 c Variables locales1346 INTEGER i1347 REAL zdu2, zdphi, ztsolv, ztvd, zscf, zucf, zcr1348 REAL friv, frih1349 REAL zcfm1(klon), zcfm2(klon)1350 REAL zcfh1(klon), zcfh2(klon)1351 c1352 c Fonctions thermodynamiques et fonctions d'instabilite1353 REAL fsta, fins, x1354 LOGICAL zxli1355 fsta(x) = 1.0 / (1.0+10.0*x*(1+8.0*x))1356 fins(x) = SQRT(1.0-18.0*x)1357 c1358 c Calculer le frottement au sol (Cdrag)1359 c1360 DO i = 1, knon1361 zdu2=max(cepdu2,u(i,1)**2+v(i,1)**2)1362 zdphi=zgeop(i,1)1363 c ztsolv = ts(i) * (1.0+RETV*q(i,1)) ! qsol approx = q(i,1)1364 ztsolv = ts(i) * (1.0+RETV*qsurf(i))1365 ztvd=(t(i,1)+zdphi/RCPD/(1.+RVTMP2*q(i,1)))1366 . *(1.+RETV*q(i,1))1367 zri(i)=zgeop(i,1)*(ztvd-ztsolv)/(zdu2*ztvd)1368 zcdn = (ckap/log(1.+zgeop(i,1)/(RG*rugos(i))))**21369 IF (zri(i) .ge. 0.) THEN ! situation stable1370 IF (.NOT.zxli) THEN1371 zscf=SQRT(1.+cd*ABS(zri(i)))1372 FRIV = AMAX1(1. / (1.+2.*CB*zri(i)/ZSCF), 0.1)1373 ! zcfm1(i) = zcdn/(1.+2.0*cb*zri(i)/ zscf)1374 zcfm1(i) = zcdn * FRIV1375 FRIH = AMAX1(1./ (1.+3.*CB*zri(i)*ZSCF), 0.1 )1376 ! zcfh1(i) = zcdn/(1.+3.0*cb*zri(i)*zscf)1377 zcfh1(i) = zcdn * FRIH1378 pcfm(i,1) = zcfm1(i)1379 pcfh(i,1) = zcfh1(i)1380 ELSE1381 pcfm(i,1) = zcdn* fsta(zri(i))1382 pcfh(i,1) = zcdn* fsta(zri(i))1383 ENDIF1384 ELSE ! situation instable1385 IF (.NOT.zxli) THEN1386 zucf=1./(1.+3.0*cb*cc*zcdn*SQRT(ABS(zri(i))1387 . *(1.0+zgeop(i,1)/(RG*rugos(i)))))1388 zcfm2(i) = zcdn*amax1((1.-2.0*cb*zri(i)*zucf),0.1)1389 zcfh2(i) = zcdn*amax1((1.-3.0*cb*zri(i)*zucf),0.1)1390 pcfm(i,1) = zcfm2(i)1391 pcfh(i,1) = zcfh2(i)1392 ELSE1393 pcfm(i,1) = zcdn* fins(zri(i))1394 pcfh(i,1) = zcdn* fins(zri(i))1395 ENDIF1396 zcr = (0.0016/(zcdn*SQRT(zdu2)))*ABS(ztvd-ztsolv)**(1./3.)1397 IF(nsrf.EQ.is_oce)pcfh(i,1)=zcdn*(1.0+zcr**1.25)**(1./1.25)1398 ENDIF1399 END DO1400 RETURN1401 END1402 1403 1442 SUBROUTINE coefkz2(nsrf, knon, paprs, pplay,t, 1404 1443 . pcfm, pcfh) … … 1611 1650 PARAMETER (isommet=klev) 1612 1651 REAL vk 1613 PARAMETER (vk=0. 35)1652 PARAMETER (vk=0.40) 1614 1653 REAL ricr 1615 1654 PARAMETER (ricr=0.4)
Note: See TracChangeset
for help on using the changeset viewer.
