Changeset 1901
- Timestamp:
- Nov 19, 2013, 10:29:06 AM (11 years ago)
- File:
-
- 1 edited
-
LMDZ5/trunk/libf/phylmd/fisrtilp.F90 (modified) (7 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ5/trunk/libf/phylmd/fisrtilp.F90
r1894 r1901 100 100 ! 101 101 INTEGER i, k, n, kk 102 REAL zqs(klon), zdqs(klon), zdelta, zcor, zcvm5 102 REAL zqs(klon), zdqs(klon), zdelta, zcor, zcvm5 103 REAL Tbef(klon),qlbef(klon),DT(klon),num,denom 104 LOGICAL convergence(klon) 105 REAL DDT0 106 PARAMETER (DDT0=.01) 107 INTEGER n_i, iter 108 103 109 REAL zrfl(klon), zrfln(klon), zqev, zqevt 104 110 REAL zifl(klon), zifln(klon), zqev0,zqevi, zqevti … … 467 473 ! de l'eau condensee: 468 474 ! 475 !verification de la valeur de iflag_fisrtilp_qsat pour iflag_ice_thermo=1 476 if ((iflag_ice_thermo.eq.1).and.(iflag_fisrtilp_qsat.ne.0)) then 477 write(*,*) " iflag_ice_thermo==1 requires iflag_fisrtilp_qsat==0", & 478 " but iflag_fisrtilp_qsat=",iflag_fisrtilp_qsat, ". Might as well stop here." 479 stop 480 endif 469 481 470 482 IF (cpartiel) THEN … … 526 538 527 539 do i=1,klon 540 Tbef(i)=zt(i) 541 qlbef(i)=0. 528 542 if (lognormale(i)) then 529 543 zpdf_sig(i)=ratqs(i,k)*zq(i) … … 538 552 zpdf_e2(i)=sign(min(abs(zpdf_e2(i)),5.),zpdf_e2(i)) 539 553 zpdf_e2(i)=1.-erf(zpdf_e2(i)) 540 endif 541 enddo 542 543 do i=1,klon 544 if (lognormale(i)) then 554 545 555 if (zpdf_e1(i).lt.1.e-10) then 546 556 rneb(i,k)=0. … … 550 560 zqn(i)=zq(i)*zpdf_e2(i)/zpdf_e1(i) 551 561 endif 562 563 qlbef(i)=max(0.,zqn(i)-zqs(i)) 564 num=-Tbef(i)+zt(i)+rneb(i,k)*RLVTT/RCPD/(1.0+RVTMP2*zq(i))*qlbef(i) 565 denom=1.+rneb(i,k)*zdqs(i) 566 DT(i)=num/denom 552 567 endif 553 554 568 enddo 569 570 n_i=1 571 ! do while ((abs(DT(i)).gt.DDT0).and.((zqn(i)-zqs(i)).gt.0.)) 572 ! iflag_fisrtilp_qsat=0: qsat ne varie pas avec T 573 ! iflag_fisrtilp_qsat > 1 : nombre d iterations max pour converger sur le calcul de qsat(T) 574 ! do while (n_i.le.iflag_fisrtilp_qsat) 575 if (iflag_fisrtilp_qsat.ge.1) then 576 do iter=1,iflag_fisrtilp_qsat 577 do i=1,klon 578 convergence(i)=abs(DT(i)).gt.DDT0 579 if (convergence(i).and.lognormale(i)) then 580 Tbef(i)=Tbef(i)+DT(i) 581 zdelta=MAX(0.,SIGN(1.,RTT-Tbef(i))) 582 zcvm5 = R5LES*RLVTT*(1.-zdelta) + R5IES*RLSTT*zdelta 583 zcvm5 = zcvm5 /RCPD/(1.0+RVTMP2*zq(i)) 584 zqs(i)= R2ES * FOEEW(Tbef(i),zdelta)/pplay(i,k) 585 zqs(i)=MIN(0.5,zqs(i)) 586 zcor=1./(1.-retv*zqs(i)) 587 zqs(i)=zqs(i)*zcor 588 zdqs(i)=FOEDE(Tbef(i),zdelta,zcvm5,zqs(i),zcor) 589 590 zpdf_sig(i)=ratqs(i,k)*zq(i) 591 zpdf_k(i)=-sqrt(log(1.+(zpdf_sig(i)/zq(i))**2)) 592 zpdf_delta(i)=log(zq(i)/zqs(i)) 593 zpdf_a(i)=zpdf_delta(i)/(zpdf_k(i)*sqrt(2.)) 594 zpdf_b(i)=zpdf_k(i)/(2.*sqrt(2.)) 595 zpdf_e1(i)=zpdf_a(i)-zpdf_b(i) 596 zpdf_e1(i)=sign(min(abs(zpdf_e1(i)),5.),zpdf_e1(i)) 597 zpdf_e1(i)=1.-erf(zpdf_e1(i)) 598 zpdf_e2(i)=zpdf_a(i)+zpdf_b(i) 599 zpdf_e2(i)=sign(min(abs(zpdf_e2(i)),5.),zpdf_e2(i)) 600 zpdf_e2(i)=1.-erf(zpdf_e2(i)) 601 602 if (zpdf_e1(i).lt.1.e-10) then 603 rneb(i,k)=0. 604 zqn(i)=zqs(i) 605 else 606 rneb(i,k)=0.5*zpdf_e1(i) 607 zqn(i)=zq(i)*zpdf_e2(i)/zpdf_e1(i) 608 endif 609 610 qlbef(i)=max(0.,zqn(i)-zqs(i)) 611 num=-Tbef(i)+zt(i)+rneb(i,k)*RLVTT/RCPD/(1.0+RVTMP2*zq(i))*qlbef(i) 612 denom=1.+rneb(i,k)*zdqs(i) 613 DT(i)=num/denom 614 n_i=n_i+1 615 endif 616 enddo 617 enddo 618 619 endif 555 620 556 621 557 622 endif ! iflag_pdf 558 623 624 if (iflag_fisrtilp_qsat.eq.-1) then 625 !CR: ATTENTION option fausse mais a existe 559 626 DO i=1,klon 560 627 IF (rneb(i,k) .LE. 0.0) THEN … … 565 632 ELSE IF (rneb(i,k) .GE. 1.0) THEN 566 633 zqn(i) = zq(i) 567 rneb(i,k) = 1.0 568 zcond(i) = MAX(0.0,zqn(i)-zqs(i))/(1 .+iflag_fisrtilp_qsat*zdqs(i))634 rneb(i,k) = 1.0 635 zcond(i) = MAX(0.0,zqn(i)-zqs(i))/(1+zdqs(i)) 569 636 rhcl(i,k)=1.0 570 637 ELSE 571 zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1 .+iflag_fisrtilp_qsat*zdqs(i))638 zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1+zdqs(i)) 572 639 rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i) 573 640 ENDIF 574 641 ENDDO 642 643 ELSE 644 645 DO i=1,klon 646 IF (rneb(i,k) .LE. 0.0) THEN 647 zqn(i) = 0.0 648 rneb(i,k) = 0.0 649 zcond(i) = 0.0 650 rhcl(i,k)=zq(i)/zqs(i) 651 ELSE IF (rneb(i,k) .GE. 1.0) THEN 652 zqn(i) = zq(i) 653 rneb(i,k) = 1.0 654 zcond(i) = MAX(0.0,zqn(i)-zqs(i)) 655 rhcl(i,k)=1.0 656 ELSE 657 zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k) 658 rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i) 659 ENDIF 660 ENDDO 661 662 663 ENDIF 664 575 665 ! do i=1,klon 576 666 ! IF (rneb(i,k) .LE. 0.0) zqn(i) = 0.0 … … 606 696 !AJ< 607 697 IF (.NOT. ice_thermo) THEN 608 DO i = 1, klon 609 zt(i) = zt(i) + zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i)) 610 ENDDO 698 if (iflag_fisrtilp_qsat.lt.1) then 699 DO i = 1, klon 700 zt(i) = zt(i) + zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i)) 701 ENDDO 702 else if (iflag_fisrtilp_qsat.gt.0) then 703 DO i= 1, klon 704 if (lognormale(i)) then 705 zt(i)=Tbef(i) 706 else 707 zt(i) = zt(i) + zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*(zq(i)+zcond(i))) 708 endif 709 ENDDO 710 endif 611 711 ELSE 612 DO i = 1, klon 613 zfice(i) = 1.0 - (zt(i)-ztglace) / (273.15-ztglace) 614 zfice(i) = MIN(MAX(zfice(i),0.0),1.0) 615 zfice(i) = zfice(i)**nexpo 616 zt(i) = zt(i) + (1.-zfice(i))*zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i)) & 712 if (iflag_fisrtilp_qsat.lt.1) then 713 DO i = 1, klon 714 zfice(i) = 1.0 - (zt(i)-ztglace) / (273.15-ztglace) 715 zfice(i) = MIN(MAX(zfice(i),0.0),1.0) 716 zfice(i) = zfice(i)**nexpo 717 zt(i) = zt(i) + (1.-zfice(i))*zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i)) & 617 718 +zfice(i)*zcond(i) * RLSTT/RCPD/(1.0+RVTMP2*zq(i)) 719 ENDDO 720 else 721 DO i=1, klon 722 zfice(i) = 1.0 - (zt(i)-ztglace) / (273.15-ztglace) 723 zfice(i) = MIN(MAX(zfice(i),0.0),1.0) 724 zfice(i) = zfice(i)**nexpo 725 !CR: ATTENTION zt different de Tbef: à corriger 726 zt(i) = zt(i) + (1.-zfice(i))*zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*(zq(i)+zcond(i))) & 727 +zfice(i)*zcond(i) * RLSTT/RCPD/(1.0+RVTMP2*(zq(i)+zcond(i))) 728 ENDDO 729 endif 618 730 ! print*,zt(i),zrfl(i),zifl(i),'temp1' 619 ENDDO620 731 ENDIF 621 732 !>AJ
Note: See TracChangeset
for help on using the changeset viewer.
