Context Navigation

hdo_surfex_mod.F @ 2616

Last change on this file since 2616 was 2593, checked in by mvals, 3 years ago
Mars GCM: hdo_surfex_mod.F, vdifc_mod.F: addings to account for the effect of kinetics in the fractionation by condensation of HDO at the surface MV
File size: 6.7 KB

Line
1	MODULE hdo_surfex_mod
2
3	IMPLICIT NONE
4
5	CONTAINS
6
7	subroutine hdo_surfex(ngrid,nlay,nq,ptimestep,
8	& zt,pplay,zq,pqsurf,
9	& old_h2o_vap,qsat,pdqsdif,dwatercap_dif,
10	& hdoflux)
11
12	use tracer_mod, only: igcm_h2o_vap, igcm_h2o_ice,
13	& igcm_hdo_vap, igcm_hdo_ice,
14	& qperemin
15	use surfdat_h, only: watercaptag
16	use geometry_mod, only: longitude_deg,latitude_deg
17	use comcstfi_h, only: pi
18
19	implicit none
20	c------------------------------------------------------------------
21	c Routine to compute the fluxes between air and surface
22	c for HDO, based of the fluxes for H2O
23	c L. Rossi.; M. Vals 2019
24	c------------------------------------------------------------------
25	include "callkeys.h"
26	include "microphys.h"
27	c------------------------------------------------------------------
28	c Arguments:
29	c ---------
30	c Inputs:
31	INTEGER, INTENT(IN) :: ngrid,nlay
32	INTEGER, INTENT(IN) :: nq ! nombre de traceurs
33
34	REAL, INTENT(IN) :: ptimestep ! pas de temps physique (s)
35	REAL, INTENT(IN) :: zt(ngrid,nlay) ! local value of temperature
36	REAL, INTENT(IN) :: pplay(ngrid,nlay) ! pression au milieu des couches (Pa)
37	REAL, INTENT(IN) :: zq(ngrid,nlay,nq) ! local value of tracers
38	REAL, INTENT(IN) :: pqsurf(ngrid,nq)
39	REAL, INTENT(IN) :: old_h2o_vap(ngrid) ! traceur d'eau avant
40	!traitement de l'eau (kg/kg)
41	REAL, INTENT(IN) :: qsat(ngrid) ! saturation mixing ratio
42	REAL, INTENT(IN) :: dwatercap_dif(ngrid) ! trend related to permanent ice
43	REAL, INTENT(INOUT) :: pdqsdif(ngrid,nq) ! tendance towards surface
44	! (kg/kg.s-1)
45
46	c Output:
47	REAL, INTENT(OUT) :: hdoflux(ngrid) ! value of vapour flux of HDO
48
49	c------------------------------------------------------------------
50	c Local variables:
51
52	REAL alpha(ngrid) ! equilibrium fractionation factor
53	REAL alpha_c(ngrid) ! real fractionation factor
54	REAL extrasublim ! sublimation in excess of surface ice
55	REAL tmpratio(ngrid) ! D/H ratio in flux to surf
56	REAL h2oflux(ngrid) ! value of vapour flux of H2O
57	! same sign as pdqsdif
58	REAL*8 satu(ngrid) ! Water vapor saturation ratio over ice
59	REAL zt1(ngrid),pplay1(ngrid)
60	REAL Dv,Dv_hdo ! Water/HDO vapor diffusion coefficient
61	INTEGER ig,l
62
63	REAL DoH_vap(ngrid)
64
65	c-----------------------------------------------------------------------
66	c Calculation of the fluxes for HDO
67	!! Calculation of the saturation ratio in the layer above the surface
68	satu(1:ngrid)=old_h2o_vap(1:ngrid) / qsat(1:ngrid)
69	!! Initialisation of the fractionation coefficient
70	alpha(1:ngrid)=1.
71	alpha_c(1:ngrid)=1.
72
73	DO ig=1,ngrid
74
75	h2oflux(ig) = pdqsdif(ig,igcm_h2o_ice) +
76	& dwatercap_dif(ig)
77
78	!! IF Sublimation
79	if (h2oflux(ig).le.0.) then
80
81	if (pqsurf(ig,igcm_h2o_ice).gt.qperemin) then
82	pdqsdif(ig,igcm_hdo_ice) =
83	& pdqsdif(ig,igcm_h2o_ice)*
84	& (pqsurf(ig,igcm_hdo_ice)/
85	& pqsurf(ig,igcm_h2o_ice) )
86	else
87	pdqsdif(ig,igcm_hdo_ice) = 0.
88	endif
89
90	pdqsdif(ig,igcm_hdo_ice)=
91	& max(pdqsdif(ig,igcm_hdo_ice),
92	& -pqsurf(ig,igcm_hdo_ice)/ptimestep)
93
94	hdoflux(ig) = pdqsdif(ig,igcm_hdo_ice)
95
96	if(watercaptag(ig)) then
97
98	!if we sublimate more than qsurf
99	if ((-h2oflux(ig)*ptimestep)
100	& .gt.pqsurf(ig,igcm_h2o_ice)) then
101
102	C dwatercap_dif is how much we sublimate in excess of
103	C pqsurf for H2O
104	C hdoflux(ig) is the flux of HDO from atm. to surf.
105	c The D/H of the old ice is supposed to be 5 SMOW
106	c We need D/H of the flux to be 5, so we need
107	c dwatercap_dif* 5 * 2 * 155.76e-6 (=1 SMOW)
108	hdoflux(ig)= hdoflux(ig)
109	& +(dwatercap_dif(ig)(2.155.76e-6)*5.)
110	endif
111	endif ! watercap
112
113	else ! condensation
114
115	if (hdofrac) then !do we use fractionation?
116	!! Calculation of the H2O vapor diffusion coefficient
117	Dv = 1./3. * sqrt( 8kbzzt(ig,1)/(pi*mh2o/nav) )
118	& * kbz * zt(ig,1) /
119	& ( pi * pplay(ig,1) * (molco2+molh2o)*(molco2+molh2o)
120	& * sqrt(1.+mh2o/mco2) )
121	!! Calculation of the HDO vapor diffusion coefficient
122	Dv_hdo = 1./3. * sqrt( 8kbzzt(ig,1)/(pi*mhdo/nav) )
123	& * kbz * zt(ig,1) /
124	& ( pi * pplay(ig,1) * (molco2+molhdo)*(molco2+molhdo)
125	& * sqrt(1.+mhdo/mco2) )
126	!! Calculation of the "equilibrium" fractionation coefficient
127	c alpha(ig) = exp(16288./zt(ig,1)**2.-9.34e-2)
128	alpha(ig) = exp(13525./zt(ig,1)**2.-5.59e-2) !Lamb
129	!! Calculation of the 'real' fractionnation coefficient (effect of kinetics, see Jouzel and Merlivat, 1984)
130	alpha_c(ig) = (alpha(ig)*satu(ig))/
131	& ( (alpha(ig)(Dv/Dv_hdo)(satu(ig)-1.)) + 1.)
132	else
133	alpha_c(ig) = 1.
134	endif
135
136	if (old_h2o_vap(ig).gt.qperemin) then
137	pdqsdif(ig,igcm_hdo_ice)=
138	& alpha_c(ig)pdqsdif(ig,igcm_h2o_ice)
139	& (zq(ig,1,igcm_hdo_vap)/
140	& old_h2o_vap(ig))
141	else
142	pdqsdif(ig,igcm_hdo_ice)= 0.
143	endif
144
145	if (hdofrac) then !do we use fractionation?
146	pdqsdif(ig,igcm_hdo_ice)=
147	& min( pdqsdif(ig,igcm_hdo_ice),
148	& (zq(ig,1,igcm_hdo_vap)/ptimestep) )
149	endif
150
151	hdoflux(ig)=pdqsdif(ig,igcm_hdo_ice)
152
153	endif !sublimation
154
155	ENDDO ! of DO ig=1,ngrid
156
157	c CALL WRITEDIAGFI(ngrid,'extrasublim',
158	c & 'extrasublimation',
159	c & ' ',2,tmpratio)
160	c CALL WRITEDIAGFI(ngrid,'alpha_c_s',
161	c & 'alpha_c_s',
162	c & ' ',2,alpha_c)
163
164	return
165	end subroutine hdo_surfex
166	c------------------------------------------------------------------
167
168	end module hdo_surfex_mod

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

Download in other formats:

Original Format