Context Navigation

phytrac_chimie.F @ 1667

Last change on this file since 1667 was 1665, checked in by slebonnois, 8 years ago
SL: oubli commit precedent + modif calcul N2 + XIOS
Property svn:executable set to ``*
File size: 12.0 KB

Rev	Line
[1305]	1	!
	2	! $Header: /home/cvsroot/LMDZ4/libf/phylmd/phytrac.F,v 1.16 2006/03/24 15:06:23 lmdzadmin Exp $
	3	!
	4	c
	5	c
	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,
[1661]	16	I pplay,
[1442]	17	O trac)
[1305]	18	c======================================================================
	19	c Auteur(s) FH
	20	c Objet: Moniteur general des tendances traceurs
	21	c
	22	cAA Remarques en vrac:
	23	cAA--------------------
	24	cAA 1/ le call phytrac se fait avec nqmax
	25	c======================================================================
	26	USE chemparam_mod
[1442]	27	use conc, only: mmean
[1305]	28	IMPLICIT none
	29
	30	#include "clesphys.h"
	31	#include "YOMCST.h"
	32	c======================================================================
	33	c Arguments:
	34	c
	35	c
	36	c EN ENTREE:
	37	c ==========
	38	c
	39	c divers:
	40	c -------
	41	c
	42
	43	REAL sza_local
	44	REAL gmtime
	45	c INTEGER, SAVE :: cpt_cloudIO !un compteur pour fichier sortie cloud_parameter en 1D
	46	INTEGER iq
	47	INTEGER i
	48	INTEGER ilon, ilev
	49	integer n_lon ! nombre de points horizontaux
	50	INTEGER n_lev ! nombre de couches verticales
	51	INTEGER nqmax ! nombre de traceurs auxquels on applique la physique
[1442]	52
[1305]	53	real pdtphys ! pas d'integration pour la physique (seconde)
	54	real lat(n_lon), lat_local(n_lon)
	55	real lon(n_lon), lon_local(n_lon)
	56	real temp(n_lon,n_lev) ! temp
[1442]	57	real trac(n_lon,n_lev,nqmax) ! traceur
	58	real trac_sav(n_lon,n_lev,nqmax)
	59	real trac_sum(n_lon,n_lev)
[1661]	60	real pplay(n_lon,n_lev) ! pression pour le mileu de chaque couche (en Pa)
[1305]	61	real lon_sun
[1442]	62
[1305]	63	logical debutphy ! le flag de l'initialisation de la physique
	64
	65	C Auxilary variables:
	66
[1442]	67	REAL, DIMENSION(n_lon,n_lev) :: mrtwv,mrtsa,
	68	+ mrwv,mrsa
	69
[1305]	70	C ps_sa: satur pressure pure SA
	71	C satps_sa: satur pres over mixture in dyne/cm2=Pa/10
	72	C----------------------------------------------------------------------------
	73
	74	C Variables liees a l'ecriture de la bande histoire : phytrac.nc
	75
	76	logical ok_sync
	77	parameter (ok_sync = .true.)
	78
	79	c modname = 'phytrac'
	80
	81	c PRINT*,'DEBUT subroutine PHYTRAC'
	82
	83	c----------------------------------
	84	c debutphy: Initiation des traceurs
	85	c----------------------------------
	86
	87	if (debutphy) then
	88
	89	if (n_lon .EQ. 1) then
	90	PRINT*,'n_lon 1D: ',n_lon
	91	end if
[1661]	92
	93	c=============================================================
	94	c == CASE INIT PHOTOCHEM ==
	95	c=============================================================
	96
[1665]	97	IF (reinit_trac.and.ok_chem) THEN
[1442]	98	PRINT*,'REINIT MIXING RATIO TRACEURS'
[1305]	99
[1661]	100	c Passage de Rm a Rv
[1305]	101	c =============================================================
[1442]	102	c Necessaire si on reprend les start.nc qui sont en MMR
[1661]	103	c SEULEMENT LES GAZ
	104	DO iq=1,nqmax-nmicro
[1442]	105	trac(:,:,iq)=trac(:,:,iq)*mmean(:,:)/M_tr(iq)
	106	END DO
