Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

traccoag_mod.f90 @ 5310

Last change on this file since 5310 was 5268, checked in by abarral, 2 weeks ago
.f90 <-> .F90 depending on cpp key use
Property svn:keywords set to `Id`
File size: 19.9 KB

Line
1	MODULE traccoag_mod
2	!
3	! This module calculates the concentration of aerosol particles in certain size bins
4	! considering coagulation and sedimentation.
5	!
6	CONTAINS
7
8	SUBROUTINE traccoag(pdtphys, gmtime, debutphy, julien, &
9	presnivs, xlat, xlon, pphis, pphi, &
10	t_seri, pplay, paprs, sh, rh, tr_seri)
11
12	USE phys_local_var_mod, ONLY: mdw, R2SO4, DENSO4, f_r_wet, surf_PM25_sulf, &
13	& budg_emi_ocs, budg_emi_so2, budg_emi_h2so4, budg_emi_part, &
14	& R2SO4B, DENSO4B, f_r_wetB, sulfmmr, SAD_sulfate, sulfmmr_mode, nd_mode, reff_sulfate
15
16	USE dimphy
17	USE infotrac_phy, ONLY : nbtr_bin, nbtr_sulgas, nbtr, id_SO2_strat
18	USE aerophys
19	USE geometry_mod, ONLY : cell_area, boundslat
20	USE mod_grid_phy_lmdz
21	USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root
22	USE mod_phys_lmdz_para, only: gather, scatter
23	USE phys_cal_mod, ONLY : year_len, year_cur, mth_cur, day_cur, hour
24	USE sulfate_aer_mod
25	USE phys_local_var_mod, ONLY: stratomask
26	USE yomcst_mod_h
27	USE print_control_mod, ONLY: lunout
28	USE strataer_local_var_mod
29
30	IMPLICIT NONE
31
32	! Input argument
33	!---------------
34	REAL,INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
35	REAL,INTENT(IN) :: gmtime ! Heure courante
36	LOGICAL,INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
37	INTEGER,INTENT(IN) :: julien ! Jour julien
38
39	REAL,DIMENSION(klev),INTENT(IN) :: presnivs! pressions approximat. des milieux couches (en PA)
40	REAL,DIMENSION(klon),INTENT(IN) :: xlat ! latitudes pour chaque point
41	REAL,DIMENSION(klon),INTENT(IN) :: xlon ! longitudes pour chaque point
42	REAL,DIMENSION(klon),INTENT(IN) :: pphis ! geopotentiel du sol
43	REAL,DIMENSION(klon,klev),INTENT(IN) :: pphi ! geopotentiel de chaque couche
44
45	REAL,DIMENSION(klon,klev),INTENT(IN) :: t_seri ! Temperature
46	REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression pour le mileu de chaque couche (en Pa)
47	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression pour chaque inter-couche (en Pa)
48	REAL,DIMENSION(klon,klev),INTENT(IN) :: sh ! humidite specifique
49	REAL,DIMENSION(klon,klev),INTENT(IN) :: rh ! humidite relative
50
51	! Output argument
52	!----------------
53	REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
54
55	! Local variables
56	!----------------
57	REAL :: m_aer_emiss_vol_daily ! daily injection mass emission
58	REAL :: m_aer ! aerosol mass
59	INTEGER :: it, k, i, j, ilon, ilev, itime, i_int, ieru
60	LOGICAL,DIMENSION(klon,klev) :: is_strato ! true = above tropopause, false = below
61	REAL,DIMENSION(klon,klev) :: m_air_gridbox ! mass of air in every grid box [kg]
62	REAL,DIMENSION(klon_glo,klev,nbtr) :: tr_seri_glo ! Concentration Traceur [U/KgA]
63	REAL,DIMENSION(klev+1) :: altLMDz ! altitude of layer interfaces in m
64	REAL,DIMENSION(klev) :: f_lay_emiss ! fraction of emission for every vertical layer
65	REAL :: f_lay_sum ! sum of layer emission fractions
66	REAL :: alt ! altitude for integral calculation
67	INTEGER,PARAMETER :: n_int_alt=10 ! number of subintervals for integration over Gaussian emission profile
68	REAL,DIMENSION(nbtr_bin) :: r_bin ! particle radius in size bin [m]
69	REAL,DIMENSION(nbtr_bin) :: r_lower ! particle radius at lower bin boundary [m]
70	REAL,DIMENSION(nbtr_bin) :: r_upper ! particle radius at upper bin boundary [m]
71	REAL,DIMENSION(nbtr_bin) :: m_part_dry ! mass of one dry particle in size bin [kg]
72	REAL :: zrho ! Density of air [kg/m3]
73	REAL :: zdz ! thickness of atm. model layer in m
74	REAL,DIMENSION(klev) :: zdm ! mass of atm. model layer in kg
75	REAL,DIMENSION(klon,klev) :: dens_aer ! density of aerosol particles [kg/m3 aerosol] with default H2SO4 mass fraction
76	REAL :: emission ! emission
77	REAL :: theta_min, theta_max ! for SAI computation between two latitudes
78	REAL :: dlat_loc
79	REAL :: latmin,latmax,lonmin,lonmax ! lat/lon min/max for injection
80	REAL :: sigma_alt, altemiss ! injection altitude + sigma for distrib
81	REAL :: pdt,stretchlong ! physic timestep, stretch emission over one day
82
83	INTEGER :: injdur_sai ! injection duration for SAI case [days]
84	INTEGER :: yr, is_bissext
85	REAL :: samoment2, samoment3! 2nd and 3rd order moments of size distribution
86
87	IF (is_mpi_root .AND. flag_verbose_strataer) THEN
88	WRITE(lunout,*) 'in traccoag: date from phys_cal_mod =',year_cur,'-',mth_cur,'-',day_cur,'-',hour
89	WRITE(lunout,*) 'IN traccoag flag_emit: ',flag_emit
90	ENDIF
91
92	! radius [m]
93	DO it=1, nbtr_bin
94	r_bin(it)=mdw(it)/2.
95	ENDDO
96
97	!--set boundaries of size bins
98	DO it=1, nbtr_bin
99	IF (it.EQ.1) THEN
100	r_upper(it)=sqrt(r_bin(it+1)*r_bin(it))
101	r_lower(it)=r_bin(it)**2./r_upper(it)
102	ELSEIF (it.EQ.nbtr_bin) THEN
103	r_lower(it)=sqrt(r_bin(it)*r_bin(it-1))
104	r_upper(it)=r_bin(it)**2./r_lower(it)
105	ELSE
106	r_lower(it)=sqrt(r_bin(it)*r_bin(it-1))
107	r_upper(it)=sqrt(r_bin(it+1)*r_bin(it))
108	ENDIF
109	ENDDO
110
111	IF (debutphy .and. is_mpi_root) THEN
112	DO it=1, nbtr_bin
113	WRITE(lunout,*) 'radius bin', it, ':', r_bin(it), '(from', r_lower(it), 'to', r_upper(it), ')'
114	ENDDO
115	ENDIF
116
117	!--initialising logical is_strato from stratomask
118	is_strato(:,:)=.FALSE.
119	WHERE (stratomask.GT.0.5) is_strato=.TRUE.
120
121	IF(flag_new_strat_compo) THEN
122	IF(debutphy) WRITE(lunout,*) 'traccoag: COMPO/DENSITY (Tabazadeh 97) + H2O kelvin effect', flag_new_strat_compo
123	! STRACOMP (H2O, P, t_seri, R -> R2SO4 + Kelvin effect) : Taba97, Socol, etc...
124	CALL stracomp_kelvin(sh,t_seri,pplay)
125	ELSE
126	IF(debutphy) WRITE(lunout,*) 'traccoag: COMPO from Bekki 2D model', flag_new_strat_compo
127	! STRACOMP (H2O, P, t_seri -> aerosol composition (R2SO4))
128	! H2SO4 mass fraction in aerosol (%)
129	CALL stracomp(sh,t_seri,pplay)
130
131	! aerosol density (gr/cm3)
132	CALL denh2sa(t_seri)
133
134	! compute factor for converting dry to wet radius (for every grid box)
135	f_r_wet(:,:) = (dens_aer_dry/(DENSO4(:,:)1000.)/(R2SO4(:,:)/100.))*(1./3.)
136	ENDIF
137
138	!--calculate mass of air in every grid box
139	DO ilon=1, klon
140	DO ilev=1, klev
141	m_air_gridbox(ilon,ilev)=(paprs(ilon,ilev)-paprs(ilon,ilev+1))/RG*cell_area(ilon)
142	ENDDO
143	ENDDO
144
145	!--initialise emission diagnostics
146	if (nErupt > 0 .and. (flag_emit == 1 .or. flag_emit == 4)) budg_emi(:,1)=0.0
147	budg_emi_ocs(:)=0.0
148	budg_emi_so2(:)=0.0
149	budg_emi_h2so4(:)=0.0
150	budg_emi_part(:)=0.0
151
152	!--sulfur emission, depending on chosen scenario (flag_emit)
153	SELECT CASE(flag_emit)
154
155	CASE(0) ! background aerosol
156	! do nothing (no emission)
157
158	CASE(1) ! volcanic eruption
159	!--only emit on day of eruption
160	! stretch emission over one day of Pinatubo eruption
161	DO ieru=1, nErupt
162	IF (year_cur==year_emit_vol(ieru).AND.mth_cur==mth_emit_vol(ieru).AND.&
163	day_cur>=day_emit_vol(ieru).AND.day_cur<(day_emit_vol(ieru)+injdur)) THEN
164
165	! daily injection mass emission
166	m_aer=m_aer_emiss_vol(ieru,1)/(REAL(injdur)*REAL(ponde_lonlat_vol(ieru)))
167	!emission as SO2 gas (with m(SO2)=64/32*m_aer_emiss)
168	m_aer=m_aer*(mSO2mol/mSatom)
169
170	WRITE(lunout,*) 'IN traccoag m_aer_emiss_vol(ieru)=',m_aer_emiss_vol(ieru,1), &
171	'ponde_lonlat_vol(ieru)=',ponde_lonlat_vol(ieru),'(injdur*ponde_lonlat_vol(ieru))', &
172	(injdur*ponde_lonlat_vol(ieru)),'m_aer_emiss_vol_daily=',m_aer,'ieru=',ieru
173	WRITE(lunout,*) 'IN traccoag, dlon=',dlon
174
175	latmin=xlat_min_vol(ieru)
176	latmax=xlat_max_vol(ieru)
177	lonmin=xlon_min_vol(ieru)
178	lonmax=xlon_max_vol(ieru)
179	altemiss = altemiss_vol(ieru)
180	sigma_alt = sigma_alt_vol(ieru)
181	pdt=pdtphys
182	! stretch emission over one day of eruption
183	stretchlong = 1.
184
185	CALL STRATEMIT(pdtphys,pdt,xlat,xlon,t_seri,pplay,paprs,tr_seri,&
186	m_aer,latmin,latmax,lonmin,lonmax,altemiss,sigma_alt,id_SO2_strat, &
187	stretchlong,1,0)
188
189	ENDIF ! emission period
190	ENDDO ! eruption number
191
192	CASE(2) ! stratospheric aerosol injections (SAI)
193	!
194	! Computing duration of SAI in days...
195	! ... starting from 0...
196	injdur_sai = 0
197	! ... then adding whole years from first to (n-1)th...
198	DO yr = year_emit_sai_start, year_emit_sai_end-1
199	! (n % 4 == 0) and (n % 100 != 0 or n % 400 == 0)
200	is_bissext = (MOD(yr,4)==0) .AND. (MOD(yr,100) /= 0 .OR. MOD(yr,400) == 0)
201	injdur_sai = injdur_sai+365+is_bissext
202	ENDDO
203	! ... then subtracting part of the first year where no injection yet...
204	is_bissext = (MOD(year_emit_sai_start,4)==0) .AND. (MOD(year_emit_sai_start,100) /= 0 .OR. MOD(year_emit_sai_start,400) == 0)
205	SELECT CASE(mth_emit_sai_start)
206	CASE(2)
207	injdur_sai = injdur_sai-31
208	CASE(3)
209	injdur_sai = injdur_sai-31-28-is_bissext
210	CASE(4)
211	injdur_sai = injdur_sai-31-28-is_bissext-31
212	CASE(5)
213	injdur_sai = injdur_sai-31-28-is_bissext-31-30
214	CASE(6)
215	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31
216	CASE(7)
217	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30
218	CASE(8)
219	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30-31
220	CASE(9)
221	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30-31-31
222	CASE(10)
223	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30-31-31-30
224	CASE(11)
225	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30-31-31-30-31
226	CASE(12)
227	injdur_sai = injdur_sai-31-28-is_bissext-31-30-31-30-31-31-30-31-30
228	END SELECT
229	injdur_sai = injdur_sai-day_emit_sai_start+1
230	! ... then adding part of the n-th year
231	is_bissext = (MOD(year_emit_sai_end,4)==0) .AND. (MOD(year_emit_sai_end,100) /= 0 .OR. MOD(year_emit_sai_end,400) == 0)
232	SELECT CASE(mth_emit_sai_end)
233	CASE(2)
234	injdur_sai = injdur_sai+31
235	CASE(3)
236	injdur_sai = injdur_sai+31+28+is_bissext
237	CASE(4)
238	injdur_sai = injdur_sai+31+28+is_bissext+31
239	CASE(5)
240	injdur_sai = injdur_sai+31+28+is_bissext+31+30
241	CASE(6)
242	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31
243	CASE(7)
244	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30
245	CASE(8)
246	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30+31
247	CASE(9)
248	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30+31+31
249	CASE(10)
250	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30+31+31+30
251	CASE(11)
252	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30+31+31+30+31
253	CASE(12)
254	injdur_sai = injdur_sai+31+28+is_bissext+31+30+31+30+31+31+30+31+30
255	END SELECT
256	injdur_sai = injdur_sai+day_emit_sai_end
257	! A security: are SAI dates of injection consistent?
258	IF (injdur_sai <= 0) THEN
259	CALL abort_physic('traccoag_mod', 'Pb in SAI dates of injection.',1)
260	ENDIF
261	! Injection in itself
262	IF (( year_emit_sai_start <= year_cur ) &
263	.AND. ( year_cur <= year_emit_sai_end ) &
264	.AND. ( mth_emit_sai_start <= mth_cur .OR. year_emit_sai_start < year_cur ) &
265	.AND. ( mth_cur <= mth_emit_sai_end .OR. year_cur < year_emit_sai_end ) &
266	.AND. ( day_emit_sai_start <= day_cur .OR. mth_emit_sai_start < mth_cur .OR. year_emit_sai_start < year_cur ) &
267	.AND. ( day_cur <= day_emit_sai_end .OR. mth_cur < mth_emit_sai_end .OR. year_cur < year_emit_sai_end )) THEN
268
269	m_aer=m_aer_emiss_sai
270	!emission as SO2 gas (with m(SO2)=64/32*m_aer_emiss)
271	m_aer=m_aer*(mSO2mol/mSatom)
272
273	latmin=xlat_sai
274	latmax=xlat_sai
275	lonmin=xlon_sai
276	lonmax=xlon_sai
277	altemiss = altemiss_sai
278	sigma_alt = sigma_alt_sai
279	pdt=0.
280	! stretch emission over whole year (360d)
281	stretchlong=FLOAT(year_len)
282
283	CALL STRATEMIT(pdtphys,pdt,xlat,xlon,t_seri,pplay,paprs,m_air_gridbox,tr_seri,&
284	m_aer,latmin,latmax,lonmin,lonmax,altemiss,sigma_alt,id_SO2_strat, &
285	stretchlong,1,0)
286
287	budg_emi_so2(:) = budg_emi(:,1)*mSatom/mSO2mol
288	ENDIF ! Condition over injection dates
289
290	CASE(3) ! --- SAI injection over a single band of longitude and between
291	! lat_min and lat_max
292
293	m_aer=m_aer_emiss_sai
294	!emission as SO2 gas (with m(SO2)=64/32*m_aer_emiss)
295	m_aer=m_aer*(mSO2mol/mSatom)
296
297	latmin=xlat_min_sai
298	latmax=xlat_max_sai
299	lonmin=xlon_sai
300	lonmax=xlon_sai
301	altemiss = altemiss_sai
302	sigma_alt = sigma_alt_sai
303	pdt=0.
304	! stretch emission over whole year (360d)
305	stretchlong=FLOAT(year_len)
306
307	CALL STRATEMIT(pdtphys,pdt,xlat,xlon,t_seri,pplay,paprs,m_air_gridbox,tr_seri,&
308	m_aer,latmin,latmax,lonmin,lonmax,altemiss,sigma_alt,id_SO2_strat, &
309	stretchlong,1,0)
310
311	budg_emi_so2(:) = budg_emi(:,1)*mSatom/mSO2mol
312
313	END SELECT ! emission scenario (flag_emit)
314
315	!--read background concentrations of OCS and SO2 and lifetimes from input file
316	!--update the variables defined in phys_local_var_mod
317	CALL interp_sulf_input(debutphy,pdtphys,paprs,tr_seri)
318
319	!--convert OCS to SO2 in the stratosphere
320	CALL ocs_to_so2(pdtphys,tr_seri,t_seri,pplay,paprs,is_strato)
321
322	!--convert SO2 to H2SO4
323	CALL so2_to_h2so4(pdtphys,tr_seri,t_seri,pplay,paprs,is_strato)
324
325	!--common routine for nucleation and condensation/evaporation with adaptive timestep
326	CALL micphy_tstep(pdtphys,tr_seri,t_seri,pplay,paprs,rh,is_strato)
327
328	!--call coagulation routine
329	CALL coagulate(pdtphys,mdw,tr_seri,t_seri,pplay,dens_aer,is_strato)
330
331	!--call sedimentation routine
332	CALL aer_sedimnt(pdtphys, t_seri, pplay, paprs, tr_seri, dens_aer)
333
334	!--compute mass concentration of PM2.5 sulfate particles (wet diameter and mass) at the surface for health studies
335	surf_PM25_sulf(:)=0.0
336	DO i=1,klon
337	DO it=1, nbtr_bin
338	IF (mdw(it) .LT. 2.5e-6) THEN
339	!surf_PM25_sulf(i)=surf_PM25_sulf(i)+tr_seri(i,1,it+nbtr_sulgas)*m_part(i,1,it) &
340	!assume that particles consist of ammonium sulfate at the surface (132g/mol)
341	!and are dry at T = 20 deg. C and 50 perc. humidity
342	surf_PM25_sulf(i)=surf_PM25_sulf(i)+tr_seri(i,1,it+nbtr_sulgas) &
343	& 132./98.dens_aer_dry4./3.RPI(mdw(it)/2.)*3 &
344	& pplay(i,1)/t_seri(i,1)/RD1.e9
345	ENDIF
346	ENDDO
347	ENDDO
348
349	!--compute
350	! sulfmmr: Sulfate aerosol concentration (dry mixing ratio) (condensed H2SO4 mmr)
351	! SAD_sulfate: SAD all aerosols (cm2/cm3) (must be WET)
352	! sulfmmr_mode: sulfate(=H2SO4 if dry) MMR in different modes (ambiguous but based on sulfmmr, it mus be DRY(?) mmr)
353	! nd_mode: DRY(?) particle concentration in different modes (part/m3)
354	sulfmmr(:,:)=0.0
355	SAD_sulfate(:,:)=0.0
356	sulfmmr_mode(:,:,:)=0.0
357	nd_mode(:,:,:)=0.0
358	reff_sulfate(:,:)=0.0
359
360	DO i=1,klon
361	DO j=1,klev
362	samoment2=0.0
363	samoment3=0.0
364	DO it=1, nbtr_bin
365	!surf_PM25_sulf(i)=surf_PM25_sulf(i)+tr_seri(i,1,it+nbtr_sulgas)*m_part(i,1,it) &
366	!assume that particles consist of ammonium sulfate at the surface (132g/mol)
367	!and are dry at T = 20 deg. C and 50 perc. humidity
368
369	! sulfmmr_mode: sulfate(=H2SO4 if dry) MMR in different modes (based on sulfmmr, it must be DRY mmr)
370	! equivalent to condensed H2SO4 mmr= H2SO4 kg / kgA in bin it
371	sulfmmr_mode(i,j,it) = tr_seri(i,j,it+nbtr_sulgas) & ! [DRY part/kgA in bin it]
372	& (4./3.)RPI(mdw(it)/2.)*3. & ! [mdw: dry diameter in m]
373	& *dens_aer_dry ! [dry aerosol mass density in kg/m3]
374
375	! sulfmmr: Sulfate aerosol concentration (dry mass mixing ratio)
376	! equivalent to total condensed H2SO4 mmr (H2SO4 kg / kgA
377	sulfmmr(i,j) = sulfmmr(i,j) + sulfmmr_mode(i,j,it)
378
379	! nd_mode: particle concentration in different modes (DRY part/m3)
380	nd_mode(i,j,it) = tr_seri(i,j,it+nbtr_sulgas) & ! [DRY part/kgA in bin it]
381	& *pplay(i,j)/t_seri(i,j)/RD ! [air mass concentration in kg air /m3A]
382
383	IF(flag_new_strat_compo) THEN
384	! SAD_sulfate: SAD WET sulfate aerosols (cm2/cm3)
385	SAD_sulfate(i,j) = SAD_sulfate(i,j) + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
386	& 4.RPI( mdw(it)f_r_wetB(i,j,it)/2. )**2. & ! [WET SA of part it in m2]
387	& *1.e-2 ! conversion from m2/m3 to cm2/cm3A
388	! samoment2 : 2nd order moment of WET sulfate aerosols (m2/m3)
389	samoment2 = samoment2 + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
390	& ( mdw(it)f_r_wetB(i,j,it)/2. )**2. ! [WET SA of part it in m2]
391	! samoment3 : 3nd order moment of WET sulfate aerosols (cm2/cm3)
392	samoment3 = samoment3 + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
393	& ( mdw(it)f_r_wetB(i,j,it)/2. )**3. ! [WET SA of part it in m2]
394	ELSE
395	! SAD_sulfate: SAD WET sulfate aerosols (cm2/cm3)
396	SAD_sulfate(i,j) = SAD_sulfate(i,j) + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
397	& 4.RPI( mdw(it)f_r_wet(i,j)/2. )**2. & ! [WET SA of part it in m2]
398	& *1.e-2 ! conversion from m2/m3 to cm2/cm3A
399	! samoment2 : 2nd order moment of WET sulfate aerosols (m2/m3)
400	samoment2 = samoment2 + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
401	& ( mdw(it)f_r_wet(i,j)/2. )**2. ! [WET SA of part it in m2]
402	! samoment3 : 3nd order moment of WET sulfate aerosols (cm2/cm3)
403	samoment3 = samoment3 + nd_mode(i,j,it) & ! [DRY part/m3A (in bin it)]
404	& ( mdw(it)f_r_wet(i,j)/2. )**3. ! [WET SA of part it in m2]
405	ENDIF
406	ENDDO
407	! reff_sulfate: effective radius of WET sulfate aerosols (cm)
408	IF(samoment2 > 1.e-30) THEN
409	reff_sulfate(i,j) = (samoment3 / samoment2) &
410	& *1.e2 ! conversion from m to cm
411
412	! Sanity check
413	IF(reff_sulfate(i,j) > 5.e-4) reff_sulfate(i,j) = 5.e-4 ! reff_sulfate max = 5 micron
414	IF(reff_sulfate(i,j) < 1.e-6) reff_sulfate(i,j) = 1.e-6 ! reff_sulfate min = 10 nm
415	ELSE
416	reff_sulfate(i,j) = 1.e-5 ! reff_sulfate N part = nul (ref = 100 nm)
417	ENDIF
418	ENDDO
419	ENDDO
420
421	END SUBROUTINE traccoag
422
423	END MODULE traccoag_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format