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

Timestamp:

Apr 6, 2017, 4:15:56 PM (8 years ago)

Author:

slebonnois

Message:

SL: bugs corrections outputs/clouds_AStol/upper_atmosphere/start2archive

Location:

trunk/LMDZ.VENUS/libf/phyvenus/cloudvenus

Files:

• new_cloud_sedim.F (modified) (3 diffs)
• new_cloud_venus.F (modified) (1 diff)

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

@@ -62 +62 @@
       real,parameter :: molrad=2.2e-10   ! CO2
       
-c     Cloud density (kg.m-3)
-c     ~~~~~~~~~~~~~~~~~~~~~~
-c      real, DIMENSION(n_lon,n_lev) ::  rho_droplet
+c     Ratio radius shell model du mode 3
+c     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+c       Ce ratio correspond aux mesures effectuées par J. Cimino (1982), Icarus
+c     Fixer ce parametre a 0 revient a une gouttelette pure en liquide acide sulfurique
+c     ATTENTION ! DOIT ETRE COHERENT AVEC new_cloud_venus !
+      REAL, PARAMETER :: qrad = 0.97
+      REAL :: qmass
+c       masse volumique du coeur (kg.m-3)
+c     ATTENTION ! DOIT ETRE COHERENT AVEC new_cloud_venus !
+      REAL, PARAMETER :: rho_core = 2500.0
 
       REAL, DIMENSION(n_lon,n_lev+1) ::
