source: trunk/LMDZ.VENUS/libf/phyvenus/phytrac_chimie.F @ 2200

Last change on this file since 2200 was 2200, checked in by flefevre, 5 years ago

correction de bugs dans le supercycling de la chimie

  • Property svn:executable set to *
File size: 10.7 KB
Line 
1      subroutine phytrac_chimie (
2     $                    debutphy,
3     $                    gmtime,
4     $                    nqmax,
5     $                    nlon,
6     $                    lat,
7     $                    lon,
8     $                    nlev,
9     $                    pdtphys,
10     $                    temp,
11     $                    pplay,
12     $                    trac,
13     $                    d_tr_chem,
14     $                    iter)
15
16      use chemparam_mod
17      use conc, only: mmean
18
19      implicit none
20     
21#include "clesphys.h"
22#include "YOMCST.h"
23
24!===================================================================
25!     input
26!===================================================================
27
28      integer :: nlon, nlev                     ! number of gridpoints and levels
29      integer :: nqmax                          ! number of tracers
30
31      real :: gmtime                            ! day fraction
32      real :: pdtphys                           ! phytrac_chimie timestep (s)
33      real, dimension(nlon,nlev) :: temp        ! temperature (k)
34      real, dimension(nlon,nlev) :: pplay       ! pressure (pa)
35      real, dimension(nlon,nlev,nqmax) :: trac  ! tracer mass mixing ratio
36
37      logical :: debutphy                       ! first call flag
38
39!===================================================================
40!     output
41!===================================================================
42
43      real, dimension(nlon,nlev,nqmax) :: d_tr_chem  ! chemical tendency for each tracer
44      integer, dimension(nlon,nlev) :: iter          ! chemical iterations
45
46!===================================================================
47!     local
48!===================================================================
49
50      real :: sza_local   ! solar zenith angle (deg)
51      real :: lon_sun
52
53      integer :: i, iq
54      integer :: ilon, ilev
55
56      real  lat(nlon), lat_local(nlon)
57      real  lon(nlon), lon_local(nlon)
58
59      real, dimension(nlon,nlev) :: mrtwv, mrtsa ! total water and total sulfuric acid
60      real, dimension(nlon,nlev) :: mrwv, mrsa   ! gas-phase water and gas-phase sulfuric acid
61      real, dimension(nlon,nlev) :: trac_sum
62      real, dimension(nlon,nlev,nqmax) :: ztrac  ! local tracer mixing ratio
63     
64!===================================================================
65!     first call : initialisations
66!===================================================================
67
68      if (debutphy) then
69     
70!-------------------------------------------------------------------
71!        case of tracers re-initialisation with chemistry
72!-------------------------------------------------------------------
73         if (reinit_trac .and. ok_chem) then
74
75            print*, "Tracers are re-initialised"
76            trac(:,:,:) = 1.0e-30
77
78            if ((i_ocs /= 0) .and. (i_co /= 0) .and. (i_hcl /= 0)
79     $           .and. (i_so2 /= 0) .and. (i_h2o /= 0) .and. (i_n2/ = 0)
80     $           .and. (i_co2 /= 0)) then
81
82               trac(:,1:22,i_ocs) = 3.e-6
83               trac(:,1:22,i_co)  = 25.e-6
84               trac(:,:,i_hcl)    = 0.4e-6
85               trac(:,1:22,i_so2) = 10.e-6
86               trac(:,1:22,i_h2o) = 30.e-6
87               trac(:,:,i_n2)     = 0.35e-1
88   
89!          ensure that sum of mixing ratios = 1
90
91               trac_sum(:,:) = 0.
92
93               do iq = 1,nqmax - nmicro
94                  if (iq /= i_co2) then
95                     trac_sum(:,:) = trac_sum(:,:) + trac(:,:,iq)
96                  end if
97               end do
98
99!          initialise co2
100
101               trac(:,:,i_co2) = 1. - trac_sum(:,:)
102
103            else
104               write(*,*) "at least one tracer is missing : stop"
105               stop
106            end if
107       
108!           convert volume to mass mixing ratio
109
110            do iq = 1,nqmax - nmicro
111               trac(:,:,iq) = trac(:,:,iq)*m_tr(iq)/mmean(:,:)
112            end do
113   
114         end if  ! reinit_trac
115
116!-------------------------------------------------------------------
117!        case of detailed microphysics without chemistry
118!-------------------------------------------------------------------
119         if (.not. ok_chem .and. ok_cloud .and. cl_scheme == 2) then
120
121!           convert mass to volume mixing ratio
122
123            do iq = 1,nqmax - nmicro
124               ztrac(:,:,iq) = trac(:,:,iq)*mmean(:,:)/m_tr(iq)
125            end do
126         
127!           initialise microphysics
128 
129            call vapors4muphy_ini(nlon,nlev,ztrac)
130
131!           convert volume to mass mixing ratio
132
133            do iq = 1,nqmax - nmicro
134               trac(:,:,iq) = ztrac(:,:,iq)*m_tr(iq)/mmean(:,:)
135            end do
136   
137         end if
138
139      end if  ! debutphy
140
141!===================================================================
142!     convert mass to volume mixing ratio : gas phase
143!===================================================================
144
145      do iq = 1,nqmax - nmicro
146         ztrac(:,:,iq) = max(trac(:,:,iq)*mmean(:,:)/m_tr(iq), 1.e-30)
147      end do
148
149!===================================================================
150!     microphysics: simplified scheme (phd aurelien stolzenbach)
151!===================================================================
152
153      if (ok_cloud .and. cl_scheme == 1) then
154
155!     convert mass to volume mixing ratio : liquid phase
156
157         ztrac(:,:,i_h2so4liq) = max(trac(:,:,i_h2so4liq)
158     $                             *mmean(:,:)/m_tr(i_h2so4liq), 1.e-30)
159         ztrac(:,:,i_h2oliq) = max(trac(:,:,i_h2oliq)
160     $                             *mmean(:,:)/m_tr(i_h2oliq), 1.e-30)
161             
162!     total water and sulfuric acid (gas + liquid)
163
164         mrtwv(:,:) = ztrac(:,:,i_h2o) + ztrac(:,:,i_h2oliq)
165         mrtsa(:,:) = ztrac(:,:,i_h2so4) + ztrac(:,:,i_h2so4liq)
166
167!     all water and sulfuric acid is put in the gas-phase
168
169         mrwv(:,:) = mrtwv(:,:)
170         mrsa(:,:) = mrtsa(:,:)
171
172!     call microphysics
173
174         call new_cloud_venus(nlev, nlon, temp, pplay,
175     $                        mrtwv, mrtsa, mrwv, mrsa)
176
177!     update water vapour and sulfuric acid
178
179         ztrac(:,:,i_h2o) = mrwv(:,:)
180         ztrac(:,:,i_h2oliq) = mrtwv(:,:) - ztrac(:,:,i_h2o)
181     
182         ztrac(:,:,i_h2so4) = mrsa(:,:)
183         ztrac(:,:,i_h2so4liq) = mrtsa(:,:) - ztrac(:,:,i_h2so4)
184
185      end if  ! simplified scheme
186
187!===================================================================
188!     microphysics: detailed scheme (phd sabrina guilbon)
189!     !!! to be confirmed whether mad_muphy expects mmr or vmr for h2o and h2so4
190!===================================================================
191
192      if (ok_cloud .and. cl_scheme == 2) then
193
194         do iq = nqmax-nmicro+1,nqmax
195            ztrac(:,:,iq) = trac(:,:,iq)
196         end do
197
198         do ilon = 1,nlon       
199            do ilev = 1, nlev
200               if (temp(ilon,ilev) < 500.) then
201                  call mad_muphy(pdtphys,                               ! timestep
202     $                           temp(ilon,ilev),pplay(ilon,ilev),      ! temperature and pressure
203     $                           ztrac(ilon,ilev,i_h2o),
204     $                           ztrac(ilon,ilev,i_h2so4),     
205     $                           ztrac(ilon,ilev,i_m0_aer),
206     $                           ztrac(ilon,ilev,i_m3_aer),   
207     $                           ztrac(ilon,ilev,i_m0_mode1drop),
208     $                           ztrac(ilon,ilev,i_m0_mode1ccn),
209     $                           ztrac(ilon,ilev,i_m3_mode1sa),
210     $                           ztrac(ilon,ilev,i_m3_mode1w),   
211     $                           ztrac(ilon,ilev,i_m3_mode1ccn),   
212     $                           ztrac(ilon,ilev,i_m0_mode2drop),
213     $                           ztrac(ilon,ilev,i_m0_mode2ccn),
214     $                           ztrac(ilon,ilev,i_m3_mode2sa),
215     $                           ztrac(ilon,ilev,i_m3_mode2w),
216     $                           ztrac(ilon,ilev,i_m3_mode2ccn))
217               else
218                  ztrac(ilon,ilev,i_m0_aer)       = 0.
219                  ztrac(ilon,ilev,i_m3_aer)       = 0.
220                  ztrac(ilon,ilev,i_m0_mode1drop) = 0.
221                  ztrac(ilon,ilev,i_m0_mode1ccn)  = 0.
222                  ztrac(ilon,ilev,i_m3_mode1sa)   = 0.
223                  ztrac(ilon,ilev,i_m3_mode1w)    = 0.
224                  ztrac(ilon,ilev,i_m3_mode1ccn)  = 0.
225                  ztrac(ilon,ilev,i_m0_mode2drop) = 0.
226                  ztrac(ilon,ilev,i_m0_mode2ccn)  = 0.
227                  ztrac(ilon,ilev,i_m3_mode2sa)   = 0.
228                  ztrac(ilon,ilev,i_m3_mode2w)    = 0.
229                  ztrac(ilon,ilev,i_m3_mode2ccn)  = 0.
230               end if
231            end do
232         end do
233
234      end if  ! detailed scheme
235           
236!===================================================================
237!     photochemistry
238!===================================================================
239
240      if (ok_chem) then
241
242         lon_sun = (0.5 - gmtime)*2.*rpi
243         lon_local(:) = lon(:)*rpi/180.
244         lat_local(:) = lat(:)*rpi/180.
245
246         do ilon = 1,nlon
247
248!           solar zenith angle
249
250            sza_local = acos(cos(lat_local(ilon))*cos(lon_local(ilon))
251     $                 *cos(lon_sun) + cos(lat_local(ilon))
252     $                 *sin(lon_local(ilon))*sin(lon_sun))*180./rpi
253     
254            call photochemistry_venus(nlev, nlon, pdtphys,
255     $                                pplay(ilon,:)/100.,
256     $                                temp(ilon,:),
257     $                                ztrac(ilon,:,:),
258     $                                mmean(ilon,:),
259     $                                sza_local, nqmax, iter(ilon,:))
260
261         end do
262
263      end if  ! ok_chem
264
265!===================================================================
266!     compute tendencies
267!===================================================================
268
269!     gas phase
270
271      do iq = 1,nqmax - nmicro
272         ztrac(:,:,iq) = max(ztrac(:,:,iq)*m_tr(iq)/mmean(:,:),
273     $                       1.e-30)
274         d_tr_chem(:,:,iq) = (ztrac(:,:,iq) - trac(:,:,iq))/pdtphys
275      end do
276
277!     liquid phase or moments
278
279      if (ok_cloud .and. cl_scheme == 1) then
280         ztrac(:,:,i_h2so4liq) = max(ztrac(:,:,i_h2so4liq)
281     $                               *m_tr(i_h2so4liq)/mmean(:,:),
282     $                               1.e-30)
283         ztrac(:,:,i_h2oliq)   = max(ztrac(:,:,i_h2oliq)
284     $                               *m_tr(i_h2oliq)/mmean(:,:),
285     $                               1.e-30)
286         d_tr_chem(:,:,i_h2so4liq) = (ztrac(:,:,i_h2so4liq)
287     $                              - trac(:,:,i_h2so4liq))/pdtphys
288         d_tr_chem(:,:,i_h2oliq) = (ztrac(:,:,i_h2oliq)
289     $                            - trac(:,:,i_h2oliq))/pdtphys
290      end if
291
292
293      if (ok_cloud .and. cl_scheme == 2) then
294         do iq = nqmax-nmicro+1,nqmax
295            d_tr_chem(:,:,iq) = (ztrac(:,:,iq) - trac(:,:,iq))/pdtphys
296         end do
297      end if
298
299      end subroutine phytrac_chimie
Note: See TracBrowser for help on using the repository browser.