[1305]	107	c =============================================================
	108
	109	c Initialisation des profils traceurs en Rv
	110	c=============================================================
[1442]	111	c initialisation sert a mettre les valeurs voulues par utilisateur pour
	112	c chaque traceur
	113	c exemple: trac(ilon,ilev,q)=xx
[1305]	114
[1442]	115	c trac_sav sert a sauver les valeurs initiales du start.nc
	116	trac_sav=trac
[1305]	117
[1442]	118	c On initialise les traceurs a zero obligatoire pour la chimie
	119	trac(:,:,:)=1.0E-30
[1305]	120
[1665]	121	c !!! Verif traceurs !!!
	122	if ( (i_ocs.ne.0).and.(i_co.ne.0).and.(i_hcl.ne.0)
	123	& .and.(i_so2.ne.0).and.(i_h2o.ne.0).and.(i_n2.ne.0)
	124	& .and.(i_co2.ne.0) ) then
	125	trac(:,1:22,i_ocs)=3.E-6
	126	trac(:,1:22,i_co)=25.E-6
	127	trac(:,:,i_hcl)=0.4E-6
	128	trac(:,1:22,i_so2)=9.E-6
	129	trac(:,1:22,i_h2o)=30.0E-6
[1305]	130
[1442]	131	c remettre tous les traceurs du start => trac(:,:,:)=trac_sav(:,:,:)
[1305]	132
[1442]	133	c N2 n est pas encore une espece chimique du modele chimique
	134	c traceur passif pour la chimie-transport
[1665]	135	trac(:,:,i_n2)=0.35d-1
[1442]	136
	137	!!!! GG: Initialization CO2 = 1 - qtot
	138	!! It assures that vmr_tot = 1
	139	c On a donc le CO2 qui est le restant d atmosphere Venus
[1665]	140	trac_sum(:,:)=0.0
[1661]	141	c SEULEMENT LES GAZ
	142	DO iq=2,nqmax-nmicro
[1442]	143	trac_sum(:,:)= trac_sum(:,:) + trac(:,:,iq)
	144	END DO
[1305]	145
[1442]	146	trac(:,:,i_co2)= 1-trac_sum(:,:)
[1665]	147	else
	148	write(,) "Réinitialisation des traceurs avec chimie "
	149	write(,) "IL manque un traceur pour les options choisies"
	150	stop
	151	endif ! verif traceurs
[1442]	152
[1305]	153	c=============================================================
	154
[1661]	155	c Passage de Rv a Rm
[1305]	156	c =============================================================
[1661]	157	c SEULEMENT LES GAZ
	158	DO iq=1,nqmax-nmicro
[1442]	159	trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/mmean(:,:)
[1305]	160	END DO
	161	c =============================================================
[1442]	162
	163	ENDIF !FIN REINIT TRAC
[1305]	164
[1661]	165	c=============================================================
	166	c == CASE INIT GAZ when muphy without chemistry ==
	167	c=============================================================
	168
	169	if (.not.ok_chem.and.ok_cloud.and.(cl_scheme.eq.2)) then
	170
	171	c Passage de Rm a Rv
	172	c =============================================================
	173	c Necessaire si on reprend les start.nc qui sont en MMR
	174	c SEULEMENT LES GAZ
	175	DO iq=1,nqmax-nmicro
	176	trac(:,:,iq)=trac(:,:,iq)*mmean(:,:)/M_tr(iq)
	177	END DO
	178	c =============================================================
	179
	180	call vapors4muphy_ini(n_lon,n_lev,trac)
	181
	182	c Passage de Rv a Rm
	183	c =============================================================
	184	c SEULEMENT LES GAZ
	185	DO iq=1,nqmax-nmicro
	186	trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/mmean(:,:)
	187	END DO
	188	c =============================================================
	189
	190	endif
[1305]	191
	192	c-------------
	193	c fin debutphy
	194	c-------------
[1442]	195
[1305]	196	ENDIF ! fin debutphy
	197
	198	c =============================================================
[1661]	199	c Passage de Rm a Rv
	200	c DES GAZ
