Changeset 1667 for trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus

Timestamp:

Feb 17, 2017, 3:04:05 PM (8 years ago)

Author:

slebonnois

Message:

SL: an other round of bug corrections for Venus clouds

Location:

trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus

Files:

• mad_muphy.F90 (modified) (11 diffs)
• new_cloud_sedim.F (moved) (moved from trunk/LMDZ.VENUS/libf/phyvenus/new_cloud_sedim.F)
• new_cloud_venus.F (moved) (moved from trunk/LMDZ.VENUS/libf/phyvenus/new_cloud_venus.F) (2 diffs)
• wsa_new.F90 (modified) (7 diffs)

trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus/mad_muphy.F90

@@ -122 +122 @@
   IF (SEDIM .and. MERGE) THEN
      redge = (r1*r2)**0.5D0 ! Edge radius of mode 1
-     write(*,*)'The intersection size or virtual bonduary', redge
+!     write(*,*)'The intersection size or virtual bonduary', redge
 
      IF (r1 .EQ. redge) THEN
-        write(*,*)'Mode 1 merge to mode 2', redge, r1
+!        write(*,*)'Mode 1 merge to mode 2', redge, r1
         CALL MERGING(r1,N1,sigma1,sigma2,dM0_merge_m1,dM3_merge_m1, &
              &   dM0_merge_m2,dM3_merge_m2)
@@ -143 +143 @@
 
      IF (r2 .EQ. redge) THEN
-        write(*,*)'Mode 2 merge to mode 1',redge,r2
+!        write(*,*)'Mode 2 merge to mode 1',redge,r2
         CALL MERGING(r2,N2,sigma2,sigma1,dM0_merge_m2,dM3_merge_m2, &
              &   dM0_merge_m1,dM3_merge_m1)
@@ -325 +325 @@
   M3_m2(3) = M3_m2(3) + dM3_m2_ccn
   !------------------------------------!
-  write(*,*)'l.653: after up', M3_m1(1),M3_m1(2),M3_m1(3),dM3_m1_SA,dM3_m1_WV,dM3_m1_ccn
-  write(*,*)'l.654: after up', M0_m1(1),M0_m1(2),dM0_m1_drop
+!  write(*,*)'l.653: after up', M3_m1(1),M3_m1(2),M3_m1(3),dM3_m1_SA,dM3_m1_WV,dM3_m1_ccn
+!  write(*,*)'l.654: after up', M0_m1(1),M0_m1(2),dM0_m1_drop
 !-----------------------------------------------------------------------
   ! To controle/check the mass conservation !
@@ -375 +375 @@
         N1 = 0.D0
      ENDIF
-     write(*,*)'rmean not ZERO',r1,N1,M0_m1(1),M3_m1(1),WSA
+!     write(*,*)'rmean not ZERO',r1,N1,M0_m1(1),M3_m1(1),WSA
   ENDIF
   
@@ -408 +408 @@
   frac_ccn_m2 = M0_m2(2)/(M0_m2(1)+M0_m2(2))
 
-  write(*,*)'rmean',r1,'NTOT',N1, 'SH2SO4', SH2SO4
+!  write(*,*)'rmean',r1,'NTOT',N1, 'SH2SO4', SH2SO4
 
 
@@ -416 +416 @@
   IF (MERGE) THEN
      redge = (r1*r2)**0.5D0 ! Edge radius of mode 1
-     write(*,*)'The intersection size or virtual bonduary', redge
+!     write(*,*)'The intersection size or virtual bonduary', redge
 
      IF (r1 .EQ. redge) THEN
-        write(*,*)'Mode 1 merge to mode 2', redge, r1
+!        write(*,*)'Mode 1 merge to mode 2', redge, r1
         CALL MERGING(r1,N1,sigma1,sigma2,dM0_merge_m1,dM3_merge_m1, &
              &   dM0_merge_m2,dM3_merge_m2)
@@ -437 +437 @@
 
      IF (r2 .EQ. redge) THEN
-        write(*,*)'Mode 2 merge to mode 1',redge,r2
+!        write(*,*)'Mode 2 merge to mode 1',redge,r2
         CALL MERGING(r2,N2,sigma2,sigma1,dM0_merge_m2,dM3_merge_m2, &
              &   dM0_merge_m1,dM3_merge_m1)
@@ -476 +476 @@
 !-----------------------------------------------------------------------
   if (dt .eq. 1) then
-     write(*,*) 'Ntot, rc mode 1: ', N1,r1
+!     write(*,*) 'Ntot, rc mode 1: ', N1,r1
      ! Log-normal distribution function
      CALL build_radius_grid(nbin,rmin,rmax,vratio)
      CALL logdist(sigma1,N1,r1,rad_cld,n_g1)
-     DO i = 1, nbin
-        open(888,file='ND_t3600_nuc-het-flux_dt1_1step_sssat')
-        write(888,'(I10,2X,3(ES15.7,1X))') i, rad_cld(i), n_g1(i)
-     ENDDO
+!     DO i = 1, nbin
+!        open(888,file='ND_t3600_nuc-het-flux_dt1_1step_sssat')
+!        write(888,'(I10,2X,3(ES15.7,1X))') i, rad_cld(i), n_g1(i)
+!     ENDDO
   endif
 
@@ -507 +507 @@
 
   MFPAIR = FPLAIR(TAIR,PAIR) ! Mean free path of air molecules
