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

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

Print tendencies in diagfi.nc and small correction

File size: 9.0 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 ! fractionnal molecular mass (m/mugaz)
22    real,save :: RLVTT !Latent heat of vaporization (J/kg)
23    real,save :: metallicity_coeff ! Coefficient to take into account the metallicity
24    contains
25
26    subroutine specie_parameters(specname)
27       
28        implicit none
29
30    !============================================================================
31    !   Load the adequate set of parameters for specname
32    !   We use Clausius Clapeyron.
33    !   Values are based on Vissher (2010) and Morley (2012).
34    !
35    !   Also set up few others useful parameters, such as epsi=m/mugaz, RLVTT and
36    !   the metallicity_coeff.
37    !   Authors
38    !   --------
39    !   Lucas Teinturier (2022)
40    !============================================================================
41        ! Inputs
42        character(*), intent(in) :: specname
43        print*,"entree dans specie_parameters avec specname = ",specname
44        if (trim(specname) .eq. "Mg") then
45            print*,"Loading data for Mg"
46            delta_vapH = 521.7E3 !J/mol
47            Tref = 2303.0 !K
48            Pref =   1.0E5 !in Pa
49            m = 140.6931
50            RLVTT = delta_vapH/(m/1000.)
51            metallicity_coeff=1.0*log(10.0)
52           
53        else if (trim(specname) .eq. "Na") then
54            print*,"Loading data for Na"
55            delta_vapH = 265.9E3
56            Tref = 1624.0
57            Pref =  1.0E5
58            m = 78.04454 !Na2S
59            RLVTT = delta_vapH/(m/1000.)
60            metallicity_coeff=0.5*log(10.0)
61           
62        else if (trim(specname) .eq. "Fe") then
63            print*,"Loading data for Fe"
64            delta_vapH = 401.9E3
65            Tref = 2903.0
66            Pref  = 1.0E5
67            m = 55.8450
68            RLVTT = delta_vapH/(m/1000.)
69            metallicity_coeff=0.0*log(10.0)
70           
71        else if (trim(specname) .eq. "Cr") then
72            print*,"Loading data for Cr"
73            delta_vapH = 394.2E3
74            Tref = 2749.0
75            Pref  = 1.0E5
76            m = 51.9961
77            RLVTT = delta_vapH/(m/1000.)
78            metallicity_coeff=0.0*log(10.0)
79           
80        else if (trim(specname) .eq. "KCl") then
81            print*,"Loading data for KCl"
82            delta_vapH = 217.9E3
83            Tref = 1495.0
84            Pref = 1.0E5
85            metallicity_coeff=0.0*log(10.0)
86            m = 74.5498
87            RLVTT = delta_vapH/(m/1000.)
88        else if (trim(specname) .eq. "Mn") then
89            print*,"Loading data for Mn"
90            delta_vapH = 455.8E3
91            Tref = 2064.0
92            Pref = 1.0E5
93            metallicity_coeff=1.0*log(10.0)
94            m = 87.003049
95            RLVTT = delta_vapH/(m/1000.)
96        else if (trim(specname) .eq. "Zn") then
97            print*,"Loading data for Zn"
98            delta_vapH = 303.9E3
99            Tref = 1238.0
100            Pref = 1.0E5
101            metallicity_coeff=1.0*log(10.0)
102            m = 97.445
103            RLVTT = delta_vapH/(m/1000.)
104        else
105            print*,"Unknow species (not in Mg, Fe, Na, KCl, Cr, Mn or Zn)"
106        endif
107        epsi = m/mugaz
108    end subroutine specie_parameters
109
110
111    subroutine specie_parameters_table(specname)
112        use datafile_mod, only: datadir
113        implicit none
114    !============================================================================
115    !   Load the adequate set of parameters for specname
116    !   From a table of traceurs
117    !
118    !   Also set up few others useful parameters, such as epsi=m/mugaz, RLVTT and
119    !   the metallicity_coeff.
120    !   Authors
121    !   --------
122    !   Noé Clément (2022)
123    !============================================================================
124
125        character(*), intent(in) :: specname
126        integer :: ios
127        character(len=500):: table_traceurs_line ! table_traceurs_line lines with parameters
128
129        open(117,file=trim(datadir)//'/table_tracers_condensable',form='formatted',status='old')
130
131        read(117,'(A)') table_traceurs_line
132
133        do
134            read(117,'(A)', iostat=ios) table_traceurs_line
135           
136            if (index(table_traceurs_line,specname) /= 0) then
137               
138                if (index(table_traceurs_line,'delta_vapH='   ) /= 0) then
139                    read(table_traceurs_line(index(table_traceurs_line,'delta_vapH=')+len('delta_vapH='):),*) delta_vapH
140                    print*, 'delta_vapH ', delta_vapH
141                end if
142                if (index(table_traceurs_line,'Tref='   ) /= 0) then
143                    read(table_traceurs_line(index(table_traceurs_line,'Tref=')+len('Tref='):),*) Tref
144                    print*, 'Tref', Tref
145                end if
146                if (index(table_traceurs_line,'Pref='   ) /= 0) then
147                    read(table_traceurs_line(index(table_traceurs_line,'Pref=')+len('Pref='):),*) Pref
148                    print*, 'Pref', Pref
149                end if
150                if (index(table_traceurs_line,'mass='   ) /= 0) then
151                    read(table_traceurs_line(index(table_traceurs_line,'mass=')+len('mass='):),*) m
152                    print*, 'mass', m
153                end if
154                if (index(table_traceurs_line,'metallicity_coeff='   ) /= 0) then
155                   read(table_traceurs_line(index(table_traceurs_line,'metallicity_coeff=')+len('metallicity_coeff='):),*) metallicity_coeff
156                   print*, 'metallicity_coeff', metallicity_coeff
157                end if
158
159                ios=1
160            end if
161
162            if (ios /= 0) exit
163           
164        end do
165
166        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
167
168            write(*,*), "ERROR : the following traceur is not referenced in the table_tracers_condensable : ", specname
169
170        end if
171       
172        close(117)
173
174        RLVTT=delta_vapH/(m/1000.)
175
176        write(*,*) 'RLVTT', RLVTT
177
178        epsi = m/mugaz
179
180    end subroutine specie_parameters_table
181
182    subroutine Psat_generic(T,p,metallicity,psat,qsat)
183        IMPLICIT NONE
184    !============================================================================
185    !   Clausius-Clapeyron relation :
186    !   d(Ln(psat))/dT = delta_vapH/(RT^2)
187    !   ->Psat = Pref * exp(-delta_vapH/R*(1/T-1/Tref))
188    !
189    !   Authors
190    !   --------
191    !   Lucas Teinturier (2022) adapted from Psat_water of Jeremy Leconte (2012)
192    !============================================================================
193        !Inputs
194        real, intent(in) :: T, p ! Temperature and pressure of the layer (in K and Pas)
195        real, intent(in) :: metallicity ! metallycity (log10)
196        !Outputs
197        real, intent(out) :: psat,qsat !saturation vapor pressure (Pa) and mass mixing ratio at saturation (kg/kg) of the layer
198
199        psat = pref*exp(-delta_vapH/(r*mugaz/1000.)*(1/T-1/Tref)-metallicity_coeff*metallicity) ! in Pa (because pref in Pa)
200
201        if (psat .gt. p) then
202            qsat = 1.
203        else
204            qsat = epsi *psat/(p-(1-epsi)*psat)
205        endif
206    end subroutine Psat_generic
207
208    subroutine Lcpdqsat_generic(T,p,psat,qsat,dqsat,dlnpsat)
209        implicit none
210
211    !===============================================================================
212    !   Compute dqsat = L/cp* d(qsat)/dT and d(Ln(psat)) = L/cp*d(Ln(psat))/dT
213    !   we have d(ln(psat))/dT =  delta_vapH/R*(1/T^2) for clausius-clapeyron
214    !   r*mugaz/1000. is the perfect gaz constant in the computation of "dummy"
215    !   Authors
216    !   --------
217    !   Lucas Teinturier (2022) adapted from Lcpdqsat_water of Jeremy Leconte (2012)
218    !===============================================================================
219        ! Inputs
220        real, intent(in) :: T ! Temperature (K)
221        real, intent(in) :: p ! Pressure (Pa)
222        real, intent(in) :: psat ! Saturation vapor pressure (Pa)
223        real, intent(in) :: qsat ! Mass mixing ratio at saturation (kg/kg)
224
225        ! Outputs
226        real, intent(out) :: dqsat,dlnpsat
227
228        ! Variables
229        real :: dummy ! used to store d(ln(psat))/dT
230
231        dummy = delta_vapH/((r*mugaz/1000.)*(T**2))
232
233        if (psat .gt. p) then
234            dqsat =0.
235        else
236            dqsat = (RLVTT/cpp) *qsat*(p/(p-(1-epsi)*psat))*dummy
237            dlnpsat = (RLVTT/cpp) * dummy
238        endif
239
240    end subroutine Lcpdqsat_generic
241
242end module generic_cloud_common_h
Note: See TracBrowser for help on using the repository browser.