source: LMDZ6/branches/contrails/libf/phylmd/Dust/sediment_mod.f90 @ 5445

Last change on this file since 5445 was 5337, checked in by Laurent Fairhead, 5 weeks ago

Getting rid of dependance to dynamics
File size: 9.5 KB
Line 
1!----- This subroutine calculates the sedimentation flux of Tracers
2!
3SUBROUTINE sediment_mod(t_seri,pplay,zrho,paprs,time_step,RHcl, &
4        id_coss,id_codu,id_scdu, &
5        ok_chimeredust, &
6        sed_ss,sed_dust,sed_dustsco, &
7        sed_ss3D,sed_dust3D,sed_dustsco3D,tr_seri)
8  !nhl     .                                       xlon,xlat,
9  !
10USE chem_mod_h
11  USE yoecumf_mod_h
12     USE dimphy
13  USE infotrac_phy, ONLY: nbtr
14USE yomcst_mod_h
15IMPLICIT NONE
16  !
17
18    ! INCLUDE "dimphy.h"
19
20  !
21   REAL :: RHcl(klon,klev)     ! humidite relative ciel clair
22   REAL :: tr_seri(klon, klev,nbtr) !conc of tracers
23   REAL :: sed_ss(klon) !sedimentation flux of Sea Salt (g/m2/s)
24   REAL :: sed_dust(klon) !sedimentation flux of dust (g/m2/s)
25   REAL :: sed_dustsco(klon) !sedimentation flux of scoarse  dust (g/m2/s)
26   REAL :: sed_ss3D(klon,klev) !sedimentation flux of Sea Salt (g/m2/s)
27   REAL :: sed_dust3D(klon,klev) !sedimentation flux of dust (g/m2/s)
28   REAL :: sed_dustsco3D(klon,klev) !sedimentation flux of scoarse  dust (g/m2/s)
29   REAL :: t_seri(klon, klev)   !Temperature at mid points of Z (K)
30   REAL :: v_dep_ss(klon,klev)  ! sed. velocity for SS m/s
31   REAL :: v_dep_dust(klon,klev)  ! sed. velocity for dust m/s
32   REAL :: v_dep_dustsco(klon,klev)  ! sed. velocity for dust m/s
33   REAL :: pplay(klon, klev)    !pressure at mid points of Z (Pa)
34   REAL :: zrho(klon, klev)     !Density of air at mid points of Z (kg/m3)
35   REAL :: paprs(klon, klev+1)    !pressure at interface of layers Z (Pa)
36   REAL :: time_step            !time step (sec)
37   LOGICAL :: ok_chimeredust
38   REAL :: xlat(klon)       ! latitudes pour chaque point
39   REAL :: xlon(klon)       ! longitudes pour chaque point
40   INTEGER :: id_coss,id_codu,id_scdu
41  !
42  !------local variables
43  !
44   INTEGER :: i, k, nbre_RH
45   PARAMETER(nbre_RH=12)
46  !
47   REAL :: lambda, ss_g
48   REAL :: mmd_ss      !mass median diameter of SS (um)
49   REAL :: mmd_dust          !mass median diameter of dust (um)
50   REAL :: mmd_dustsco          !mass median diameter of scoarse dust (um)
51   REAL :: rho_ss(nbre_RH),rho_ss1 !density of sea salt (kg/m3)
52   REAL :: rho_dust          !density of dust(kg/m3)
53   REAL :: v_stokes, CC, v_sed, ss_growth_f(nbre_RH)
54   REAL :: sed_flux(klon,klev)  ! sedimentation flux g/m2/s
55   REAL :: air_visco(klon,klev)
56   REAL :: zdz(klon,klev)       ! layers height (m)
57   REAL :: temp                 ! temperature in degree Celius
58  !
59   INTEGER :: RH_num
60   REAL :: RH_MAX, DELTA, rh, RH_tab(nbre_RH)
61   PARAMETER (RH_MAX=95.)
62  !
63   DATA RH_tab/0.,10.,20.,30.,40.,50.,60.,70.,80.,85.,90.,95./
64  !
65  !
66   DATA rho_ss/2160. ,2160. ,2160.,  2160,  1451.6, 1367.9, &
67         1302.9,1243.2,1182.7, 1149.5,1111.6, 1063.1/
68  !
69   DATA ss_growth_f/0.503, 0.503, 0.503, 0.503, 0.724, 0.782, &
70         0.838, 0.905, 1.000, 1.072, 1.188, 1.447/
71  !
72  !
73   mmd_ss=12.7   !dia -um at 80% for bin 0.5-20 um but 90% of real mmd
74            ! obsolete      mmd_dust=2.8  !micrometer for bin 0.5-20 and 0.5-10 um
75  ! 4tracer SPLA:       mmd_dust=11.0  !micrometer for bin 0.5-20 and 0.5-10 um
76  !3days       mmd_dust=3.333464  !micrometer for bin 0.5-20 and 0.5-10 um
77  !3days       mmd_dustsco=12.91315  !micrometer for bin 0.5-20 and 0.5-10 um
78  !JE20140911       mmd_dust=3.002283  !micrometer for bin 0.5-20 and 0.5-10 um
79  !JE20140911       mmd_dustsco=13.09771  !micrometer for bin 0.5-20 and 0.5-10 um
80  !JE20140911        mmd_dust=5.156346  !micrometer for bin 0.5-20 and 0.5-10 um
81  !JE20140911        mmd_dustsco=15.56554  !micrometer for bin 0.5-20 and 0.5-10 um
82    IF (ok_chimeredust) THEN
83  !JE20150212<< : changes in ustar in dustmod changes emission distribution
84     ! mmd_dust=3.761212  !micrometer for bin 0.5-3 and 0.5-10 um
85     ! mmd_dustsco=15.06167  !micrometer for bin 3-20 and 0.5-10 um
86  !JE20150212>>
87  !JE20150618: Change in div3 of dustmod changes distribution. now is div3=6
88  !div=3        mmd_dust=3.983763
89  !div=3        mmd_dustsco=15.10854
90    mmd_dust=3.898047
91    mmd_dustsco=15.06167
92    ELSE
93    mmd_dust=11.0  !micrometer for bin 0.5-20 and 0.5-10 um
94    mmd_dustsco=100. ! absurd value, bin not used in this scheme
95    ENDIF
96
97
98   rho_dust=2600. !kg/m3
99  !
100  !--------- Air viscosity (poise=0.1 kg/m-sec)-----------
101  !
102   DO k=1, klev
103   DO i=1, klon
104  !
105   zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
106  !
107   temp=t_seri(i,k)-RTT
108  !
109   IF (temp.LT.0.) THEN
110     air_visco(i,k)=(1.718+0.0049*temp-1.2e-5*temp*temp)*1.e-4
111   ELSE
112     air_visco(i,k)=(1.718+0.0049*temp)*1.e-4
113   ENDIF
114  !
115   ENDDO
116   ENDDO
117  !
118  !--------- for Sea Salt -------------------
119  !
120  !
121  !
122   IF(id_coss>0) THEN
123   DO k=1, klev
124   DO i=1,klon
125  !
126  !---cal. correction factor hygroscopic growth of aerosols
127  !
128    rh=MIN(RHcl(i,k)*100.,RH_MAX)
129    RH_num = INT( rh/10. + 1.)
130    IF (rh.gt.85.) RH_num=10
131    IF (rh.gt.90.) RH_num=11
132    DELTA=(rh-RH_tab(RH_num))/(RH_tab(RH_num+1)-RH_tab(RH_num))
133  !
134    ss_g=ss_growth_f(rh_num) + &
135          DELTA*(ss_growth_f(RH_num+1)-ss_growth_f(RH_num))
136
137    rho_ss1=rho_ss(rh_num) + &
138          DELTA*(rho_ss(RH_num+1)-rho_ss(RH_num))
139  !
140    v_stokes=RG*(rho_ss1-zrho(i,k))* & !m/sec
141          (mmd_ss*ss_g)*(mmd_ss*ss_g)* &
142          1.e-12/(18.0*air_visco(i,k)/10.)
143  !
144   lambda=6.6*1.e-8*(103125/pplay(i,k))*(t_seri(i,k)/293.15)
145  !
146   CC=1.0+1.257*lambda/(mmd_ss*ss_g)/1.e6  ! C-correction factor
147  !
148   v_sed=v_stokes*CC                       ! m/sec !orig
149  !
150  !---------check for v_sed*dt<zdz
151  !
152   IF (v_sed*time_step.GT.zdz(i,k)) THEN
153     v_sed=zdz(i,k)/time_step
154   ENDIF
155  !
156   v_dep_ss(i,k)= v_sed
157   sed_flux(i,k)= tr_seri(i,k,id_coss)*v_sed !g/cm3*m/sec
158   ! !sed_ss3D(i,k)= -sed_flux(i,k)/zdz(i,k)      !g/cm3*sec !!!!!!!
159  ! ! conc_sed_ss3D(i,k)=sed_flux(i,k)*1.e6      !g/m3*sec !!!!!!!
160  !
161   ENDDO          !klon
162   ENDDO          !klev
163  !
164  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
165   sed_ss3D(:,:)=0.0  ! initialisation
166
167   DO k=1, klev
168   DO i=1, klon
169   sed_ss3D(i,k)=sed_ss3D(i,k)- &
170         sed_flux(i,k)/zdz(i,k) !!!!!!!!!!!!!!!!!!!!!!
171   ENDDO          !klon
172   ENDDO          !klev
173  !
174   DO k=1, klev-1
175   DO i=1, klon
176    sed_ss3D(i,k)=sed_ss3D(i,k)+ &
177          sed_flux(i,k+1)/zdz(i,k) !!!!!!!!
178
179   ENDDO          !klon
180   ENDDO          !klev
181
182  DO k = 1, klev
183  DO i = 1, klon
184      tr_seri(i,k,id_coss)=tr_seri(i,k,id_coss)+ &
185            sed_ss3D(i,k)*time_step
186  ENDDO
187  ENDDO
188
189  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
190  !
191   DO i=1, klon
192     sed_ss(i)=sed_flux(i,1)*1.e6*1.e3    !--unit mg/m2/s
193   ENDDO          !klon
194   ELSE
195    DO i=1, klon
196      sed_ss(i)=0.
197    ENDDO
198   ENDIF
199  !
200  !
201
202  !--------- For dust ------------------
203  !
204  !
205   IF(id_codu>0) THEN
206   DO k=1, klev
207   DO i=1,klon
208  !
209    v_stokes=RG*(rho_dust-zrho(i,k))* & !m/sec
210          mmd_dust*mmd_dust* &
211          1.e-12/(18.0*air_visco(i,k)/10.)
212  !
213   lambda=6.6*1.e-8*(103125/pplay(i,k))*(t_seri(i,k)/293.15)
214   CC=1.0+1.257*lambda/(mmd_dust)/1.e6        !dimensionless
215   v_sed=v_stokes*CC                       !m/sec
216  !
217  !---------check for v_sed*dt<zdz
218  !
219   IF (v_sed*time_step.GT.zdz(i,k)) THEN
220     v_sed=zdz(i,k)/time_step
221   ENDIF
222
223  !
224   v_dep_dust(i,k)= v_sed
225   sed_flux(i,k)  = tr_seri(i,k,id_codu)*v_sed !g/cm3.m/sec
226   ! !sed_dust3D(i,k)= -sed_flux(i,k)/zdz(i,k)      !g/cm3*sec !!!!!!!
227  !
228   ENDDO          !klon
229   ENDDO          !klev
230
231  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
232   sed_dust3D(:,:)=0.0  ! initialisation
233
234   DO k=1, klev
235   DO i=1, klon
236   sed_dust3D(i,k)=sed_dust3D(i,k)- &
237         sed_flux(i,k)/zdz(i,k)
238   ENDDO          !klon
239   ENDDO          !klev
240
241  !
242  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
243
244   DO k=1, klev-1
245   DO i=1, klon
246    sed_dust3D(i,k)=sed_dust3D(i,k) + &
247          sed_flux(i,k+1)/zdz(i,k)
248   ENDDO          !klon
249   ENDDO          !klev
250  !
251  DO k = 1, klev
252  DO i = 1, klon
253     tr_seri(i,k,id_codu)=tr_seri(i,k,id_codu)+ &
254           sed_dust3D(i,k)*time_step
255  ENDDO
256  ENDDO
257
258
259   DO i=1, klon
260     sed_dust(i)=sed_flux(i,1)*1.e6*1.e3    !--unit mg/m2/s
261   ENDDO          !klon
262   ELSE
263    DO i=1, klon
264      sed_dust(i)=0.
265    ENDDO
266   ENDIF
267  !
268
269
270  !--------- For scoarse  dust ------------------
271  !
272  !
273   IF(id_scdu>0) THEN
274   DO k=1, klev
275   DO i=1,klon
276  !
277    v_stokes=RG*(rho_dust-zrho(i,k))* & !m/sec
278          mmd_dustsco*mmd_dustsco* &
279          1.e-12/(18.0*air_visco(i,k)/10.)
280  !
281   lambda=6.6*1.e-8*(103125/pplay(i,k))*(t_seri(i,k)/293.15)
282   CC=1.0+1.257*lambda/(mmd_dustsco)/1.e6        !dimensionless
283   v_sed=v_stokes*CC                       !m/sec
284  !
285  !---------check for v_sed*dt<zdz
286
287
288   IF (v_sed*time_step.GT.zdz(i,k)) THEN
289     v_sed=zdz(i,k)/time_step
290   ENDIF
291
292  !
293   v_dep_dustsco(i,k)= v_sed
294   sed_flux(i,k)     = tr_seri(i,k,id_scdu)*v_sed !g/cm3.m/sec
295   ! !sed_dustsco3D(i,k)= -sed_flux(i,k)/zdz(i,k)      !g/cm3*sec !!!!!!!
296  !
297   ENDDO          !klon
298   ENDDO          !klev
299
300  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
301   sed_dustsco3D(:,:)=0.0  ! initialisation
302
303   DO k=1, klev
304   DO i=1, klon
305   sed_dustsco3D(i,k)=sed_dustsco3D(i,k)- &
306         sed_flux(i,k)/zdz(i,k)
307   ENDDO          !klon
308   ENDDO          !klev
309  !
310   DO k=1, klev-1
311   DO i=1, klon
312    sed_dustsco3D(i,k)=sed_dustsco3D(i,k) + &
313          sed_flux(i,k+1)/zdz(i,k)
314   ENDDO          !klon
315   ENDDO          !klev
316
317  DO k = 1, klev
318  DO i = 1, klon
319   tr_seri(i,k,id_scdu)=tr_seri(i,k,id_scdu)+ &
320         sed_dustsco3D(i,k)*time_step
321  ENDDO
322  ENDDO
323  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
324
325
326  !
327   DO i=1, klon
328     sed_dustsco(i)=sed_flux(i,1)*1.e6*1.e3    !--unit mg/m2/s
329   ENDDO          !klon
330   ELSE
331    DO i=1, klon
332      sed_dustsco(i)=0.
333    ENDDO
334   ENDIF
335  !
336
337
338
339
340  !
341   RETURN
342END SUBROUTINE sediment_mod
Note: See TracBrowser for help on using the repository browser.