source: trunk/LMDZ.GENERIC/libf/phystd/generic_cloud_common_h.F90 @ 2725

Last change on this file since 2725 was 2725, checked in by aslmd, 2 years ago

switch to epsi_generic & RLVTT_generic

File size: 10.3 KB
Line 
1module generic_cloud_common_h
2
3    !============================================================================
4    ! Purpose
5    ! -------
6    ! Set up relevant constants and parameters for generic condensable tracers cycles,
7    ! and generic condensable tracers cloud properties
8    !
9    ! Authors
10    ! -------
11    ! Adapted from watercommon_h.F90
12    ! by Lucas Teinturier & Noé Clément (2022)
13    !============================================================================
14
15    use comcstfi_mod, only: r, cpp, mugaz
16    implicit none
17    real, save :: m  ! molecular mass of the specie (g/mol)
18    real, save :: delta_vapH ! Enthalpy of vaporization (J/mol)
19    real,save :: Tref ! Ref temperature for Clausis-Clapeyron (K)
20    real,save :: Pref ! Reference pressure for Clausius Clapeyron (Pa)
21    real, save :: epsi_generic ! fractionnal molecular mass (m/mugaz)
22    real,save :: RLVTT_generic !Latent heat of vaporization (J/kg)
23    real,save :: metallicity_coeff ! Coefficient to take into account the metallicity
24
25   
26    contains
27
28    subroutine specie_parameters(specname)
29       
30        implicit none
31
32    !============================================================================
33    !   Load the adequate set of parameters for specname
34    !   We use Clausius Clapeyron.
35    !   Values are based on Vissher (2010) and Morley (2012).
36    !
37    !   Also set up few others useful parameters, such as epsi_generic=m/mugaz, RLVTT_generic and
38    !   the metallicity_coeff.
39    !   Authors
40    !   --------
41    !   Lucas Teinturier (2022)
42    !============================================================================
43        ! Inputs
44        character(*), intent(in) :: specname
45        print*,"entree dans specie_parameters avec specname = ",specname
46        if (trim(specname) .eq. "Mg") then
47            print*,"Loading data for Mg"
48            delta_vapH = 521.7E3 !J/mol
49            Tref = 2303.0 !K
50            Pref =   1.0E5 !in Pa
51            m = 140.6931
52            RLVTT_generic = delta_vapH/(m/1000.)
53            metallicity_coeff=1.0*log(10.0)
54           
55        else if (trim(specname) .eq. "Na") then
56            print*,"Loading data for Na"
57            delta_vapH = 265.9E3
58            Tref = 1624.0
59            Pref =  1.0E5
60            m = 78.04454 !Na2S
61            RLVTT_generic = delta_vapH/(m/1000.)
62            metallicity_coeff=0.5*log(10.0)
63           
64        else if (trim(specname) .eq. "Fe") then
65            print*,"Loading data for Fe"
66            delta_vapH = 401.9E3
67            Tref = 2903.0
68            Pref  = 1.0E5
69            m = 55.8450
70            RLVTT_generic = delta_vapH/(m/1000.)
71            metallicity_coeff=0.0*log(10.0)
72           
73        else if (trim(specname) .eq. "Cr") then
74            print*,"Loading data for Cr"
75            delta_vapH = 394.2E3
76            Tref = 2749.0
77            Pref  = 1.0E5
78            m = 51.9961
79            RLVTT_generic = delta_vapH/(m/1000.)
80            metallicity_coeff=0.0*log(10.0)
81           
82        else if (trim(specname) .eq. "KCl") then
83            print*,"Loading data for KCl"
84            delta_vapH = 217.9E3
85            Tref = 1495.0
86            Pref = 1.0E5
87            metallicity_coeff=0.0*log(10.0)
88            m = 74.5498
89            RLVTT_generic = delta_vapH/(m/1000.)
90        else if (trim(specname) .eq. "Mn") then
91            print*,"Loading data for Mn"
92            delta_vapH = 455.8E3
93            Tref = 2064.0
94            Pref = 1.0E5
95            metallicity_coeff=1.0*log(10.0)
96            m = 87.003049
97            RLVTT_generic = delta_vapH/(m/1000.)
98        else if (trim(specname) .eq. "Zn") then
99            print*,"Loading data for Zn"
100            delta_vapH = 303.9E3
101            Tref = 1238.0
102            Pref = 1.0E5
103            metallicity_coeff=1.0*log(10.0)
104            m = 97.445
105            RLVTT_generic = delta_vapH/(m/1000.)
106        else
107            print*,"Unknow species (not in Mg, Fe, Na, KCl, Cr, Mn or Zn)"
108        endif
109        epsi_generic = m/mugaz
110    end subroutine specie_parameters
111
112
113    subroutine specie_parameters_table(specname)
114        use datafile_mod, only: datadir
115        implicit none
116    !============================================================================
117    !   Load the adequate set of parameters for specname
118    !   From a table of traceurs
119    !
120    !   Also set up few others useful parameters, such as epsi_generic=m/mugaz, RLVTT_generic and
121    !   the metallicity_coeff.
122    !   Authors
123    !   --------
124    !   Noé Clément (2022)
125    !============================================================================
126
127        character(*), intent(in) :: specname
128        integer :: ios
129        character(len=500):: table_traceurs_line ! table_traceurs_line lines with parameters
130        !Check if we're using newstart.F or the gcm
131        if (index(trim(datadir),"/u/lmdz/WWW/planets/") .ne. 0) then !this is the default path writter in datafile_mod
132            print*,"We're using newstart.F, we don't need to load the table_tracers_condensable"
133            return
134        endif
135        open(117,file=trim(datadir)//'/table_tracers_condensable',form='formatted',status='old')
136
137        read(117,'(A)') table_traceurs_line
138
139        do
140            read(117,'(A)', iostat=ios) table_traceurs_line
141           
142            if (index(table_traceurs_line,specname) /= 0) then
143               
144                if (index(table_traceurs_line,'delta_vapH='   ) /= 0) then
145                    read(table_traceurs_line(index(table_traceurs_line,'delta_vapH=')+len('delta_vapH='):),*) delta_vapH
146                    print*, 'delta_vapH ', delta_vapH
147                end if
148                if (index(table_traceurs_line,'Tref='   ) /= 0) then
149                    read(table_traceurs_line(index(table_traceurs_line,'Tref=')+len('Tref='):),*) Tref
150                    print*, 'Tref', Tref
151                end if
152                if (index(table_traceurs_line,'Pref='   ) /= 0) then
153                    read(table_traceurs_line(index(table_traceurs_line,'Pref=')+len('Pref='):),*) Pref
154                    print*, 'Pref', Pref
155                end if
156                if (index(table_traceurs_line,'mass='   ) /= 0) then
157                    read(table_traceurs_line(index(table_traceurs_line,'mass=')+len('mass='):),*) m
158                    print*, 'mass', m
159                end if
160                if (index(table_traceurs_line,'metallicity_coeff='   ) /= 0) then
161                   read(table_traceurs_line(index(table_traceurs_line,'metallicity_coeff=')+len('metallicity_coeff='):),*) metallicity_coeff
162                   print*, 'metallicity_coeff', metallicity_coeff
163                end if
164
165                ios=1
166            end if
167
168            if (ios /= 0) exit
169           
170        end do
171
172        if (ios .eq. -1) then ! If the file has been read but there is no data on the speci, then we have ios equal to -1
173
174            write(*,*), "ERROR : the following traceur is not referenced in the table_tracers_condensable : ", specname
175
176        end if
177       
178        close(117)
179
180        RLVTT_generic=delta_vapH/(m/1000.)
181
182        write(*,*) 'RLVTT_generic', RLVTT_generic
183
184        epsi_generic = m/mugaz
185
186    end subroutine specie_parameters_table
187
188    subroutine Psat_generic(T,p,metallicity,psat,qsat)
189        IMPLICIT NONE
190    !============================================================================
191    !   Clausius-Clapeyron relation :
192    !   d(Ln(psat))/dT = delta_vapH/(RT^2)
193    !   ->Psat = Pref * exp(-delta_vapH/R*(1/T-1/Tref))
194    !
195    !   Authors
196    !   --------
197    !   Lucas Teinturier (2022) adapted from Psat_water of Jeremy Leconte (2012)
198    !============================================================================
199        !Inputs
200        real, intent(in) :: T, p ! Temperature and pressure of the layer (in K and Pas)
201        real, intent(in) :: metallicity ! metallycity (log10)
202        !Outputs
203        real, intent(out) :: psat,qsat !saturation vapor pressure (Pa) and mass mixing ratio at saturation (kg/kg) of the layer
204
205        psat = pref*exp(-delta_vapH/(r*mugaz/1000.)*(1/T-1/Tref)-metallicity_coeff*metallicity) ! in Pa (because pref in Pa)
206
207        if (psat .gt. p) then
208            qsat = 1. ! is the maximum amount of vapor that a parcel can hold without condensation, in specific concentration.
209        else
210            qsat = epsi_generic *psat/(p-(1-epsi_generic)*psat)
211        endif
212    end subroutine Psat_generic
213
214    subroutine Lcpdqsat_generic(T,p,psat,qsat,dqsat,dlnpsat)
215        implicit none
216
217    !===============================================================================
218    !   Compute dqsat = L/cp* d(qsat)/dT and d(Ln(psat)) = L/cp*d(Ln(psat))/dT
219    !   we have d(ln(psat))/dT =  delta_vapH/R*(1/T^2) for clausius-clapeyron
220    !   r*mugaz/1000. is the perfect gaz constant in the computation of "dummy"
221    !   Authors
222    !   --------
223    !   Lucas Teinturier (2022) adapted from Lcpdqsat_water of Jeremy Leconte (2012)
224    !===============================================================================
225        ! Inputs
226        real, intent(in) :: T ! Temperature (K)
227        real, intent(in) :: p ! Pressure (Pa)
228        real, intent(in) :: psat ! Saturation vapor pressure (Pa)
229        real, intent(in) :: qsat ! Mass mixing ratio at saturation (kg/kg)
230
231        ! Outputs
232        real, intent(out) :: dqsat,dlnpsat
233
234        ! Variables
235        real :: dummy ! used to store d(ln(psat))/dT
236
237        dummy = delta_vapH/((r*mugaz/1000.)*(T**2))
238
239        if (psat .gt. p) then
240            dqsat =0.
241        else
242            dqsat = (RLVTT_generic/cpp) *qsat*(p/(p-(1-epsi_generic)*psat))*dummy
243            dlnpsat = (RLVTT_generic/cpp) * dummy
244        endif
245
246    end subroutine Lcpdqsat_generic
247
248    !==================================================================
249    subroutine Tsat_generic(p,Tsat)
250
251        implicit none
252
253    !==================================================================
254    !     Purpose
255    !     -------
256    !     Compute the saturation temperature
257    !     for a given pressure (Pa)
258    !
259    !     Authors
260    !     -------
261    !     Noé Clément (2022)
262    !
263    !==================================================================
264
265        ! input
266        real p
267   
268        ! output
269        real Tsat
270
271
272        Tsat = 1/(1/Tref - (r*mugaz/1000.)/delta_vapH*Log(p/Pref))
273
274    end subroutine Tsat_generic
275    !==================================================================
276
277end module generic_cloud_common_h
Note: See TracBrowser for help on using the repository browser.