[1305]	201	c =============================================================
[1661]	202	DO iq=1,nqmax-nmicro
	203	trac(:,:,iq)=MAX(trac(:,:,iq)*mmean(:,:)/M_tr(iq),1.E-30)
	204	END DO
[1305]	205	c =============================================================
	206
	207
	208	c =============================================================
[1661]	209	c Appel Microphysique (schema simple, these Aurelien Stolzenbach)
[1305]	210	c =============================================================
	211
[1661]	212	IF (ok_cloud .AND. cl_scheme.eq.1) THEN
	213
	214	c Passage de Rm a Rv pour les liq
	215	trac(:,:,i_h2so4liq)=MAX(trac(:,:,i_h2so4liq)
	216	& *mmean(:,:)/M_tr(i_h2so4liq),1.E-30)
	217	trac(:,:,i_h2oliq)=MAX(trac(:,:,i_h2oliq)
	218	& *mmean(:,:)/M_tr(i_h2oliq),1.E-30)
	219
[1442]	220	c PRINT*,'DEBUT CLOUD'
[1305]	221	c On remet tout le RM liq dans la partie gaz
[1661]	222	c !!! On reforme un nuage a chaque fois !!!
[1305]	223
[1442]	224	DO ilev=1, n_lev
	225	DO ilon=1, n_lon
	226	mrtwv(ilon,ilev)=trac(ilon,ilev,i_h2o) +
	227	& trac(ilon,ilev,i_h2oliq)
	228	mrtsa(ilon,ilev)=trac(ilon,ilev,i_h2so4) +
	229	& trac(ilon,ilev,i_h2so4liq)
	230	mrwv(ilon,ilev)=mrtwv(ilon,ilev)
	231	mrsa(ilon,ilev)=mrtsa(ilon,ilev)
	232	ENDDO
	233	ENDDO
	234
	235	CALL new_cloud_venus(n_lev, n_lon,
[1661]	236	e temp,pplay,
[1442]	237	e mrtwv,mrtsa,
	238	e mrwv,mrsa)
[1305]	239
	240	c =========================================
	241	c Actualisation des mixing ratio liq et gaz
	242	c =========================================
[1661]	243	c Si la routine new_cloud_venus n'a pas actualise mrwv et mrsa
[1305]	244	c on a alors bien mr=mrt pour sa et wv, donc les parties liq sont=0 hors du nuage
	245	c ou si on ne condense pas
	246
	247	c PRINT*,'DEBUT ACTUALISATION OUTPUT CLOUD'
	248	c si tout se passe bien, mrtwv et mrtsa ne changent pas
[1442]	249	DO ilev=1, n_lev
	250	DO ilon=1, n_lon
	251	trac(ilon,ilev,i_h2o) = mrwv(ilon,ilev)
	252	trac(ilon,ilev,i_h2oliq) = mrtwv(ilon,ilev) -
	253	& trac(ilon,ilev,i_h2o)
[1305]	254
[1442]	255	trac(ilon,ilev,i_h2so4) = mrsa(ilon,ilev)
	256	trac(ilon,ilev,i_h2so4liq) = mrtsa(ilon,ilev) -
	257	& trac(ilon,ilev,i_h2so4)
	258	ENDDO
	259	ENDDO
[1305]	260
	261	c =============================================================
