Changeset 3536 for LMDZ6/trunk/libf


Ignore:
Timestamp:
Jun 13, 2019, 7:19:07 PM (6 years ago)
Author:
oboucher
Message:

Setting up the new nucleation routine (compiles but not yet tested)

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/StratAer/nucleation_tstep_mod.F90

    r3534 r3536  
    66CONTAINS
    77
    8 SUBROUTINE nucleation_rate(rhoa,t_seri,pplay,rh,a_xm,b_xm,c_xm,nucl_rate,ntot,x)
     8SUBROUTINE nucleation_rate(rhoa,t_seri,pplay,rh,a_xm,b_xm,c_xm,nucl_rate,ntot_n,x_n)
    99
    1010  USE aerophys
    1111  USE infotrac
    12   USE YOMCST, ONLY : RPI, RD, RKBOL
     12  USE YOMCST, ONLY : RPI, RD, RMD, RKBOL, RNAVO
    1313
    1414  IMPLICIT NONE
    1515
    16   ! input variable
     16  ! input variables
    1717  LOGICAL, PARAMETER :: flag_new_nucl=.TRUE.
    18   REAL rhoa !H2SO4 number density [molecules/cm3]
    19   REAL t_seri
    20   REAL pplay
    21   REAL rh
     18  REAL rhoa    ! H2SO4 number density [molecules/cm3]
     19  REAL t_seri  ! temperature (K)
     20  REAL pplay   ! pressure (Pa)
     21  REAL rh      ! relative humidity (between 0 and 1)
    2222  REAL a_xm, b_xm, c_xm
    2323
    24   ! output variable
     24  ! output variables
    2525  REAL nucl_rate
    26   REAL ntot ! total number of molecules in the critical cluster
    27   REAL x    ! mole fraction of H2SO4 in the critical cluster
    28 
    29   ! local variable
    30   REAL                                          :: jnuc !nucleation rate in 1/cm3s (10^-7-10^10 1/cm3s)
    31   REAL                                          :: rc   !radius of the critical cluster in nm
     26  REAL ntot_n ! total number of molecules in the critical cluster for neutral nucleation
     27  REAL x_n    ! mole fraction of H2SO4 in the critical cluster for neutral nucleation
     28
     29  ! local variables
     30  REAL jnuc_n ! nucleation rate in 1/cm3s (10^-7-10^10 1/cm3s) for neutral nucleation
     31  REAL rc_n   ! radius of the critical cluster in nm for neutral nucleation
     32  REAL na_n   ! sulfuric acid molecules in the neutral critical cluster (NOT IN USE)
     33  LOGICAL kinetic_n ! true if kinetic neutral nucleation (NOT IN USE)
     34  LOGICAL kinetic_i ! true if kinetic ion-induced nucleation (NOT IN USE)
     35  REAL rhoatres ! threshold concentration of H2SO4 (1/cm^3) for neutral kinetic nucleation (NOT IN USE)
    3236  REAL VH2SO4mol
     37  REAL ntot_i, x_i, jnuc_i, rc_i, na_i, n_i ! quantities for charged nucleation (NOT IN USE)
     38  REAL csi     ! Ion condensation sink (s-1) NOT IN USE
     39  REAL airn    ! Air molecule concentration in (cm-3) NOT IN USE
     40  REAL ipr     ! Ion pair production rate (cm-3 s-1) NOT IN USE
    3341
    3442  ! call nucleation routine
    35 !  IF (.NOT.flag_new_nucl) THEN
    36     CALL binapara(t_seri,rh,rhoa,jnuc,x,ntot,rc)
    37 !  ELSE
    38 !    CALL newbinapara(t_seri,rh,rhoa,jnuc,x,ntot,rc)
    39 !  ENDIF
    40 
    41   IF (ntot < 4.0) THEN
    42     !set jnuc to collision rate of two H2SO4 molecules (following personal communication of Ulrike Niemeier and Hanna Vehkamäki)
    43     VH2SO4mol=mH2SO4mol/(1.E-3*(a_xm+t_seri*(b_xm+t_seri*c_xm))) !cm3
    44     jnuc = rhoa**2. *(3./4.*RPI)**(1./6.) *(12.*RKBOL*t_seri/mH2SO4mol)**0.5 &
    45          & *100.*(2.*VH2SO4mol**(1./3.))**2. !1/(cm3s)
    46     ntot=2.0
    47     x=1.0
    48   ENDIF
    49 
    50   ! convert jnuc from particles/cm3/s to kg(H2SO4)/kgA/s
    51   nucl_rate=jnuc*ntot*x*mH2SO4mol/(pplay/t_seri/RD/1.E6)
     43  IF (.NOT.flag_new_nucl) THEN
     44    ! Use older routine from Hanna Vehkamäki (FMI)
     45    CALL binapara(t_seri,rh,rhoa,jnuc_n,x_n,ntot_n,rc_n)
     46    ! when total number of molecules is too small
     47    ! then set jnuc_n to collision rate of two H2SO4 molecules (following personal communication of Ulrike Niemeier and Hanna Vehkamäki)
     48    IF (ntot_n < 4.0) THEN
     49      VH2SO4mol=mH2SO4mol/(1.E-3*(a_xm+t_seri*(b_xm+t_seri*c_xm))) !cm3
     50      jnuc_n = rhoa**2. *(3./4.*RPI)**(1./6.) *(12.*RKBOL*t_seri/mH2SO4mol)**0.5 &
     51           & *100.*(2.*VH2SO4mol**(1./3.))**2. !1/(cm3s)
     52      ntot_n=2.0
     53      x_n=1.0
     54    ENDIF
     55  ELSE
     56    ! Use new routine from Anni Maattanen (LATMOS)
     57    csi=0.0   ! no charged nucleation for now
     58    ipr=-1.0  ! dummy value to make sure charged nucleation does not occur
     59    airn=0.0  ! NOT IN USE
     60!   airn=pplay/t_seri/RD/1.E3*RNAVO/RMD ! molec cm-3 (for future use, to be confirmed)
     61    CALL newbinapara(t_seri,rh,rhoa,csi,airn,ipr,jnuc_n,ntot_n,jnuc_i,ntot_i, &
     62                   & x_n,x_i,na_n,na_i,rc_n,rc_i,n_i,kinetic_n,kinetic_i,rhoatres)
     63  ENDIF
     64
     65  ! convert jnuc_n from particles/cm3/s to kg(H2SO4)/kgA/s
     66  nucl_rate=jnuc_n*ntot_n*x_n*mH2SO4mol/(pplay/t_seri/RD/1.E6)
    5267
    5368END SUBROUTINE nucleation_rate
    … …  
    155170
    156171  IF (t < 190.15) THEN
    157 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): temperature < 190.15 K, using 190.15 K (T=',t,'K)'
     172!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): temperature < 190.15 K, using 190.15 K (T=',t,'K)'
    158173     t=190.15
    159174  ENDIF
    160175
    161176  IF (t > 300.15) THEN
    162 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): temperature > 300.15 K, using 300.15 K'
     177!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): temperature > 300.15 K, using 300.15 K'
    163178     t=300.15
    164179  ENDIF
    165180
    166181  IF (rh < 0.0001) THEN
    167 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): saturation ratio of water < 0.0001, using 0.0001'
     182!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): saturation ratio of water < 0.0001, using 0.0001'
    168183     rh=0.0001
    169184  ENDIF
    170185
    171186  IF (rh > 1.0) THEN
    172 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): saturation ratio of water > 1 using 1'
     187!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): saturation ratio of water > 1 using 1'
    173188!     print *, 'rh=',rh
    174189     rh=1.0
    … …  
    176191
    177192  IF (rhoa < 1.e4) THEN
    178 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): sulfuric acid concentration < 1e4 1/cm3, using 1e4 1/cm3'
     193!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): sulfuric acid concentration < 1e4 1/cm3, using 1e4 1/cm3'
    179194     rhoa=1.e4
    180195  ENDIF
    181196
    182197  IF (rhoa > 1.e11) THEN
    183 !     print *,'Warning (ilon=',ilon,'ilat=',ilat,'): sulfuric acid concentration > 1e11 1/cm3, using 1e11 1/cm3'
     198!     print *,'Warning (ilon=',ilon,'ilev=',ilev,'): sulfuric acid concentration > 1e11 1/cm3, using 1e11 1/cm3'
    184199     rhoa=1.e11
    185200  ENDIF
    … …  
    266281
    267282  IF (jnuc < 1.E-7) THEN
    268 !     print *,'Warning (ilon=',ilon,'ilat=',ilat'): nucleation rate < 1E-7/cm3s, using 0.0/cm3s,'
     283!     print *,'Warning (ilon=',ilon,'ilev=',ilev'): nucleation rate < 1E-7/cm3s, using 0.0/cm3s,'
    269284     jnuc=0.0
    270285  ENDIF
    271286
    272287  IF (jnuc > 1.e10) THEN
    273 !     print *,'Warning (ilon=',ilon,'ilat=',ilat, &
     288!     print *,'Warning (ilon=',ilon,'ilev=',ilev, &
    274289!        & '): nucleation rate > 1e10/cm3s, using 1e10/cm3s'
    275290     jnuc=1.e10
    … …  
    288303!---------------------------------------------------------------------------------------------------
    289304
    290 SUBROUTINE newbinapara(t,satrat,rhoa,csi,airn,ipr,jnuc_n,ntot_n,jnuc_i,ntot_i,        &
    291                    &   x_n,x_i,na_n,na_i,rc_n,rc_i,n_i,kinetic_n,kinetic_i,rhoatres)
     305SUBROUTINE newbinapara(t,satrat,rhoa,csi,airn,ipr,jnuc_n_real,ntot_n_real,jnuc_i_real,ntot_i_real,        &
     306                   &   x_n_real,x_i_real,na_n_real,na_i_real,rc_n_real,rc_i_real,n_i_real,                &
     307                   &   kinetic_n,kinetic_i,rhoatres_real)
    292308
    293309  !    Fortran 90 subroutine newbinapara
    … …  
    330346 
    331347  !Global intent in
    332   DOUBLE PRECISION,INTENT(IN) :: t         ! temperature in K
    333   DOUBLE PRECISION,INTENT(IN) :: satrat    ! saturatio ratio of water (between zero and 1)
    334   DOUBLE PRECISION,INTENT(IN) :: rhoa      ! sulfuric acid concentration in 1/cm3
    335   DOUBLE PRECISION,INTENT(IN) :: csi       ! Ion condensation sink (s-1)
    336   DOUBLE PRECISION,INTENT(IN) :: airn      ! Air molecule concentration in (cm-3)
    337   DOUBLE PRECISION,INTENT(IN) :: ipr       ! Ion pair production rate (cm-3 s-1)
     348  REAL,INTENT(IN) :: t         ! temperature in K
     349  REAL,INTENT(IN) :: satrat    ! saturatio ratio of water (between zero and 1)
     350  REAL,INTENT(IN) :: rhoa      ! sulfuric acid concentration in 1/cm3
     351  REAL,INTENT(IN) :: csi       ! Ion condensation sink (s-1)
     352  REAL,INTENT(IN) :: airn      ! Air molecule concentration in (cm-3)
     353  REAL,INTENT(IN) :: ipr       ! Ion pair production rate (cm-3 s-1)
    338354  !Global intent out
    339   DOUBLE PRECISION,INTENT(OUT) :: jnuc_n   ! Neutral nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
    340   DOUBLE PRECISION,INTENT(OUT) :: ntot_n   ! total number of molecules in the neutral critical cluster
    341   DOUBLE PRECISION,INTENT(OUT) :: jnuc_i   ! Charged nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
    342   DOUBLE PRECISION,INTENT(OUT) :: ntot_i   ! total number of molecules in the charged critical cluster
    343   DOUBLE PRECISION,INTENT(OUT) :: x_n      ! mole fraction of H2SO4 in the neutral critical cluster
    344   DOUBLE PRECISION,INTENT(OUT) :: x_i      ! mole fraction of H2SO4 in the charged critical cluster
     355  REAL,INTENT(OUT) :: jnuc_n_real   ! Neutral nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
     356  REAL,INTENT(OUT) :: ntot_n_real   ! total number of molecules in the neutral critical cluster
     357  REAL,INTENT(OUT) :: jnuc_i_real   ! Charged nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
     358  REAL,INTENT(OUT) :: ntot_i_real   ! total number of molecules in the charged critical cluster
     359  REAL,INTENT(OUT) :: x_n_real      ! mole fraction of H2SO4 in the neutral critical cluster
     360  REAL,INTENT(OUT) :: x_i_real      ! mole fraction of H2SO4 in the charged critical cluster
    345361                                           ! (note that x_n=x_i in nucleation regime)
    346   DOUBLE PRECISION,INTENT(OUT) :: na_n     ! sulfuric acid molecules in the neutral critical cluster
    347   DOUBLE PRECISION,INTENT(OUT) :: na_i     ! sulfuric molecules in the charged critical cluster
    348   DOUBLE PRECISION,INTENT(OUT) :: rc_n     ! radius of the charged critical cluster in nm
    349   DOUBLE PRECISION,INTENT(OUT) :: rc_i     ! radius of the charged critical cluster in nm
    350   DOUBLE PRECISION,INTENT(OUT) :: n_i      ! number of ion pairs in air (cm-3)
     362  REAL,INTENT(OUT) :: na_n_real     ! sulfuric acid molecules in the neutral critical cluster
     363  REAL,INTENT(OUT) :: na_i_real     ! sulfuric molecules in the charged critical cluster
     364  REAL,INTENT(OUT) :: rc_n_real     ! radius of the charged critical cluster in nm
     365  REAL,INTENT(OUT) :: rc_i_real     ! radius of the charged critical cluster in nm
     366  REAL,INTENT(OUT) :: n_i_real      ! number of ion pairs in air (cm-3)
     367  REAL,INTENT(OUT) :: rhoatres_real ! threshold concentration of H2SO4 (1/cm^3) for neutral kinetic nucleation
    351368  LOGICAL,INTENT(OUT)  :: kinetic_n        ! true if kinetic neutral nucleation
    352369  LOGICAL,INTENT(OUT)  :: kinetic_i        ! true if kinetic ion-induced nucleation
    353   DOUBLE PRECISION,INTENT(OUT) :: rhoatres ! threshold concentration of H2SO4 (1/cm^3) for neutral kinetic nucleation
    354370
    355371  ! Local
     372  DOUBLE PRECISION :: jnuc_n      ! Neutral nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
     373  DOUBLE PRECISION :: ntot_n      ! total number of molecules in the neutral critical cluster
     374  DOUBLE PRECISION :: jnuc_i      ! Charged nucleation rate in 1/cm3s (J>10^-7 1/cm3s)
     375  DOUBLE PRECISION :: ntot_i      ! total number of molecules in the charged critical cluster
     376  DOUBLE PRECISION :: x_n         ! mole fraction of H2SO4 in the neutral critical cluster
     377  DOUBLE PRECISION :: x_i         ! mole fraction of H2SO4 in the charged critical cluster
     378                                              ! (note that x_n=x_i in nucleation regime)
     379  DOUBLE PRECISION :: na_n        ! sulfuric acid molecules in the neutral critical cluster
     380  DOUBLE PRECISION :: na_i        ! sulfuric molecules in the charged critical cluster
     381  DOUBLE PRECISION :: rc_n        ! radius of the charged critical cluster in nm
     382  DOUBLE PRECISION :: rc_i        ! radius of the charged critical cluster in nm
     383  DOUBLE PRECISION :: n_i         ! number of ion pairs in air (cm-3)
     384  DOUBLE PRECISION :: rhoatres    ! threshold concentration of H2SO4 (1/cm^3) for neutral kinetic nucleation
    356385  DOUBLE PRECISION :: x           ! mole fraction of H2SO4 in the critical cluster
    357386  DOUBLE PRECISION :: satratln    ! bounded water saturation ratio for neutral case (between 5.E-6 - 1.0)
    … …  
    812841
    813842  ENDIF
     843
     844!--conversion from double precision to float in case the model is run in single precision
     845  jnuc_n_real   = jnuc_n
     846  ntot_n_real   = ntot_n
     847  jnuc_i_real   = jnuc_i
     848  ntot_i_real   = ntot_i
     849  x_n_real      = x_n
     850  x_i_real      = x_i
     851  na_n_real     = na_n
     852  na_i_real     = na_i
     853  rc_n_real     = rc_n
     854  rc_i_real     = rc_i
     855  n_i_real      = n_i
     856  rhoatres_real = rhoatres
    814857 
    815858END SUBROUTINE newbinapara
Note: See TracChangeset for help on using the changeset viewer.