source: trunk/LMDZ.MARS/libf/phymars/nuclea.F @ 358

Last change on this file since 358 was 358, checked in by aslmd, 13 years ago

LMDZ.MARS:

A FIRST VERSION WITH SCAVENGING
TRANSPARENT TO CASUAL USER
FOR WATER CYCLE FOLKS, SEE BELOW

[TN and AS solved conflicts due to r330 >> r357]

07/11/11 == JBM

Changed watercloud.F to call two separate routines,

simpleclouds.F or improvedclouds.F, which are a simplified and
full microphysical scheme for cloud formation, respectively.
Removed the tag called "improved" in watercloud.F, and added
another tag called "microphys" which is defined in callphys.F
instead. Changed routines: callkeys, inifis, physiq, watercloud.

Reimplemented the use of typical profiles of dust particle sizes

and CCNs in simpleclouds.F. Corrected the previously used
analytical CCN profile. Moved ccn_factor to simpleclouds.F,
which won't be used in the new microphysical scheme. Changed
routines: aeropacity, physiq, simpleclouds, watercloud.

Computed rdust at the end of callsedim instead of updatereffrad,

except at firstcall. Renamed rfall into rsedcloud and moved it
in simpleclouds. Moved nuice_sed to initracer.F and added it to
tracer.h. Changed routines: callsedim, physiq, tracer.h,
watercloud, initracer, simpleclouds, updatereffrad.

Added two tracers for the CCNs. Added the different tests in

initracer.F to make sure that, for example, microphys is not
used without doubleq, etc. Corrected an inconsistency in
callsedim.F, and changed the way r0 is updated. Changes
routines: callsedim, inifis, initracer, physiq, testphys1d,
tracer.h.

Added the ability to have a spatially variable density in

newsedim (same method as that used for the radius of
sedimentation). Required the addition of one input to newsedim,
which is the size of the density variable "rho". Changed
routines: callsedim, newsedim.

Added an output to aeropacity called "tauscaling", which is a

factor to convert dust_mass and dust_number to true quantities,
based on tauvis. Removed ccn and qdust from aeropacity, which
are obsolete now.

Wrote improvedcloud.F which includes all the microphysical

scheme, and connected it to the sedimentation scheme. Added and
changed routines: callsedim, physiq, growthrate, nuclea,
improvedclouds, initracer, watercloud, microphys.h.

07/11/11 == TN

Added CCN & dust tracers updates in physiq.F

Corrected a bug that can give negative CCN numbers, due to the
use of single precision values (only 6 decimals) whereas up to 10E+10
CCN can be formed or disappears...

Corrected physical bug that causes h2o_vap or h2o_ice

to be negative in improvedclouds.F.

Corrected physical bug that causes CCN & dust tracers

to be negative in improvedclouds.F when all ice sublimates and
releases dust

Added parameter contact mteta in callphys.def

Default value is still 0.95, see inifis.F

Changed tendancy computation for dust_number in improvedclouds.F

that was not the right one. Indeed, the scavenged dust_number tracer
is computed from the dust_mass one, and its tendancy before scavenging
must be taken into account to compute its scavenging's tendancy.
File size: 3.7 KB
Line 
1*******************************************************
2*                                                     *
3      subroutine nuclea(ph2o,temp,sat,n_ccn,nucrate)
4      implicit none
5*                                                     *
6*   This subroutine computes the nucleation rate      *
7*   as given in Pruppacher & Klett (1978) in the      *
8*   case of water ice forming on a solid substrate.   *
9*     Definition refined by Keese (jgr,1989)          *
10*   Authors: F. Montmessin                            *
11*     Adapted for the LMD/GCM by J.-B. Madeleine      *
12*     (October 2011)                                  *
13*******************************************************
14
15#include "dimensions.h"
16#include "dimphys.h"
17#include "comcstfi.h"
18#include "tracer.h"
19#include "microphys.h"
20
21c     Inputs
22      DOUBLE PRECISION ph2o,sat
23      DOUBLE PRECISION n_ccn(nbin_cld)
24      REAL temp
25
26c     Output
27      DOUBLE PRECISION nucrate(nbin_cld)
28
29c     Local variables
30      DOUBLE PRECISION nh2o
31      DOUBLE PRECISION sig      ! Water-ice/air surface tension  (N.m)
32      external sig
33      DOUBLE PRECISION rstar    ! Radius of the critical germ (m)
34      DOUBLE PRECISION gstar    ! # of molecules forming a critical embryo
35      DOUBLE PRECISION fistar   ! Activation energy required to form a critical embryo (J)
36      DOUBLE PRECISION zeldov   ! Zeldovitch factor (no dim)
37      DOUBLE PRECISION fshape   ! function defined at the end of the file
38      DOUBLE PRECISION deltaf
39
40c     Ratio rstar/radius of the nucleating dust particle
41      double precision xratio
42     
43      double precision mtetalocal ! local mteta in double precision
44
45      integer i
46
47c     *************************************************
48
49      mtetalocal = mteta
50
51      if (sat .gt. 1.) then    ! minimum condition to activate nucleation
52
53        nh2o   = ph2o / kbz / temp
54        rstar  = 2. * sig(temp) * vo1 / (rgp*temp*dlog(sat))
55        gstar  = 4. * nav * pi * (rstar**3) / (3.*vo1)
56
57c       Loop over size bins
58        do 200 i=1,nbin_cld
59
60          if ( n_ccn(i) .lt. 1e-10 ) then
61c           no dust, no need to compute nucleation!
62            nucrate(i)=0.
63            goto 200
64          endif
65
66          xratio   = rad_cld(i) / rstar
67          fistar = (4./3.*pi) * sig(temp) * (rstar**2.) *
68     &             fshape(mtetalocal,xratio)
69          deltaf = (2.*desorp-surfdif-fistar)/
70     &             (kbz*temp)
71          deltaf = min( max(deltaf, -100.), 100.)
72
73          if (deltaf.eq.-100.) then
74            nucrate(i) = 0.
75          else
76            zeldov = sqrt ( fistar /
77     &               (3.*pi*kbz*temp*(gstar**2.)) )
78            nucrate(i)= zeldov * kbz * temp * rstar
79     &                  * rstar * 4. * pi
80     &                  * ( nh2o*rad_cld(i) )**2.
81     &                  / ( fshape(mtetalocal,xratio) * nus * m0 )
82     &                  * dexp (deltaf)
83          endif
84
85200     continue
86
87      else
88
89        do i=1,nbin_cld
90          nucrate(i) = 0.
91        enddo
92
93      endif
94
95      return
96      end
97
98*********************************************************
99      double precision function fshape(cost,rap)
100      implicit none
101*        function computing the f(m,x) factor           *
102* related to energy required to form a critical embryo  *
103*********************************************************
104
105      double precision cost,rap
106      double precision phi
107      double precision a,b,c
108
109      phi = sqrt( 1. - 2.*cost*rap + rap**2 )
110      a = 1. + ( (1.-cost*rap)/phi )**3
111      b = (rap**3) * (2.-3.*(rap-cost)/phi+((rap-cost)/phi)**3)
112      c = 3. * cost * (rap**2) * ((rap-cost)/phi-1.)
113
114      fshape = 0.5*(a+b+c)
115
116      if (rap.gt.3000.) fshape = ((2.+cost)*(1.-cost)**2)/4.
117
118      return
119      end
Note: See TracBrowser for help on using the repository browser.