-  write(*,*) 'MFPAIR',MFPAIR, TAIR, PAIR
+!  write(*,*) 'MFPAIR',MFPAIR, TAIR, PAIR
   Kn = MFPAIR/ri1             ! Knudsen number 
 
@@ -531 +531 @@
   real, intent(in) :: TAIR, PAIR, MRSA_loc
   real :: ROSAS, PSASAS_ZELE,  STSAS
-  write(*,*)'checkpoint change_wsa', ST, RHOSA, RHOsasat, PSASAS_ZELE(TAIR,WSA), TAIR, WSA
+!  write(*,*)'checkpoint change_wsa', ST, RHOSA, RHOsasat, PSASAS_ZELE(TAIR,WSA), TAIR, WSA
   ST = STSAS(TAIR,WSA)      ! Surface tension of sulfuric acid solution/vapor (N/m)
 
@@ -539 +539 @@
   SH2SO4 = (PAIR*MRSA_loc)/RHOsasat    ! Saturation Ratio
 
-  write(*,*)'checkpoint change_wsa', PSASAS_ZELE(TAIR,WSA)
+!  write(*,*)'checkpoint change_wsa', PSASAS_ZELE(TAIR,WSA)
 
 END SUBROUTINE change_wsa

trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus/new_cloud_venus.F

@@ -119 +119 @@
 !       Elle fait appel à la fct/ssrtine ITERWV()
 
-        CALL BRACWV(WVMIN,WVMAX,NBRAC,RMODE,
-     &  mrt_wv(ilon,ilev),mrt_sa(ilon,ilev),TT(ilon,ilev),
-     &  PP(ilon,ilev),FLAG,WSAFLAG,NBROOT)
+        CALL BRACWV(TT(ilon,ilev),PP(ilon,ilev),WVMIN,WVMAX,NBRAC,
+     &  RMODE,mrt_wv(ilon,ilev),mrt_sa(ilon,ilev),FLAG,WSAFLAG,NBROOT)
 
         SELECT CASE(FLAG)
@@ -130 +129 @@
 !       Elle fait appel ˆ la fct/ssrtine ITERWV()
 
-          WH2SO4(ilon,ilev)=IRFRMWV(WVMIN,WVMAX,WVACC,MAXITE,RMODE,
-     &    TT(ilon,ilev),PP(ilon,ilev),
+          WH2SO4(ilon,ilev)=IRFRMWV(TT(ilon,ilev),PP(ilon,ilev),
+     &    WVMIN,WVMAX,WVACC,MAXITE,RMODE,
      &    mrt_wv(ilon,ilev),mrt_sa(ilon,ilev),NBROOT)
 

trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus/wsa_new.F90

@@ -204 +204 @@
   WSAMIN = 0.1D0
   WSAMAX = 1.0D0
-  LPPWV = LOG(PAIR*WV)
+  LPPWV = DLOG(PAIR*WV)
 
   ! Appel Bracket de KEEQ         
@@ -293 +293 @@
         ! Si NROOT>1 on place la borne sup output ˆ la borne min du calcul en i             
         IF (NROOT.GT.1) THEN
-           XB=(1.-LOG(DBLE(I+1))/LOG(DBLE(N)))*XMAX
+           XB=(1.-DLOG(DBLE(I+1))/DLOG(DBLE(N)))*XMAX
         ENDIF
 
@@ -314 +314 @@
 
      DO J=N-1,1,-1
-        X=(1.-LOG(DBLE(N-J))/LOG(DBLE(N)))*XMAX
+        X=(1.-DLOG(DBLE(N-J))/DLOG(DBLE(N)))*XMAX
         !             write(*,*) 'BRACWV, bf 4th ITERWV (cas 2) '
         CALL ITERWV(TAIR,PAIR,X,WVLIQ,WVEQ,WVTOT,WSAOUT,SATOT,RADIUS)
@@ -365 +365 @@
   DX=(XB-XA)/N
   X=XB
-  FP=KEEQ(RADIUS,X,LPPWVINP)
+  FP=KEEQ(TAIR,RADIUS,X,LPPWVINP)
 
   DO I=N,1,-1
      X=X-DX
 
-     FC=KEEQ(RADIUS,X,LPPWVINP)
+     FC=KEEQ(TAIR,RADIUS,X,LPPWVINP)
 
      IF ((FP*FC).LE.0.) THEN
@@ -752 +752 @@
 
   ! Physical constants:
-  REAL :: MH2O, KEEQ
+  REAL :: KEEQ
 
   !     External functions needed:
@@ -764 +764 @@
   REAL :: C1
 
-  MH2O= MWV
-  C1=2.0D0*MH2O/RGAS
+  C1=2.0D0*MWV/RGAS
 
   KEEQ=LPPWV-C1*STSAS(TAIR,WX)/(TAIR*RADIUS*ROSAS(TAIR,WX))- &
@@ -818 +817 @@
 
   !pure acid satur vapor pressure
-  lpar= -11.695D0 + LOG(pstand) ! Zeleznik
-  acidps = 1/360.15D0 - 1.0/T + 0.38D0/545.D0*(1.0+LOG(360.15D0/T)-360.15D0/T)
+  lpar= -11.695D0 + DLOG(pstand) ! Zeleznik
+  acidps = 1/360.15D0 - 1.0/T + 0.38D0/545.D0*(1.0+DLOG(360.15D0/T)-360.15D0/T)
   acidps = 10156.D0*acidps + lpar
-  acidps = EXP(acidps)    !Pa
+  acidps = DEXP(acidps)    !Pa
 
   !acid sat.vap.PP over mixture (flat surface):