Context Navigation

calc_ysat.F90 @ 3094

Last change on this file since 3094 was 3090, checked in by slebonnois, 15 months ago
BdeBatz? : Cleans microphysics and makes few corrections for physics
File size: 8.5 KB

Rev	Line
[1966]	1	SUBROUTINE calc_ysat(nlon,nlay,press,temp,ysat)
	2
	3	! ==============================================================================
	4	! Purpose
	5	! -------
	6	! Compute saturation profiles (mol/mol) for chemical tracers
	7	!
	8	! Authors
	9	! -------
	10	! S. Lebonnois
	11	! -> inicondens.F in old Titan, with T,P in planetary average
[2326]	12	! J. Vatant d'Ollone
	13	! -> 2017 : Adapt to new physics, move to F90 and compute every grid point
[3090]	14	! -> 2018 : Update saturation profiles from recent litterature (Vuitton 2019)
	15	! Modified
	16	! --------
	17	! B. de Batz de Trenquelléon
	18	! -> 2022 : Update saturation profiles and correct the pressure unit
[1966]	19	! ==============================================================================
	20
	21
	22	!-----------------------------------------------------------------------
	23	! Declarations:
	24	! -------------
	25
	26	USE comchem_h, only: nkim, cnames
	27
	28	IMPLICIT NONE
	29
	30	! Arguments :
	31	! -----------
	32	INTEGER, INTENT(IN) :: nlon ! # of grid points in the chunk
	33	INTEGER, INTENT(IN) :: nlay ! # of vertical layes
[3083]	34	REAL, DIMENSION(nlon,nlay), INTENT(IN) :: press ! Mid-layers pressure (!!! mbar !!!)
[1966]	35	REAL, DIMENSION(nlon,nlay), INTENT(IN) :: temp ! Mid-layers temperature (K)
	36	REAL, DIMENSION(nlon,nlay,nkim), INTENT(OUT) :: ysat ! Saturation profiles (mol/mol)
	37
	38	! Local variables :
	39	! -----------------
	40	INTEGER :: ic
	41	REAL,DIMENSION(:,:), ALLOCATABLE :: x
	42	! -------------------------------------------------------------------------
	43
	44	ALLOCATE(x(nlon,nlay))
	45
[3090]	46	! By default, ysat = 1, meaning we do not condense
[1966]	47	ysat(:,:,:) = 1.0
	48
[3090]	49	! Main loop
	50	do ic = 1, nkim
	51
	52	! CH4 :
	53	!------
	54	if(trim(cnames(ic)).eq."CH4") then
	55	!where (temp(:,:).lt.90.65)
	56	! ysat(:,:,ic) = 10.0**(4.42507e0 - ( ( ( 1165560.7e0 / temp(:,:) &
	57	! - 115352.19e0 ) / temp(:,:) + 4055.6016e0 ) / temp(:,:) &
	58	! + 453.92414e0 ) / temp(:,:) ) / press(:,:) * 1013.25e0
	59	!elsewhere
	60	! ysat(:,:,ic) = 10.0**(3.901408e0 - ( (154567.02e0 / temp(:,:) - 1598.8512e0) &
	61	! / temp(:,:) + 437.54809e0 ) / temp(:,:)) / press(:,:) * 1013.25e0
	62	!endwhere
[1966]	63
[3090]	64	! Fray and Schmidt (2009)
	65	ysat(:,:,ic) = (1.e3 / press(:,:)) * exp(1.051e1 - 1.110e3/temp(:,:) &
	66	- 4.341e3/temp(:,:)2 + 1.035e5/temp(:,:)3 - 7.910e5/temp(:,:)**4)
[1966]	67
[3090]	68	! Forcing CH4 to 1.4% minimum
	69	where (ysat(:,:,ic).lt.0.014) ysat(:,:,ic) = 0.014
[1966]	70
[3090]	71	! C2H2 :
	72	!-------
	73	else if(trim(cnames(ic)).eq."C2H2") then
	74	!ysat(:,:,ic) = 10.0**(6.09748e0-1644.1e0/temp(:,:)+7.42346e0 &
	75	! * alog10(1.0e3/temp(:,:)) ) / press(:,:)*1013.25e0/760.0