[1661]	262	c PRINT*,'FIN CLOUD, scheme 1'
	263	ENDIF
	264
	265	c =============================================================
	266	c Appel Microphysique (schema complet, these Sabrina Guilbon)
	267	c =============================================================
	268
	269	IF (ok_cloud .AND. cl_scheme.eq.2) THEN
	270
	271	c Boucle sur grille (n_lon) et niveaux (n_lev)
	272	DO ilon=1, n_lon
	273	DO ilev=1, n_lev
	274
	275	CALL MAD_MUPHY(pdtphys, ! Timestep
	276	& temp(ilon,ilev),pplay(ilon,ilev), ! Temperature and pressure
	277	& trac(ilon,ilev,i_h2o),trac(ilon,ilev,i_h2so4), ! Mixing ratio of SA and W
	278	& trac(ilon,ilev,i_m0_aer),trac(ilon,ilev,i_m3_aer), ! Moments of aerosols
	279	& trac(ilon,ilev,i_m0_mode1drop),trac(ilon,ilev,i_m0_mode1ccn), ! Moments of mode 1
	280	& trac(ilon,ilev,i_m3_mode1sa),trac(ilon,ilev,i_m3_mode1w), ! Moments of mode 1
	281	& trac(ilon,ilev,i_m3_mode1ccn), ! Moments of mode 1
	282	& trac(ilon,ilev,i_m0_mode2drop),trac(ilon,ilev,i_m0_mode2ccn), ! Moments of mode 2
	283	& trac(ilon,ilev,i_m3_mode2sa),trac(ilon,ilev,i_m3_mode2w), ! Moments of mode 2
	284	& trac(ilon,ilev,i_m3_mode2ccn)) ! Moments of mode 2
	285
	286	ENDDO
	287	ENDDO
	288
	289
	290	c =============================================================
	291	c PRINT*,'FIN CLOUD, scheme 2'
	292	ENDIF
	293
[1442]	294
	295	c=============================================================
	296	c CHIMIE: Boucle sur les lon, lat (n_lon)
	297	c=============================================================
	298
	299	c AS:
	300	c Ici, la longitude au midi local se deplace vers l'Ouest
	301	c c'est le sens terrestre
[1661]	302	c pour Venus on prend juste l'oppose de la longitude et on a la rotation
	303	c de Venus et donc le midi local qui se deplace vers l'Est
[1442]	304
	305	lon_sun = (0.5 - gmtime) * 2.0 * RPI
	306	lon_local = lon * RPI/180.0E+0
	307	lat_local = lat * RPI/180.0E+0
	308
[1661]	309	IF (ok_chem) THEN
	310
	311	c PRINT*,'DEBUT CHEMISTRY'
[1442]	312	DO ilon=1, n_lon
	313
	314	c calcul sza_local pour obtenir des sza_local > 90, utile pour la chimie
	315	sza_local = acos(cos(lat_local(ilon))cos(lon_local(ilon))
	316	& cos(lon_sun) + cos(lat_local(ilon))*sin(lon_local(ilon))
	317	& sin(lon_sun)) 180.0E+0/RPI
[1305]	318
[1442]	319	c PRINT*,'sza_local :', sza_local
	320
[1305]	321	c =============================================================
	322	c Appel Photochimie
	323	c =============================================================
[1661]	324	c Pression en hPa => pplay/100.
[1305]	325
	326	CALL new_photochemistry_venus(n_lev, n_lon, pdtphys,
[1661]	327	e pplay(ilon,:)/100.,
[1305]	328	e temp(ilon,:),
	329	e trac(ilon,:,:),
[1452]	330	e mmean(ilon,:),
[1305]	331	e sza_local, nqmax)
	332	c =============================================================
[1661]	333
	334	END DO
[1305]	335	c PRINT*,'FIN CHEMISTRY'
[1442]	336
[1661]	337	END IF
	338	c=============================================================
[1305]	339
	340	c =============================================================
[1661]	341	c Passage de Rv a Rm
[1305]	342	c =============================================================
[1661]	343	c DES GAZ
	344	DO iq=1,nqmax-nmicro
[1442]	345	trac(:,:,iq)=trac(:,:,iq)*M_tr(iq)/mmean(:,:)
	346
[1661]	347	END DO
	348
	349	IF (ok_cloud .AND. cl_scheme.eq.1) THEN
	350	c Passage de Rm a Rv pour les liq
	351	trac(:,:,i_h2so4liq)=trac(:,:,i_h2so4liq)*M_tr(i_h2so4liq)
	352	& /mmean(:,:)
	353	trac(:,:,i_h2oliq)=trac(:,:,i_h2oliq)*M_tr(i_h2oliq)
	354	& /mmean(:,:)
	355	ENDIF
	356
[1305]	357	c =============================================================
[1442]	358	C PRINT*,'FIN PHYTRAC'
[1305]	359	RETURN
	360	END

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

Download in other formats:

Original Format