source: trunk/LMDZ.MARS/libf/phymars/growthrate.F

Last change on this file was 3008, checked in by emillour, 17 months ago

Mars PCM:
Some code cleanup around microphysics. Turn microphys.h into module
microphys_h.F90, and while at it also turn nuclea.F, growthrate.F90 and
massflowrateco2.F90 into modules.
EM

File size: 3.8 KB
RevLine 
[3008]1      MODULE growthrate_mod
2     
3      IMPLICIT NONE
4     
5      CONTAINS
6
[2407]7      subroutine growthrate(temp,pmid,psat,rcrystal,res,Dv)
[38]8
[1036]9      use tracer_mod, only: rho_ice
[3008]10      use comcstfi_h, only: pi
11      use microphys_h, only: kbz, mh2o, mco2, molco2, molh2o
12      use microphys_h, only: nav, rgp, To
13     
[38]14      IMPLICIT NONE
15
16c=======================================================================
17c
18c     Determination of the water ice crystal growth rate
19c
[358]20c     Authors: F. Montmessin
21c       Adapted for the LMD/GCM by J.-B. Madeleine (October 2011)
[626]22c       Use of resistances in the analytical function
23c            instead of growth rate - T. Navarro (2012)
[358]24c     
[38]25c=======================================================================
26
27c-----------------------------------------------------------------------
28c   declarations:
29c   -------------
30
31c
32c   arguments:
33c   ----------
34
[626]35c     Input
[3008]36      REAL, INTENT(IN) :: temp     ! temperature in the middle of the layer (K)
37      REAL, INTENT(IN) :: pmid     ! pressure in the middle of the layer (K)
38      REAL, INTENT(IN) :: psat   ! water vapor saturation pressure (Pa)
39      REAL, INTENT(IN) :: rcrystal ! crystal radius before condensation (m)
[38]40
[626]41c     Output
[3008]42      REAL, INTENT(OUT) :: res      ! growth resistance (res=Rk+Rd)
43      REAL, INTENT(OUT) :: Dv       ! water vapor diffusion coefficient
[626]44
[38]45c   local:
46c   ------
47
[633]48      REAL k,Lv                 
49      REAL knudsen           ! Knudsen number (gas mean free path/particle radius)
[2407]50      REAL afactor,lambda       ! Intermediate computations for growth rate
[633]51      REAL Rk,Rd
[626]52     
53     
[38]54
55c-----------------------------------------------------------------------
56c      Ice particle growth rate by diffusion/impegement of water molecules
57c                r.dr/dt = (S-Seq) / (Seq*Rk+Rd)
58c        with r the crystal radius, Rk and Rd the resistances due to
59c        latent heat release and to vapor diffusion respectively
60c-----------------------------------------------------------------------
61
62c     - Equilibrium saturation accounting for KeLvin Effect
[358]63c      seq=exp(2*sigh2o*mh2o/(rho_ice*rgp*t*r))
64c      (already computed in improvedcloud.F)
[38]65
66c     - Thermal conductibility of CO2
[358]67      k  = (0.17913 * temp - 13.9789) * 4.184e-4
[38]68c     - Latent heat of h2o (J.kg-1)
[530]69      Lv = (2834.3
70     &        - 0.28  * (temp-To)
71     &        - 0.004 * (temp-To) * (temp-To) ) * 1.e+3
[38]72
73c     - Constant to compute gas mean free path
74c     l= (T/P)*a, with a = (  0.707*8.31/(4*pi*molrad**2 * avogadro))
[530]75      afactor = 0.707*rgp/(4 * pi * molco2 * molco2 * nav)
[38]76
77c     - Compute Dv, water vapor diffusion coefficient
78c       accounting for both kinetic and continuum regime of diffusion,
79c       the nature of which depending on the Knudsen number.
80
[358]81      Dv = 1./3. * sqrt( 8*kbz*temp/(pi*mh2o/nav) )* kbz * temp /
[530]82     &   ( pi * pmid * (molco2+molh2o)*(molco2+molh2o)
83     &        * sqrt(1.+mh2o/mco2) )
[626]84     
[358]85      knudsen = temp / pmid * afactor / rcrystal
[38]86      lambda  = (1.333+0.71/knudsen) / (1.+1./knudsen)
[626]87     
88c      Dv is not corrected. Instead, we use below coefficients coeff1, coeff2
89c      Dv      = Dv / (1. + lambda * knudsen)
[38]90
91c     - Compute Rk
[633]92      Rk = Lv*Lv* rho_ice * mh2o / (k*rgp*temp*temp)
[38]93c     - Compute Rd
[358]94      Rd = rgp * temp *rho_ice / (Dv*psat*mh2o)
[626]95     
96     
[633]97      res = Rk + Rd*(1. + lambda * knudsen)
[626]98     
[633]99      !coeff1 = real(Rk + Rd*(1. + lambda * knudsen))
100      !coeff2 = real(Rk + Rd*(1. - lambda * knudsen))
101     
[626]102c Below are growth rate used for other schemes
[358]103c     - Compute growth=rdr/dt, then r(t+1)= sqrt(r(t)**2.+2.*growth*dt)
[626]104c      growth = 1. / (seq*Rk+Rd)
105c      growth = (ph2o/psat-seq) / (seq*Rk+Rd)
[358]106c      rf   = sqrt( max( r**2.+2.*growth*timestep , 0. ) )
[38]107c      dr   = rf-r
108
[3008]109      end subroutine growthrate
110     
111      END MODULE growthrate_mod
[38]112
Note: See TracBrowser for help on using the repository browser.