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

Last change on this file since 1443 was 1442, checked in by slebonnois, 10 years ago

SL: update of the Venus GCM, + corrections on routines used for newstart/start2archive for Titan and Venus, + some modifications on tools

  • Property svn:executable set to *
File size: 9.0 KB
Line 
1!
2! $Header: /home/cvsroot/LMDZ4/libf/phylmd/phytrac.F,v 1.16 2006/03/24 15:06:23 lmdzadmin Exp $
3!
4c
5c
6      SUBROUTINE  phytrac_chimie (
7     I                    debutphy,
8     I                    gmtime,
9     I                    nqmax,
10     I                    n_lon,
11     I                    lat,
12     I                    lon,
13     I                    n_lev,
14     I                    pdtphys,
15     I                    temp,
16     I                    pplev,
17     O                    trac)
18c======================================================================
19c Auteur(s) FH
20c Objet: Moniteur general des tendances traceurs
21c
22cAA Remarques en vrac:
23cAA--------------------
24cAA 1/ le call phytrac se fait avec nqmax
25c======================================================================
26c      USE ioipsl
27c      USE infotrac
28c      USE control_mod
29c      USE dimphy
30c      USE comgeomphy
31      USE chemparam_mod
32      use conc, only: mmean
33      IMPLICIT none
34     
35c#include "dimensions.h"
36#include "clesphys.h"
37c#include "temps.h"
38c#include "paramet.h"
39c#include "comcstfi.h" !me permet de recuperer mugaz et d'autres constantes comme rad,pi etc
40#include "YOMCST.h"
41c======================================================================
42c Arguments:
43c
44c
45c   EN ENTREE:
46c   ==========
47c
48c   divers:
49c   -------
50c
51
52      REAL  sza_local
53      REAL  gmtime
54c      INTEGER, SAVE :: cpt_cloudIO  !un compteur pour fichier sortie cloud_parameter en 1D
55      INTEGER  iq
56      INTEGER  i
57      INTEGER  ilon, ilev
58      integer  n_lon  ! nombre de points horizontaux
59      INTEGER  n_lev  ! nombre de couches verticales
60      INTEGER  nqmax ! nombre de traceurs auxquels on applique la physique
61
62      real  pdtphys  ! pas d'integration pour la physique (seconde)
63      real  lat(n_lon), lat_local(n_lon)
64      real  lon(n_lon), lon_local(n_lon)
65      real  temp(n_lon,n_lev) ! temp
66      real  trac(n_lon,n_lev,nqmax) ! traceur
67      real  trac_sav(n_lon,n_lev,nqmax)
68      real  trac_sum(n_lon,n_lev)
69      real  pplev(n_lon,n_lev)  ! pression pour le mileu de chaque couche (en Pa)
70      real  lon_sun
71
72      logical debutphy       ! le flag de l'initialisation de la physique
73
74C     Auxilary variables:
75 
76      REAL, DIMENSION(n_lon,n_lev) :: mrtwv,mrtsa,
77     +                                mrwv,mrsa
78     
79C    ps_sa: satur pressure pure SA
80C    satps_sa: satur pres over mixture in dyne/cm2=Pa/10
81C----------------------------------------------------------------------------
82       
83C Variables liees a l'ecriture de la bande histoire : phytrac.nc
84 
85      logical ok_sync
86      parameter (ok_sync = .true.)
87
88c      modname = 'phytrac'
89
90c      PRINT*,'DEBUT subroutine PHYTRAC'
91
92c----------------------------------
93c debutphy: Initiation des traceurs
94c----------------------------------
95
96      if (debutphy) then
97     
98         if (n_lon .EQ. 1) then           
99         PRINT*,'n_lon 1D: ',n_lon
100         end if
101                 
102c         if ((n_lon .GT. 1) .AND. ok_chem) then
103c !!! DONC 3D !!!       
104c           CALL chemparam_ini()
105c         endif
106         
107c         if ((n_lon .GT. 1) .AND. ok_cloud) then
108c !!! DONC 3D !!!
109c         CALL cloud_ini(n_lon,n_lev)
110c         endif
111           
112       IF (reinit_trac) THEN
113       PRINT*,'REINIT MIXING RATIO TRACEURS'
114
115c       =============================================================
116c                                       Passage de Rm à Rv
117c       =============================================================
118c       Necessaire si on reprend les start.nc qui sont en MMR
119
120         DO iq=1,nqmax
121          trac(:,:,iq)=trac(:,:,iq)*mmean(:,:)/M_tr(iq)
122         END DO
123c       =============================================================
124         
125c=============================================================
126c               Initialisation des profils traceurs en Rv
127c=============================================================
128c initialisation sert a mettre les valeurs voulues par utilisateur pour
129c chaque traceur
130c exemple: trac(ilon,ilev,q)=xx
131
132c     trac_sav sert a sauver les valeurs initiales du start.nc     
133      trac_sav=trac
134
135c     On initialise les traceurs a zero obligatoire pour la chimie
136      trac(:,:,:)=1.0E-30
137
138c     !!! AVEC NUAGE !!!
139      trac(:,1:22,i_ocs)=3.E-6
140      trac(:,:,i_hcl)=0.4E-6
141      trac(:,1:22,i_so2)=10.E-6
142      trac(:,1:22,i_h2o)=30.0E-6
143
144c     remettre tous les traceurs du start => trac(:,:,:)=trac_sav(:,:,:)
145
146c     N2 n est pas encore une espece chimique du modele chimique
147c     traceur passif pour la chimie-transport
148      trac(:,:,i_n2)=0.35d-1
149   
150!!!! GG:   Initialization CO2 = 1 - qtot
151!!    It assures that vmr_tot = 1
152c     On a donc le CO2 qui est le restant d atmosphere Venus 
153         trac_sum(:,:)=0.0
154        DO iq=2,nqmax
155         trac_sum(:,:)= trac_sum(:,:) + trac(:,:,iq)
156        END DO
157
158        trac(:,:,i_co2)= 1-trac_sum(:,:)
159       
160c=============================================================
161     
162c       =============================================================
163c                                       Passage de Rv à Rm
164c       =============================================================
165         DO iq=1,nqmax
166          trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/mmean(:,:)
167         END DO
168c       =============================================================
169   
170       ENDIF  !FIN REINIT TRAC
171
172       
173c-------------
174c fin debutphy
175c-------------
176
177      ENDIF  ! fin debutphy
178
179c       =============================================================
180c                                       Passage de Rm à Rv
181c       =============================================================
182       DO iq=1,nqmax
183         trac(:,:,iq)=MAX(trac(:,:,iq)*mmean(:,:)/M_tr(iq),1.E-30)
184       END DO
185c       =============================================================
186
187
188c       =============================================================
189c                        Appel Microphysique (sans nucleation)
190c                  Volume Mixing Ratio
191c       =============================================================
192
193         IF (ok_cloud) THEN
194               
195c      PRINT*,'DEBUT CLOUD'     
196c     On remet tout le RM liq dans la partie gaz
197c     !!! On reforme un nuage à chaque fois !!!
198
199      DO ilev=1, n_lev
200      DO ilon=1, n_lon         
201      mrtwv(ilon,ilev)=trac(ilon,ilev,i_h2o) +
202     &  trac(ilon,ilev,i_h2oliq)
203      mrtsa(ilon,ilev)=trac(ilon,ilev,i_h2so4) +
204     &  trac(ilon,ilev,i_h2so4liq)
205      mrwv(ilon,ilev)=mrtwv(ilon,ilev)
206      mrsa(ilon,ilev)=mrtsa(ilon,ilev)
207      ENDDO
208      ENDDO
209                   
210      CALL new_cloud_venus(n_lev, n_lon,
211     e temp,pplev,
212     e mrtwv,mrtsa,
213     e mrwv,mrsa)
214
215c       =========================================               
216c       Actualisation des mixing ratio liq et gaz
217c       =========================================
218c       Si la routine new_cloud_venus n'a pas actualisé mrwv et mrsa
219c       on a alors bien mr=mrt pour sa et wv, donc les parties liq sont=0 hors du nuage
220c       ou si on ne condense pas
221
222c      PRINT*,'DEBUT ACTUALISATION OUTPUT CLOUD'
223c    si tout se passe bien, mrtwv et mrtsa ne changent pas
224      DO ilev=1, n_lev
225      DO ilon=1, n_lon       
226      trac(ilon,ilev,i_h2o) = mrwv(ilon,ilev)
227      trac(ilon,ilev,i_h2oliq) = mrtwv(ilon,ilev) -
228     &  trac(ilon,ilev,i_h2o)
229     
230      trac(ilon,ilev,i_h2so4) = mrsa(ilon,ilev)
231      trac(ilon,ilev,i_h2so4liq) = mrtsa(ilon,ilev) -
232     &  trac(ilon,ilev,i_h2so4)
233      ENDDO
234      ENDDO
235
236c       =============================================================
237c      PRINT*,'FIN CLOUD'
238      ENDIF
239           
240c=============================================================
241c               CHIMIE: Boucle sur les lon, lat (n_lon)
242c=============================================================
243
244c     AS:
245c     Ici, la longitude au midi local se deplace vers l'Ouest
246c     c'est le sens terrestre
247c     pour Vénus on prend juste l'opposé de la longitude et on a la rotation
248c     de Vénus et donc le midi local qui se déplace vers l'Est
249     
250      lon_sun = (0.5 - gmtime) * 2.0 * RPI
251      lon_local = lon * RPI/180.0E+0
252      lat_local = lat * RPI/180.0E+0
253       
254      DO ilon=1, n_lon
255
256c     calcul sza_local pour obtenir des sza_local > 90, utile pour la chimie
257      sza_local = acos(cos(lat_local(ilon))*cos(lon_local(ilon))*
258     & cos(lon_sun) + cos(lat_local(ilon))*sin(lon_local(ilon))
259     & *sin(lon_sun))* 180.0E+0/RPI
260     
261c      PRINT*,'sza_local :', sza_local     
262   
263      IF (ok_chem) THEN
264c      PRINT*,'DEBUT CHEMISTRY'
265c       =============================================================
266c                                       Appel Photochimie
267c       =============================================================
268c     Pression en hPa => pplev/100.
269       
270      CALL new_photochemistry_venus(n_lev, n_lon, pdtphys,
271     e                         pplev(ilon,:)/100.,
272     e                         temp(ilon,:),
273     e                         trac(ilon,:,:),
274     e                         sza_local, nqmax)
275c       =============================================================
276c      PRINT*,'FIN CHEMISTRY'
277   
278        END IF
279
280      END DO
281c       =============================================================
282c                                       Passage de Rv à Rm
283c       =============================================================
284        DO iq=1,nqmax
285c               trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/RMD
286                trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/mmean(:,:)
287
288        END DO
289c       =============================================================   
290C      PRINT*,'FIN PHYTRAC'
291      RETURN
292      END
Note: See TracBrowser for help on using the repository browser.