[1966]	76
[3090]	77	! Fray and Schmidt (2009)
	78	ysat(:,:,ic) = (1.e3 / press(:,:)) * exp(1.340e1 - 2.536e3/temp(:,:))
[3083]	79
[3090]	80	! C2H4 :
	81	!-------
	82	else if(trim(cnames(ic)).eq."C2H4") then
	83	! not far from Fray and Schmidt (2009)
	84	where (temp(:,:).lt.89.0)
	85	ysat(:,:,ic) = 10.0**(1.5477e0 + (1.0e0/temp(:,:) - 0.011e0) &
	86	* (16537.0e0*(1.0e0/temp(:,:) - 0.011e0) - 1038.1e0)) &
	87	/ press(:,:) * 1.01325e0 / 760.0
	88	elsewhere (temp(:,:).lt.104.0)
	89	ysat(:,:,ic) = 10.0**(8.724e0 - 901.6e0/(temp(:,:) - 2.555e0) ) &
	90	/ press(:,:) * 1013.25e0 / 760.0
	91	elsewhere (temp(:,:).lt.120.0)
	92	ysat(:,:,ic) = 10.0**(50.79e0 - 1703.0e0/temp(:,:) - 17.141e0 &
	93	* alog10(temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0
	94	elsewhere (temp(:,:).lt.155.0)
	95	ysat(:,:,ic) = 10.0**(6.74756e0 - 585.0e0/(temp(:,:) - 18.16e0) ) &
	96	/ press(:,:) * 1013.25e0 / 760.0
	97	endwhere
	98
	99	! C2H6 :
	100	!-------
	101	else if(trim(cnames(ic)).eq."C2H6") then
	102	where (temp(:,:).lt.90.)
	103	!ysat(:,:,ic) = 10.0**(10.01e0-1085.0e0/(temp(:,:)-0.561e0) ) &
	104	! / press(:,:) * 1013.25e0 / 760.0e0
	105	! Fray and Schmidt (2009)
	106	ysat(:,:,ic) = (1.e3 / press(:,:)) * exp ( 15.11 - 2207./temp(:,:) &
	107	- 24110./(temp(:,:)2) + 7.744E5/(temp(:,:)3) &
	108	- 1.161E7/(temp(:,:)4) + 6.763E7/(temp(:,:)5))
	109	elsewhere
	110	ysat(:,:,ic) = 10.0**(5.9366e0 - 1086.17e0/temp(:,:) + 3.83464e0 &
	111	* alog10(1.0e3/temp(:,:)) ) / press(:,:)*1013.25e0/760.0
	112	endwhere
	113
	114	! CH3CCH :
	115	!---------
	116	else if((trim(cnames(ic)).eq."CH3CCH") .or. (trim(cnames(ic)).eq."CH2CCH2")) then
	117	ysat(:,:,ic) = 10.0*(2.8808e0 - 4.5e0(249.9e0 - temp(:,:)) &
[1966]	118	/ (1.15e0temp(:,:) - 37.485e0) ) / press(:,:) 1013.25e0 / 760.0e0
	119
[3090]	120	! C3H6 :
	121	!-------
	122	else if(trim(cnames(ic)).eq."C3H6") then
	123	ysat(:,:,ic) = 10.0*(7.4463e0 - 1028.5654e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
[1966]	124
[3090]	125	! C3H8 :
	126	!-------
[1966]	127	else if(trim(cnames(ic)).eq."C3H8") then
	128	ysat(:,:,ic) = 10.0*(7.217e0 - 994.30251e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
	129
[3090]	130	! C4H2 :
	131	!-------
	132	else if((trim(cnames(ic)).eq."C4H2").or.(trim(cnames(ic)).eq."C4H2s")) then
	133	ysat(:,:,ic) = 10.0**(96.26781e0 - 4651.872e0/temp(:,:) - 31.68595e0 &
[1966]	134	* alog10(temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0e0
[3090]	135
	136	! C4H4 :
	137	!-------
	138	else if(trim(cnames(ic)).eq."C4H4") then
	139	ysat(:,:,ic)= 1.0e3 * exp(9.3898e0 - 2203.57/(temp(:,:)-43.15e0) ) / press(:,:)
[1966]	140
[3090]	141	! C4H6 :
	142	!-------
	143	else if(trim(cnames(ic)).eq."C4H6") then
	144	ysat(:,:,ic)= 10.0*(2.8808e0 - 4.6e0(262.3e0 - temp(:,:)) &
[1966]	145	/(1.15e0temp(:,:) - 39.345e0) ) / press(:,:) 1013.25e0 / 760.0e0
	146
[3090]	147	! C4H10 :
	148	!--------
	149	else if(trim(cnames(ic)).eq."C4H10") then
	150	ysat(:,:,ic)= 10.0*(8.446e0 - 1461.2e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
[1966]	151
[3090]	152	! C6H2 :
	153	!-------
	154	else if(trim(cnames(ic)).eq."C6H2") then
	155	ysat(:,:,ic)= 10.0*(4.666e0 - 4956e0/temp(:,:) + 25.845e0 &
[1966]	156	alog10(1.0e3/temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0e0
	157
[3090]	158	! C8H2 :
	159	!-------
	160	else if(trim(cnames(ic)).eq."C8H2") then
	161	ysat(:,:,ic)= 10.0*(3.95e0 - 6613e0/temp(:,:) + 35.055e0 &
[1966]	162	alog10(1.0e3/temp(:,:)) ) / press(:,:) * 1013.25e0 / 760.0e0
	163
[3090]	164	! C6H6 :
	165	!-------
	166	else if(trim(cnames(ic)).eq."AC6H6") then
	167	!x = 1.0e0 - temp(:,:) / 562.2e0
	168	!ysat(:,:,ic)= 48.9e3 * exp( ( 1.33213 * x*1.5 - 6.98273 x &
	169	! - x*3 (2.62863 + 3.33399 * x*3) ) 562.2e0/temp(:,:) ) / press(:,:)
[1966]	170
[3090]	171	! Fray and Schmidt (2009)
	172	ysat(:,:,ic) = (1.e3 / press(:,:)) * exp (17.35-5663./temp(:,:))
[1966]	173
[3090]	174	! HCN :
	175	!-------
	176	else if(trim(cnames(ic)).eq."HCN") then
	177	!ysat(:,:,ic)= 10.0*(8.6165e0 - 1516.5e0/(temp(:,:) - 26.2e0) ) / press(:,:) 1013.25e0 / 760.0e0
	178
	179	! Fray and Schmidt (2009)
	180	ysat(:,:,ic) = (1.e3 / press(:,:)) * exp (1.393e1 - 3.624e3/temp(:,:) &
	181	- 1.325e5/temp(:,:)2 + 6.314e6/temp(:,:)3 - 1.128e8/temp(:,:)**4)
[1966]	182
[3090]	183	! CH3CN :
	184	!--------
	185	else if(trim(cnames(ic)).eq."CH3CN") then
	186	ysat(:,:,ic)= 10.0*(8.458e0 - 1911.7e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
	187
	188	! C2H3CN :
	189	!-------
	190	else if(trim(cnames(ic)).eq."C2H3CN") then
	191	ysat(:,:,ic)= 10.0*(9.3051e0 - 2782.21/(temp(:,:) - 51.15e0) ) / press(:,:) 1013.25e0 / 760.0e0
[1966]	192
[3090]	193	! NCCN :
	194	!-------
	195	else if(trim(cnames(ic)).eq."NCCN") then
	196	ysat(:,:,ic)= 10.0*(7.454e0 - 1832e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
	197
	198	! HC3N :
	199	!-------
	200	else if(trim(cnames(ic)).eq."HC3N") then
	201	ysat(:,:,ic)= 10.0*(7.7446e0 - 1453.5609e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
	202
	203	! C4N2 :
	204	!-------
	205	else if(trim(cnames(ic)).eq."C4N2") then
	206	ysat(:,:,ic)= 10.0*(8.269e0 - 2155.0e0/temp(:,:) ) / press(:,:) 1013.25e0 / 760.0e0
[3083]	207
[3090]	208	endif
[1966]	209	enddo
	210
	211	END SUBROUTINE calc_ysat

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

Download in other formats:

Original Format