source: trunk/LMDZ.MARS/libf/aeronomars/deposition.F

Last change on this file was 3466, checked in by emillour, 2 months ago

Mars PCM:
More tidying in aeronomars:

  • remove unused "inv.F" and remove "dtridgl.F" which is not used here and is a duplicate of the "dtridgl" routine in phymars/swr_toon.F
  • turn aeronomars routines to modules, for those which aren't in modules yet.

EM

File size: 6.0 KB
RevLine 
[3466]1      module deposition_mod
2     
3      implicit none
4     
5      contains
6     
[1047]7      subroutine deposition(ngrid, nlayer, nq,
8     &                      ig, ig_vl1, pplay, pplev, zzlay, zzlev,
[415]9     $                      zu, zv, zt, zycol, ptimestep, co2ice)
10cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
11c
12c     dry deposition of chemical species
13c
14cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
15c
[1047]16      use surfdat_h, only: z0 ! surface roughness
17      use conc_mod, only: rnew ! specific gas constant
[415]18      implicit none
19c
20c
21c     input
22c
[1047]23      integer,intent(in) :: ngrid ! number of atmospheric columns
24      integer,intent(in) :: nlayer ! number of atmospheric layers
25      integer,intent(in) :: nq ! number of tracers
26      integer ig                     ! grid point index
27      integer ig_vl1                 ! viking 1 grid point
28      real    pplay(ngrid,nlayer)    ! pressure at the middle of the layers (pa)
29      real    pplev(ngrid,nlayer+1)  ! pressure at layer boundaries (pa)
30      real    zzlay(ngrid,nlayer)    ! altitude at the middle of the layers (m)
31      real    zzlev(ngrid,nlayer+1)  ! altitude at layer boundaries (m)
32      real    zu(ngrid,nlayer)       ! u component of the wind (m.s-1)
33      real    zv(ngrid,nlayer)       ! v component of the wind (m.s-1)
34      real    zt(ngrid,nlayer)       ! temperature (k)
35      real    zycol(nlayer,nq)       ! composition (volume mixing ratio)
36      real    ptimestep              ! physical timestep (s)
37      real    co2ice(ngrid)          ! co2 ice surface layer (kg.m-2)
[415]38c
39c     local
40c
41      real    ubar                       ! wind module (m.s-1)
42      real    ustar                      ! friction velocity (m.s-1)
43      real    karman                     ! von karman constant
44      real    ra                         ! aerodynamic resistance (s.m-1)
45      real    rb                         ! quasi-laminar layer resistance (s.m-1)
46      real    vd                         ! dry deposition velocity (cm.s-1)
47      real    deltaz                     ! thickness of first layer (m)
48      real    mu0                        ! dynamic viscosity of co2 at 293.15 k (pa.s)
49      real    mu                         ! dynamic viscosity of co2 (pa.s)
50      real    nu                         ! kinematic viscosity of co2 (cm2.s-1)
51      real    rho                        ! density
52      real    d                          ! molecular diffusivity of methane (cm2.s-1)
53      real    schmidt                    ! schmidt number
54      real    prandtl                    ! prandtl number
55      real    p0                         ! reference pressure (pa)
56      real    loss                       ! loss rate (s-1)                       
57      integer iq
58c
59      real    nuch4
60      real    tau
61      real    gravity
62      real    gam
63      real    dp
64      real    cd
65c
66      data    karman / 0.4 /
67      data    prandtl / 0.72 /
68      data    mu0 / 14.8e-6 /
69      data    p0 /1.e5/
70c
71c     deposition is only active on surface uncovered by ice
72c
73c     if ((.not. watercaptag(ig)) .and. (co2ice(ig) .eq. 0.)) then
74c
75c        wind module (m.s-1)
76c
77         ubar = (zu(ig,1)**2. + zv(ig,1)**2.)**0.5
78c
79c        friction velocity (m.s-1)
80c
[421]81         ustar = ubar*karman/log(zzlay(ig,1)/z0(ig))
[415]82c
83c        aerodynamic resistance (s.m-1)
84c
[421]85         ra = 1./(karman*ustar)*log(zzlay(ig,1)/z0(ig))
[415]86c
87c        molecular diffusivity of methane in *air*  (cm2.s-1)
88c        massman, atmospheric environment, 32, 1111-1127, 1998
89c
90         d = 0.1952*(p0/pplay(ig,1))*(zt(ig,1)/273.15)**1.81
91c
92c        dynamic viscosity: temperature dependance (pa.s)
93c        sutherland's formula
94c
95         mu = mu0*(293.15 + 240.)/(zt(ig,1) + 240.)
96     $           *(zt(ig,1)/293.15)**(3./2.)
97c
98c        density (kg.m-3)
99c
100         rho = pplay(ig,1)/(rnew(ig,1)*zt(ig,1))
101c
102c        kinematic viscosity (cm2.s-1)
103c
104         nu = mu/rho*1.e4
105c
106c        schmidt number
107c
108         schmidt = nu/d
109c
110c        quasi-laminar layer resistance (s.m-1)
111c
112         rb = (2./(karman*ustar))*(schmidt/prandtl)**2./3.
113c
114c        dry deposition velocity (m.s-1)
115c
116         vd = 1./(ra + rb)
117c
118c        thickness of first layer (m)
119c
120         deltaz = zzlev(ig,2) - zzlev(ig,1)
121c
122c        loss rate (s-1)
123c
124         loss = vd/deltaz
125c
126c        test
127c
128c        loss = 1./(3600.*6.)
129c        vd = loss*deltaz
130c
131         nuch4 = sqrt(8.*8.31*zt(ig,1)
132     $                /(3.1416*16.e-3))
133         tau = 6.*3600.
134         gravity = 3.7
135         dp = pplev(ig,1) - pplev(ig,2)
[421]136         cd = (karman/log(zzlay(ig,1)/z0(ig)))**2.
[415]137c
138         gam = (4./nuch4)*dp
139     $        /(tau*gravity*rho - 1./(cd*ubar))
140c
141c        methane index
142c
143         iq = 12
144c
145c        update methane in first layer
146c
147c        zycol(1,iq) = zycol(1,iq)*exp(-loss*ptimestep)
148c
149         if (ig .eq. ig_vl1) then
150            print*,'**** deposition ****'
[421]151            print*,'z0     = ', z0(ig), 'm'
[415]152            print*,'deltaz = ', deltaz, ' m'
153            print*,'deltap = ', dp, 'pa'
154            print*,'zzlay  = ', zzlay(ig, 1), ' m'
155            print*,'pplay  = ', pplay(ig, 1), ' pa'
156            print*,'t      = ', zt(ig,1), ' k'
157            print*,'u      = ', zu(ig,1), ' m.s-1'
158            print*,'v      = ', zv(ig,1), ' m.s-1'
159            print*,'rho    = ', rho, ' kg.m-3'
160            print*,'ubar   = ', ubar, ' m.s-1'
161            print*,'d      = ', d, ' cm2.s-1'
162            print*,'mu     = ', mu, ' pa.s'
163            print*,'nu     = ', nu, ' cm2.s-1'
164            print*,'schmi  = ', schmidt
165            print*,'Ra     = ', ra, ' s.m-1'
166            print*,'Rb     = ', rb, ' s.m-1'
167            print*,'vd     = ', vd*100., 'cm.s-1'
168            print*,'tau dep= ', 1./loss, 's'
169            print*,'R      = ', rnew(ig,1)
170            print*,'nuch4  = ', nuch4, 'm.s-1'
171            print*,'tau    = ', tau, ' s'
172            print*,'gamma  = ', gam
173            print*,'taugrho= ',tau*gravity*rho
174            print*,'1surcdu= ', 1./(cd*ubar)
175
176            print*,'********************'
177         end if
178c
179c     end if
180c
[3466]181      end subroutine deposition
182
183      end module deposition_mod
Note: See TracBrowser for help on using the repository browser.