Changeset 98 for LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F
- Timestamp:
- Jul 5, 2000, 4:58:04 PM (24 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F
r86 r98 1 1 SUBROUTINE clmain(dtime,pctsrf,t,q,u,v, 2 . soil_model,ts,soilcap,soilflux, 3 . paprs,pplay,radsol,snow,qsol, 4 . xlat, rugos, 2 . ok_veget,ts, 3 . paprs,pplay,radsol,snow,qsol,evap,albe, 4 . rain_f, snow_f, solsw, sollw, 5 . rlon, rlat, rugos, 6 . debut, lafin, 5 7 . d_t,d_q,d_u,d_v,d_ts, 6 8 . flux_t,flux_q,flux_u,flux_v,cdragh,cdragm, … … 41 43 c beta-----input-R- coefficient de l'evaporation reelle (0 a 1) 42 44 c dif_grnd-input-R- coeff. de diffusion (chaleur) vers le sol profond 43 c xlat-----input-R- latitude en degree45 c rlat-----input-R- latitude en degree 44 46 c rugos----input-R- longeur de rugosite (en m) 45 47 c … … 75 77 REAL u(klon,klev), v(klon,klev) 76 78 REAL paprs(klon,klev+1), pplay(klon,klev), radsol(klon) 77 REAL xlat(klon)79 REAL rlon(klon), rlat(klon) 78 80 REAL d_t(klon, klev), d_q(klon, klev) 79 81 REAL d_u(klon, klev), d_v(klon, klev) 80 REAL flux_t(klon,klev ), flux_q(klon,klev)82 REAL flux_t(klon,klev, nbsrf), flux_q(klon,klev, nbsrf) 81 83 REAL dflux_t(klon), dflux_q(klon) 82 REAL flux_u(klon,klev ), flux_v(klon,klev)84 REAL flux_u(klon,klev, nbsrf), flux_v(klon,klev, nbsrf) 83 85 REAL rugmer(klon) 84 86 REAL cdragh(klon), cdragm(klon) 87 LOGICAL debut, lafin, ok_veget 85 88 cAA INTEGER itr 86 89 cAA REAL tr(klon,klev,nbtr) … … 93 96 REAL snow(klon,nbsrf) 94 97 REAL qsol(klon,nbsrf) 98 REAL evap(klon,nbsrf) 99 REAL albe(klon,nbsrf) 100 real rain_f(klon), snow_f(klon) 101 REAL sollw(klon), solsw(klon) 95 102 REAL rugos(klon,nbsrf) 96 103 cAA … … 103 110 c====================================================================== 104 111 REAL yts(klon), yrugos(klon), ypct(klon) 105 REAL ycal(klon), ybeta(klon), ydif(klon) 112 REAL ycal(klon), ybeta(klon), ydif(klon), yalb(klon),yevap(klon) 106 113 REAL yu1(klon), yv1(klon) 114 real ysnow(klon), yqsol(klon) 115 real yrain_f(klon), ysnow_f(klon) 116 real ysollw(klon), ysolsw(klon), ysolswnet(klon) 107 117 REAL yrugm(klon), yrads(klon) 108 118 REAL y_d_ts(klon) … … 144 154 c====================================================================== 145 155 146 write(*,*)'CLMAIN.NEW'147 148 156 DO k = 1, klev ! epaisseur de couche 149 157 DO i = 1, klon … … 174 182 d_ts(i,nsrf) = 0.0 175 183 ENDDO 176 ENDDO 184 END DO 185 C§§§ PB 186 flux_t = 0. 187 flux_q = 0. 188 flux_u = 0. 189 flux_v = 0. 177 190 DO k = 1, klev 178 191 DO i = 1, klon 179 192 d_t(i,k) = 0.0 180 193 d_q(i,k) = 0.0 181 flux_t(i,k) = 0.0182 flux_q(i,k) = 0.0194 c$$$ flux_t(i,k) = 0.0 195 c$$$ flux_q(i,k) = 0.0 183 196 d_u(i,k) = 0.0 184 197 d_v(i,k) = 0.0 185 flux_u(i,k) = 0.0186 flux_v(i,k) = 0.0198 c$$$ flux_u(i,k) = 0.0 199 c$$$ flux_v(i,k) = 0.0 187 200 zcoefh(i,k) = 0.0 188 201 ENDDO … … 203 216 c 204 217 c prescrire les proprietes du sol: 205 CALL calbeta(dtime,nsrf,snow,qsol, beta,capsol,dif_grnd) 206 IF (.NOT.soil_model) THEN 207 DO i = 1, klon 208 cal(i) = RCPD * capsol(i) 209 totalflu(i) = radsol(i) 210 ENDDO 211 ELSE 212 DO i = 1, klon 213 totalflu(i) = soilflux(i,nsrf) + radsol(i) 214 IF (nsrf.EQ.is_oce) THEN 215 cal(i) = 0.0 216 ELSE 217 cal(i) = RCPD / soilcap(i,nsrf) 218 ENDIF 219 ENDDO 220 ENDIF 221 c 218 c CALL calbeta(dtime,nsrf,snow,qsol, beta,capsol,dif_grnd) 219 c IF (.NOT.soil_model) THEN 220 c DO i = 1, klon 221 c cal(i) = RCPD * capsol(i) 222 c totalflu(i) = radsol(i) 223 c ENDDO 224 c ELSE 225 c DO i = 1, klon 226 c totalflu(i) = soilflux(i,nsrf) + radsol(i) 227 c IF (nsrf.EQ.is_oce) THEN 228 c cal(i) = 0.0 229 c ELSE 230 c cal(i) = RCPD / soilcap(i,nsrf) 231 c ENDIF 232 c ENDDO 233 c ENDIF 234 c 235 totalflu = radsol 236 222 237 c chercher les indices: 223 238 DO j = 1, klon … … 237 252 ypct(j) = pctsrf(i,nsrf) 238 253 yts(j) = ts(i,nsrf) 254 ysnow(j) = snow(i,nsrf) 255 yevap(j) = evap(i,nsrf) 256 yqsol(j) = qsol(i,nsrf) 257 yalb(j) = albe(i,nsrf) 258 yrain_f(j) = rain_f(i) 259 ysnow_f(j) = snow_f(i) 260 ysolsw(j) = solsw(i) 261 ysollw(j) = sollw(i) 239 262 yrugos(j) = rugos(i,nsrf) 240 263 yu1(j) = u1lay(i) 241 264 yv1(j) = v1lay(i) 242 265 yrads(j) = totalflu(i) 243 ycal(j) = cal(i)244 ybeta(j) = beta(i)245 ydif(j) = dif_grnd(i)266 c ycal(j) = cal(i) 267 c ybeta(j) = beta(i) 268 c ydif(j) = dif_grnd(i) 246 269 ypaprs(j,klev+1) = paprs(i,klev+1) 247 270 ENDDO … … 259 282 ENDDO 260 283 c 261 cAA IF (itr.GE.1) THEN262 cAA DO it = 1, itr263 cAA DO k = 1, klev264 cAA DO j = 1, knon265 cAA i = ni(j)266 cAA ytr(j,k,it) = tr(i,k,it)267 cAA ENDDO268 cAA ENDDO269 cAA DO j = 1, knon270 cAA i = ni(j)271 cAA yflxsrf(j,it) = flux_surf(i,it)272 cAA ENDDO273 cAA ENDDO274 cAA ENDIF275 284 c 276 285 c calculer Cdrag et les coefficients d'echange … … 287 296 ENDDO 288 297 c 289 c parametrisation non-locale:290 IF (ok_nonloc) THEN291 DO i = 1, knon292 y_cd_h(i) = ycoefh(i,1)293 y_cd_m(i) = ycoefm(i,1)294 ENDDO295 CALL nonlocal(knon, ypaprs, ypplay,296 . yts,ybeta,yu,yv,yt,yq,297 . y_cd_h, y_cd_m, ycoefm0, ycoefh0, ygamt, ygamq)298 DO k = 1, klev299 DO i = 1, knon300 ycoefm(i,k) = MAX(ycoefm(i,k),ycoefm0(i,k))301 ycoefh(i,k) = MAX(ycoefh(i,k),ycoefh0(i,k))302 ENDDO303 ENDDO304 ELSE305 IF (.NOT.contreg) THEN306 DO k = 2, klev307 DO i = 1, knon308 ygamq(i,k) = 0.0309 ygamt(i,k) = 0.0310 ENDDO311 ENDDO312 ELSE313 DO k = 3, klev314 DO i = 1, knon315 ygamq(i,k) = 0.0316 ygamt(i,k) = -1.0E-03317 ENDDO318 ENDDO319 DO i = 1, knon320 ygamq(i,2) = 0.0321 ygamt(i,2) = -2.5E-03322 ENDDO323 ENDIF324 ENDIF325 298 c 326 299 c calculer la diffusion de "q" et de "h" 327 CALL clqh(knon, dtime, nsrf,yu1, yv1, 300 CALL clqh(knon, dtime, nsrf, ni, pctsrf, rlon, rlat, 301 e yu1, yv1, 328 302 e ycoefh,yt,yq,yts,ypaprs,ypplay,ydelp,yrads, 329 e ycal,ybeta,ydif,ygamt,ygamq, 303 e yevap,yalb, ysnow, yqsol, yrain_f, ysnow_f, 304 e ysollw, ysolsw, 330 305 s y_d_t, y_d_q, y_d_ts, 331 306 s y_flux_t, y_flux_q, y_dflux_t, y_dflux_q) … … 363 338 c 364 339 DO k = 1, klev 340 DO j = 1, knon 341 i = ni(j) 342 ycoefh(j,k) = ycoefh(j,k) * ypct(j) 343 ycoefm(j,k) = ycoefm(j,k) * ypct(j) 344 y_d_t(j,k) = y_d_t(j,k) * ypct(j) 345 y_d_q(j,k) = y_d_q(j,k) * ypct(j) 346 C§§§ PB 347 flux_t(i,k,nsrf) = y_flux_t(j,k) 348 flux_q(i,k,nsrf) = y_flux_q(j,k) 349 flux_u(i,k,nsrf) = y_flux_u(j,k) 350 flux_v(i,k,nsrf) = y_flux_v(j,k) 351 c$$$ PB y_flux_t(j,k) = y_flux_t(j,k) * ypct(j) 352 c$$$ PB y_flux_q(j,k) = y_flux_q(j,k) * ypct(j) 353 y_d_u(j,k) = y_d_u(j,k) * ypct(j) 354 y_d_v(j,k) = y_d_v(j,k) * ypct(j) 355 c$$$ PB y_flux_u(j,k) = y_flux_u(j,k) * ypct(j) 356 c$$$ PB y_flux_v(j,k) = y_flux_v(j,k) * ypct(j) 357 ENDDO 358 ENDDO 359 360 evap(:,nsrf) = - flux_q(:,1,nsrf) 361 c 365 362 DO j = 1, knon 366 ycoefh(j,k) = ycoefh(j,k) * ypct(j) 367 ycoefm(j,k) = ycoefm(j,k) * ypct(j) 368 y_d_t(j,k) = y_d_t(j,k) * ypct(j) 369 y_d_q(j,k) = y_d_q(j,k) * ypct(j) 370 y_flux_t(j,k) = y_flux_t(j,k) * ypct(j) 371 y_flux_q(j,k) = y_flux_q(j,k) * ypct(j) 372 y_d_u(j,k) = y_d_u(j,k) * ypct(j) 373 y_d_v(j,k) = y_d_v(j,k) * ypct(j) 374 y_flux_u(j,k) = y_flux_u(j,k) * ypct(j) 375 y_flux_v(j,k) = y_flux_v(j,k) * ypct(j) 376 ENDDO 377 ENDDO 378 c 379 DO j = 1, knon 380 i = ni(j) 363 i = ni(j) 381 364 d_ts(i,nsrf) = y_d_ts(j) 365 albe(i,nsrf) = yalb(j) 366 snow(i,nsrf) = ysnow(j) 367 qsol(i,nsrf) = yqsol(j) 382 368 rugmer(i) = yrugm(j) 383 369 cdragh(i) = cdragh(i) + ycoefh(j,1) … … 400 386 d_t(i,k) = d_t(i,k) + y_d_t(j,k) 401 387 d_q(i,k) = d_q(i,k) + y_d_q(j,k) 402 403 flux_q(i,k) = flux_q(i,k) + y_flux_q(j,k)388 c$$$ PB flux_t(i,k) = flux_t(i,k) + y_flux_t(j,k) 389 c$$$ flux_q(i,k) = flux_q(i,k) + y_flux_q(j,k) 404 390 d_u(i,k) = d_u(i,k) + y_d_u(j,k) 405 391 d_v(i,k) = d_v(i,k) + y_d_v(j,k) 406 407 flux_v(i,k) = flux_v(i,k) + y_flux_v(j,k)392 c$$$ PB flux_u(i,k) = flux_u(i,k) + y_flux_u(j,k) 393 c$$$ flux_v(i,k) = flux_v(i,k) + y_flux_v(j,k) 408 394 zcoefh(i,k) = zcoefh(i,k) + ycoefh(j,k) 409 395 ENDDO … … 430 416 99999 CONTINUE 431 417 c 418 432 419 RETURN 433 420 END 434 SUBROUTINE clqh(knon,dtime,nisurf,u1lay,v1lay,coef, 421 SUBROUTINE clqh(knon,dtime,nisurf,knindex,pctsrf, rlon, rlat, 422 e u1lay,v1lay,coef, 435 423 e t,q,ts,paprs,pplay, 436 e delp,radsol,cal,beta,dif_grnd, gamt,gamq, 424 e delp,radsol,evap,albedo,snow,qsol, 425 e precip_rain, precip_snow, 426 e lwdown, swdown, 437 427 s d_t, d_q, d_ts, flux_t, flux_q,dflux_s,dflux_l) 438 428 … … 446 436 #include "dimensions.h" 447 437 #include "dimphy.h" 438 #include "YOMCST.h" 439 #include "YOETHF.h" 440 #include "FCTTRE.h" 441 #include "indicesol.h" 448 442 c Arguments: 449 443 INTEGER knon … … 462 456 REAL q(klon,klev) ! humidite specifique (kg/kg) 463 457 REAL ts(klon) ! temperature du sol (K) 458 REAL evap(klon) ! evaporation au sol 464 459 REAL paprs(klon,klev+1) ! pression a inter-couche (Pa) 465 460 REAL pplay(klon,klev) ! pression au milieu de couche (Pa) 466 461 REAL delp(klon,klev) ! epaisseur de couche en pression (Pa) 467 462 REAL radsol(klon) ! ray. net au sol (Solaire+IR) W/m2 463 REAL albedo(klon) ! albedo de la surface 464 REAL snow(klon) ! hauteur de neige 465 REAL qsol(klon) ! humidite de la surface 466 real precip_rain(klon), precip_snow(klon) 467 integer knindex(klon) 468 real pctsrf(klon,nbsrf) 469 real rlon(klon), rlat(klon) 468 470 c 469 471 REAL d_t(klon,klev) ! incrementation de "t" … … 521 523 c====================================================================== 522 524 REAL zcor, zdelta, zcvm5 523 #include "YOMCST.h" 524 #include "YOETHF.h" 525 #include "FCTTRE.h" 526 #include "indicesol.h" 525 logical contreg 526 parameter (contreg=.true.) 527 527 c====================================================================== 528 528 c Rajout pour l'interface … … 530 530 integer jour 531 531 integer nisurf 532 integer knindex(klon)533 532 logical debut, lafin, ok_veget 534 real rlon(klon), rlat(klon) 535 real zlev(klon), zlflu(klon) 533 real zlev1(klon) 536 534 real temp_air(klon), spechum(klon) 537 535 real hum_air(klon), ccanopy(klon) 538 536 real tq_cdrag(klon), petAcoef(klon), peqAcoef(klon) 539 537 real petBcoef(klon), peqBcoef(klon) 540 real precip_rain(klon), precip_snow(klon)541 538 real lwdown(klon), swnet(klon), swdown(klon), ps(klon) 542 539 real p1lay(klon) … … 545 542 546 543 ! Parametres de sortie 547 real evap(klon),fluxsens(klon), fluxlat(klon)544 real fluxsens(klon), fluxlat(klon) 548 545 real tsol_rad(klon), tsurf_new(klon), alb_new(klon) 549 546 real emis_new(klon), z0_new(klon) 550 real dflux_l(klon), dflux_s(klon)551 547 real pctsrf_new(klon,nbsrf) 552 548 c 549 550 if (.not. contreg) then 551 do k = 2, klev 552 do i = 1, knon 553 gamq(i,k) = 0.0 554 gamt(i,k) = 0.0 555 enddo 556 enddo 557 else 558 do k = 3, klev 559 do i = 1, knon 560 gamq(i,k)= 0.0 561 gamt(i,k)= -1.0e-03 562 enddo 563 enddo 564 do i = 1, knon 565 gamq(i,2) = 0.0 566 gamt(i,2) = -2.5e-03 567 enddo 568 endif 553 569 554 570 DO i = 1, knon … … 650 666 spechum=q(:,1) 651 667 p1lay = pplay(:,1) 668 zlev1 = delp(:,1) 669 swnet = swdown * (1 - albedo) 670 c En attendant mieux 671 hum_air = 0. 672 ccanopy = 0. 673 tq_cdrag = 0. 652 674 653 675 CALL interfsurf(itime, dtime, jour, 654 . klon, nisurf, knon, knindex, rlon, rlat,676 . klon, iim, jjm, nisurf, knon, knindex, pctsrf, rlon, rlat, 655 677 . debut, lafin, ok_veget, 656 . zlev , zlflu,u1lay, v1lay, temp_air, spechum, hum_air, ccanopy,678 . zlev1, u1lay, v1lay, temp_air, spechum, hum_air, ccanopy, 657 679 . tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, 658 680 . precip_rain, precip_snow, lwdown, swnet, swdown, 659 . ts, p1lay, cal, beta, coef1lay, psref, radsol, dif_grnd, 681 . albedo, snow, qsol, 682 . ts, p1lay, coef1lay, psref, radsol, 660 683 . ocean, 661 684 . evap, fluxsens, fluxlat, dflux_l, dflux_s, 662 . tsol_rad, tsurf_new, alb_new, emis_new, z0_new, pctsrf_new) 685 . tsol_rad, tsurf_new, alb_new, emis_new, z0_new, pctsrf_new, 686 . zmasq) 663 687 664 688 flux_t(:,1) = fluxsens 665 flux_q(:,1) = fluxlat689 flux_q(:,1) = - evap 666 690 d_ts = tsurf_new - ts 691 albedo = alb_new 667 692 668 693 c==== une fois on a zx_h_ts, on peut faire l'iteration ======== … … 1228 1253 RETURN 1229 1254 END 1230 SUBROUTINE calbeta(dtime,indice, snow,qsol,1255 SUBROUTINE calbeta(dtime,indice,knon,snow,qsol, 1231 1256 . vbeta,vcal,vdif) 1232 1257 IMPLICIT none … … 1238 1263 c Calculer quelques parametres pour appliquer la couche limite 1239 1264 c ------------------------------------------------------------ 1240 #include "dimensions.h"1241 #include "dimphy.h"1265 !#include "dimensions.h" 1266 !#include "dimphy.h" 1242 1267 #include "YOMCST.h" 1243 1268 #include "indicesol.h" … … 1256 1281 c 1257 1282 REAL dtime 1258 REAL snow(k lon,nbsrf), qsol(klon,nbsrf)1259 INTEGER indice 1260 C 1261 REAL vbeta(k lon)1262 REAL vcal(k lon)1263 REAL vdif(k lon)1283 REAL snow(knon,nbsrf), qsol(knon,nbsrf) 1284 INTEGER indice, knon 1285 C 1286 REAL vbeta(knon) 1287 REAL vcal(knon) 1288 REAL vdif(knon) 1264 1289 C 1265 1290 IF (indice.EQ.is_oce) THEN 1266 DO i = 1, k lon1291 DO i = 1, knon 1267 1292 vcal(i) = 0.0 1268 1293 vbeta(i) = 1.0 … … 1272 1297 c 1273 1298 IF (indice.EQ.is_sic) THEN 1274 DO i = 1, k lon1299 DO i = 1, knon 1275 1300 vcal(i) = calice 1276 1301 IF (snow(i,is_sic) .GT. 0.0) vcal(i) = calsno … … 1282 1307 c 1283 1308 IF (indice.EQ.is_ter) THEN 1284 DO i = 1, k lon1309 DO i = 1, knon 1285 1310 vcal(i) = calsol 1286 1311 IF (snow(i,is_ter) .GT. 0.0) vcal(i) = calsno … … 1291 1316 c 1292 1317 IF (indice.EQ.is_lic) THEN 1293 DO i = 1, k lon1318 DO i = 1, knon 1294 1319 vcal(i) = calice 1295 1320 IF (snow(i,is_lic) .GT. 0.0) vcal(i) = calsno
Note: See TracChangeset
for help on using the changeset viewer.