@@ -131 +138 @@
 c           pbndlay(:,51)=0 (en parallèle c'est sûr), ne pas l'utiliser pour Fse
         
-        DO imode=1, nbr_mode
+c     Sedimentation pour une gouttelette mode 3 de type J. Cimino, 1982, Icarus
+c     c.a.d 97% radius due a solide 3% radius acide sulfurique  
+        DO imode=1, nbr_mode - 1
          DO l = cloudmin, cloudmax
             DO ig=1,n_lon
@@ -238 +247 @@
       ENDDO
 
+c****************************************************************
+c        On calcule le F_sed du mode 3 + coeff*(Fsed1 + Fsed2)
+c****************************************************************
+         DO l = cloudmin, cloudmax
+            DO ig=1,n_lon
+
+c       calcul de qmass
+            qmass=(rho_core*qrad**3)/
+     &    (rho_core*qrad**3+rho_droplet(ig,l)*(1.-qrad**3))
+
+c     RD=1000.*RNAVO*RKBOL/RMD avec RMD=43.44 Masse molaire atm venus en g.mol-1    
+           D_stokes=(((qmass*rho_core+(1.-qmass)*rho_droplet(ig,l))
+     &      -pmidlay(ig,l)/(RD*pt(ig,l))))
+     &  *(2./9.)*(RG/VISCOSITY_CO2(pt(ig,l)))
+      
+           l_mean=(pt(ig,l)/pmidlay(ig,l))*
+     &  (0.707*R/(4.*RPI* molrad*molrad * RNAVO))
+      
+           R_mode0=R_MEDIAN(ig,l,3)*
+     &     EXP(-LOG(STDDEV(ig,l,3))**2.)
+              IF ((l_mean/(R_mode0)).GT.10.) THEN
+              Rp_DL=R_MEDIAN(ig,l,3)*
+     &        EXP(3.*LOG(STDDEV(ig,l,3))**2.)
+              ELSE
+              Rp_DL=R_MEDIAN(ig,l,3)*
+     &        EXP(4.*LOG(STDDEV(ig,l,3))**2.)
+              ENDIF
+                
+           a=1.246*l_mean
+        
+           c=0.87/l_mean
+        
+           b_exp=0.42*l_mean*EXP(-c*Rp_DL)
+        
+           A1=a+b_exp*(1.+c*Rp_DL
+     &  +0.5*(Rp_DL*c)**2
+     &  +1./6.*(Rp_DL*c)**3)
+     
+           A2=1.-b_exp*(c
+     &  +Rp_DL*c**2
+     &  +0.5*(Rp_DL**2)*(c**3))
+        
+           A3=0.5*b_exp*(c**2+Rp_DL*c**3)
+        
+           A4=-b_exp*1./6.*c**3
+
+c     Addition des Flux de tous les modes presents      
+       F_sed(ig,l)=F_sed(ig,l)
+     &  +((1.-qmass)/(1.-qmass*K_MASS(ig,l,3)))*(
+     &  (qmass*rho_core+(1.-qmass)*rho_droplet(ig,l))*4./3.*RPI*
+     &  NBRTOT(ig,l,3)*1.0E6*D_stokes*(
+     &  A1*R_MEDIAN(ig,l,3)**4
+     &  *EXP(8.0*LOG(STDDEV(ig,l,3))**2.)
+     &  +A2*R_MEDIAN(ig,l,3)**5
+     &  *EXP(12.5*LOG(STDDEV(ig,l,3))**2.)
+     &  +A3*R_MEDIAN(ig,l,3)**6
+     &  *EXP(18.0*LOG(STDDEV(ig,l,3))**2.)
+     &  +A4*R_MEDIAN(ig,l,3)**7
+     &  *EXP(24.5*LOG(STDDEV(ig,l,3))**2.)))
+      
+c      PRINT*,' APRES dTime: F_sed=',F_sed(ig,l), ig, l
+      
+        IF (F_sed(ig,l).GT.m_lay(ig,l)) THEN
+        PRINT*,'==============================================='
+        PRINT*,'WARNING On a epuise la couche', ig, l
+        PRINT*,'On epuise pas une couche avec une espèce
+     &   minoritaire, c est pas bien maaaaaal'
+            PRINT*,'Water',zqi_wv(ig,l),'Sulfuric Acid',zqi_sa(ig,l)
+        PRINT*,'F_sed:',F_sed(ig,l),'m_lay:',m_lay(ig,l)
+        PRINT*,'F_sed/dtphy',F_sed(ig,l)/ptimestep
+        PRINT*,'Pbnd top',pbndlay(ig,l+1),'Temp',pt(ig,l),'Rho',
+     &  rho_droplet(ig,l)
+                PRINT*,'Ntot',NBRTOT(ig,l,:)
+                PRINT*,'StdDev',STDDEV(ig,l,:),'Rmed',R_MEDIAN(ig,l,:)
+                PRINT*,'K_MASS',K_MASS(ig,l,:)
+                PRINT*,'WSA',WH2SO4(ig,l),'RHO',rho_droplet(ig,l)
+                
+c               ELSE
+c               
+c               PRINT*,'~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
+c       PRINT*,'WARNING On a PAS epuise la couche', ig, l
+c       PRINT*,'F_sed:',F_sed(ig,l),'m_lay:',m_lay(ig,l)
+c       PRINT*,'F_sed/dtphy',F_sed(ig,l)/ptimestep
+c       PRINT*,'Pbnd top',pbndlay(ig,l+1),'Temp',pt(ig,l),'Rho',
+c     &         rho_droplet(ig,l)(ig,l)
+c               PRINT*,'Ntot',NBRTOT(ig,l),'Ntot m3',NBRTOT(ig,l)*1.0e6
+c               PRINT*,'StdDev',STDDEV(ig,l),'Rmed',R_MEDIAN(ig,l) 
+            STOP                
+              ENDIF
+        
+           IF (F_sed(ig,l).LT.0.0d0) THEN
+              PRINT*,"F_sed est négatif !!!"
+              PRINT*,'F_sed:',F_sed(ig,l),'m_lay:',m_lay(ig,l)
+        PRINT*,'F_sed/dtphy',F_sed(ig,l)/ptimestep
+        PRINT*,'Pbnd top',pbndlay(ig,l+1),'Pmid',pmidlay(ig,l)
+        PRINT*,'Temp',pt(ig,l),'Rho',
+     &  rho_droplet(ig,l)
+                PRINT*,'Ntot',NBRTOT(ig,l,imode),'Ntot m3',
+     &  NBRTOT(ig,l,imode)*1.0e6
+                PRINT*,'StdDev',STDDEV(ig,l,imode),'Rmed',
+     &          R_MEDIAN(ig,l,imode)
+                PRINT*,'A1',A1,'A2',A2
+                PRINT*,'A3',A1,'A4',A2
+                PRINT*,'D_stokes',D_stokes
+                STOP
+           ENDIF
+           
+              ENDDO
+              
+c           ELSE           
+c           F_sed(:,l)=0.0d0           
+c           ENDIF
+           
+         ENDDO
+c****************************************************************
+
 c     Passage du Flux au Flux pour un pas de temps (== kg.m-2)     
       F_sed(:,:)=F_sed(:,:)*ptimestep

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

@@ -62 +62 @@
 c       Ce ratio correspond aux mesures effectuées par J. Cimino (1982), Icarus
 c     Fixer ce parametre a 0 revient a une gouttelette pure en liquide acide sulfurique
+c     ATTENTION ! DOIT ETRE COHERENT AVEC new_cloud_sedim !
       REAL, PARAMETER :: qrad = 0.97
       REAL :: qmass
 c       masse volumique du coeur (kg.m-3)
+c     ATTENTION ! DOIT ETRE COHERENT AVEC new_cloud_sedim !
       REAL, PARAMETER :: rho_core = 2500.0
 !----------------------------------------------------------------------------