Changeset 2594 for LMDZ5/branches/testing/libf/phylmd/cloudth.F90
- Timestamp:
- Jul 18, 2016, 9:41:10 PM (8 years ago)
- Location:
- LMDZ5/branches/testing
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ5/branches/testing
- Property svn:mergeinfo changed
/LMDZ5/trunk merged: 2547-2567,2569,2571-2574,2576-2589
- Property svn:mergeinfo changed
-
LMDZ5/branches/testing/libf/phylmd/cloudth.F90
r2544 r2594 5 5 6 6 7 USE IOIPSL, ONLY : getin8 7 IMPLICIT NONE 9 8 … … 20 19 #include "FCTTRE.h" 21 20 #include "thermcell.h" 21 #include "nuage.h" 22 22 23 23 INTEGER itap,ind1,ind2 … … 54 54 REAL sth,senv,sigma1s,sigma2s,xth,xenv 55 55 REAL Tbef,zdelta,qsatbef,zcor 56 REAL alpha,qlbef56 REAL qlbef 57 57 REAL ratqs(ngrid,klev) ! determine la largeur de distribution de vapeur 58 58 … … 62 62 REAL erf 63 63 64 REAL, SAVE :: iflag_cloudth_vert, iflag_cloudth_vert_omp=065 66 67 LOGICAL, SAVE :: first=.true.68 69 70 64 71 65 72 66 73 67 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 74 ! Astuce pour gérer deux versions de cloudth en attendant 75 ! de converger sur une version nouvelle 68 ! Gestion de deux versions de cloudth 76 69 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 77 IF (first) THEN 78 !$OMP MASTER 79 CALL getin('iflag_cloudth_vert',iflag_cloudth_vert_omp) 80 !$OMP END MASTER 81 !$OMP BARRIER 82 iflag_cloudth_vert=iflag_cloudth_vert_omp 83 first=.false. 84 ENDIF 85 IF (iflag_cloudth_vert==1) THEN 86 CALL cloudth_vert(ngrid,klev,ind2, & 70 71 IF (iflag_cloudth_vert.GE.1) THEN 72 CALL cloudth_vert(ngrid,klev,ind2, & 87 73 & ztv,po,zqta,fraca, & 88 74 & qcloud,ctot,zpspsk,paprs,ztla,zthl, & 89 75 & ratqs,zqs,t) 90 91 76 RETURN 77 ENDIF 92 78 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 93 94 79 95 80 … … 185 170 cenv(ind1,ind2)=0.5*(1.+1.*erf(xenv)) 186 171 ctot(ind1,ind2)=fraca(ind1,ind2)*cth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv(ind1,ind2) 187 ! ctot(ind1,ind2)=alpha*cth(ind1,ind2)+(1.-1.*alpha)*cenv(ind1,ind2)188 189 190 172 191 173 qlth(ind1,ind2)=sigma2s*((exp(-1.*xth**2)/sqrt2pi)+xth*sqrt(2.)*cth(ind1,ind2)) 192 174 qlenv(ind1,ind2)=sigma1s*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt(2.)*cenv(ind1,ind2)) 193 175 qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2) 194 ! qltot(ind1,ind2)=alpha*qlth(ind1,ind2)+(1.-1.*alpha)*qlenv(ind1,ind2)195 196 197 ! print*,senv,sth,sigma1s,sigma2s,fraca(ind1,ind2),'senv et sth et sig1 et sig2 et alpha'198 176 199 177 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! … … 269 247 & ratqs,zqs,t) 270 248 271 272 IMPLICIT NONE273 274 275 249 !=========================================================================== 276 250 ! Auteur : Arnaud Octavio Jam (LMD/CNRS) … … 280 254 281 255 256 USE ioipsl_getin_p_mod, ONLY : getin_p 257 258 IMPLICIT NONE 259 282 260 #include "YOMCST.h" 283 261 #include "YOETHF.h" 284 262 #include "FCTTRE.h" 285 263 #include "thermcell.h" 286 264 #include "nuage.h" 265 287 266 INTEGER itap,ind1,ind2 288 267 INTEGER ngrid,klev,klon,l,ig … … 320 299 REAL IntJ,IntI1,IntI2,IntI3,coeffqlenv,coeffqlth 321 300 REAL Tbef,zdelta,qsatbef,zcor 322 REAL alpha,qlbef301 REAL qlbef 323 302 REAL ratqs(ngrid,klev) ! determine la largeur de distribution de vapeur 303 ! Change the width of the PDF used for vertical subgrid scale heterogeneity 304 ! (J Jouhaud, JL Dufresne, JB Madeleine) 305 REAL,SAVE :: vert_alpha 306 LOGICAL, SAVE :: firstcall = .TRUE. 324 307 325 308 REAL zpdf_sig(ngrid),zpdf_k(ngrid),zpdf_delta(ngrid) … … 327 310 REAL zqs(ngrid), qcloud(ngrid) 328 311 REAL erf 329 330 331 332 333 312 334 313 !------------------------------------------------------------------------------ … … 355 334 sqrtpi=sqrt(pi) 356 335 357 336 IF (firstcall) THEN 337 vert_alpha=0.5 338 CALL getin_p('cloudth_vert_alpha',vert_alpha) 339 WRITE(*,*) 'cloudth_vert_alpha = ', vert_alpha 340 firstcall=.FALSE. 341 ENDIF 358 342 359 343 !------------------------------------------------------------------------------- … … 425 409 cenv(ind1,ind2)=0.5*(1.+1.*erf(xenv)) 426 410 ctot(ind1,ind2)=fraca(ind1,ind2)*cth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv(ind1,ind2) 427 ! ctot(ind1,ind2)=alpha*cth(ind1,ind2)+(1.-1.*alpha)*cenv(ind1,ind2)428 429 430 411 431 412 qlth(ind1,ind2)=sigma2s*((exp(-1.*xth**2)/sqrt2pi)+xth*sqrt(2.)*cth(ind1,ind2)) 432 413 qlenv(ind1,ind2)=sigma1s*((exp(-1.*xenv**2)/sqrt2pi)+xenv*sqrt(2.)*cenv(ind1,ind2)) 433 414 qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2) 434 ! qltot(ind1,ind2)=alpha*qlth(ind1,ind2)+(1.-1.*alpha)*qlenv(ind1,ind2)435 415 436 437 ! print*,senv,sth,sigma1s,sigma2s,fraca(ind1,ind2),'senv et sth et sig1 et sig2 et alpha' 438 439 416 IF (iflag_cloudth_vert == 1) THEN 440 417 !------------------------------------------------------------------------------- 441 418 ! Version 2: Modification selon J.-Louis. On condense ?? partir de qsat-ratqs … … 479 456 qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2) 480 457 458 ELSE IF (iflag_cloudth_vert == 2) THEN 459 460 !------------------------------------------------------------------------------- 461 ! Version 3: Modification Jean Jouhaud. On condense a partir de -delta s 462 !------------------------------------------------------------------------------- 463 ! deltasenv=aenv*ratqs(ind1,ind2)*po(ind1) 464 ! deltasth=ath*ratqs(ind1,ind2)*zqta(ind1,ind2) 465 ! deltasenv=aenv*ratqs(ind1,ind2)*zqsatenv(ind1,ind2) 466 ! deltasth=ath*ratqs(ind1,ind2)*zqsatth(ind1,ind2) 467 deltasenv=aenv*vert_alpha*sigma1s 468 deltasth=ath*vert_alpha*sigma2s 469 470 xenv1=-(senv+deltasenv)/(sqrt(2.)*sigma1s) 471 xenv2=-(senv-deltasenv)/(sqrt(2.)*sigma1s) 472 xth1=-(sth+deltasth)/(sqrt(2.)*sigma2s) 473 xth2=-(sth-deltasth)/(sqrt(2.)*sigma2s) 474 ! coeffqlenv=(sigma1s)**2/(2*sqrtpi*deltasenv) 475 ! coeffqlth=(sigma2s)**2/(2*sqrtpi*deltasth) 476 477 cth(ind1,ind2)=0.5*(1.-1.*erf(xth1)) 478 cenv(ind1,ind2)=0.5*(1.-1.*erf(xenv1)) 479 ctot(ind1,ind2)=fraca(ind1,ind2)*cth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*cenv(ind1,ind2) 480 481 IntJ=0.5*senv*(1-erf(xenv2))+(sigma1s/sqrt2pi)*exp(-1.*xenv2**2) 482 IntI1=(((senv+deltasenv)**2+(sigma1s)**2)/(8*deltasenv))*(erf(xenv2)-erf(xenv1)) 483 IntI2=(sigma1s**2/(4*deltasenv*sqrtpi))*(xenv1*exp(-1.*xenv1**2)-xenv2*exp(-1.*xenv2**2)) 484 IntI3=((sqrt2*sigma1s*(senv+deltasenv))/(4*sqrtpi*deltasenv))*(exp(-1.*xenv1**2)-exp(-1.*xenv2**2)) 485 486 ! IntI1=0.5*(0.5*sqrtpi*(erf(xenv2)-erf(xenv1))+xenv1*exp(-1.*xenv1**2)-xenv2*exp(-1.*xenv2**2)) 487 ! IntI2=xenv2*(exp(-1.*xenv2**2)-exp(-1.*xenv1**2)) 488 ! IntI3=0.5*sqrtpi*xenv2**2*(erf(xenv2)-erf(xenv1)) 489 490 qlenv(ind1,ind2)=IntJ+IntI1+IntI2+IntI3 491 ! qlenv(ind1,ind2)=IntJ 492 ! print*, qlenv(ind1,ind2),'VERIF EAU' 493 494 IntJ=0.5*sth*(1-erf(xth2))+(sigma2s/sqrt2pi)*exp(-1.*xth2**2) 495 IntI1=(((sth+deltasth)**2+(sigma2s)**2)/(8*deltasth))*(erf(xth2)-erf(xth1)) 496 IntI2=(sigma2s**2/(4*deltasth*sqrtpi))*(xth1*exp(-1.*xth1**2)-xth2*exp(-1.*xth2**2)) 497 IntI3=((sqrt2*sigma2s*(sth+deltasth))/(4*sqrtpi*deltasth))*(exp(-1.*xth1**2)-exp(-1.*xth2**2)) 498 499 qlth(ind1,ind2)=IntJ+IntI1+IntI2+IntI3 500 ! qlth(ind1,ind2)=IntJ 501 ! print*, IntJ,IntI1,IntI2,IntI3,qlth(ind1,ind2),'VERIF EAU2' 502 qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2) 503 504 505 506 507 ENDIF ! of if (iflag_cloudth_vert==1 or 2) 508 481 509 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 510 482 511 if (cenv(ind1,ind2).lt.1.e-10.or.cth(ind1,ind2).lt.1.e-10) then 483 512 ctot(ind1,ind2)=0.
Note: See TracChangeset
for help on using the changeset viewer.