From 755d54135abbf333b31cf67631b9778737b8e0a2 Mon Sep 17 00:00:00 2001 From: Nick Wogan Date: Fri, 27 Oct 2023 15:57:03 -0700 Subject: [PATCH] Added supercritical saturation to gases to try to create more physical behavior --- src/adiabat/clima_adiabat_general.f90 | 12 ++---- src/clima_saturationdata.f90 | 47 ++++++++++++++++++++++- templates/EarlyMars/species.yaml | 2 + templates/runaway_greenhouse/species.yaml | 2 + 4 files changed, 53 insertions(+), 10 deletions(-) diff --git a/src/adiabat/clima_adiabat_general.f90 b/src/adiabat/clima_adiabat_general.f90 index 5f71b2a..ae39afb 100644 --- a/src/adiabat/clima_adiabat_general.f90 +++ b/src/adiabat/clima_adiabat_general.f90 @@ -352,9 +352,7 @@ subroutine make_profile(T_surf, P_i_surf, & ! compute saturation vapor pressures P_sat = huge(1.0_dp) if (allocated(d%sp%g(i)%sat)) then - if (T_surf < d%sp%g(i)%sat%T_critical) then - P_sat = d%RH(i)*d%sp%g(i)%sat%sat_pressure(T_surf) - endif + P_sat = d%RH(i)*d%sp%g(i)%sat%sat_pressure(T_surf) endif ! determine if species are condensing, or not @@ -640,9 +638,7 @@ subroutine root_fcn(d, P, T, gout) ! compute saturation vapor pressures P_sat = huge(1.0_dp) if (allocated(d%sp%g(i)%sat)) then - if (T < d%sp%g(i)%sat%T_critical) then - P_sat = d%RH(i)*d%sp%g(i)%sat%sat_pressure(T) - endif + P_sat = d%RH(i)*d%sp%g(i)%sat%sat_pressure(T) endif if (d%sp_type(i) == CondensingSpeciesType) then @@ -703,13 +699,13 @@ subroutine mixing_ratios(d, P, T, f_i_layer, f_dry) f_moist = 0.0_dp do i = 1,d%sp%ng if (d%sp_type(i) == CondensingSpeciesType) then - f_i_layer(i) = d%RH(i)*d%sp%g(i)%sat%sat_pressure(T)/P + f_i_layer(i) = min(d%RH(i)*d%sp%g(i)%sat%sat_pressure(T)/P,1.0_dp) f_moist = f_moist + f_i_layer(i) endif enddo ! fraction of the atmosphere that is dry - f_dry = 1.0_dp - f_moist + f_dry = max(1.0_dp - f_moist, 1.0e-40_dp) ! mixing ratios of dry species do i = 1,d%sp%ng if (d%sp_type(i) == DrySpeciesType) then diff --git a/src/clima_saturationdata.f90 b/src/clima_saturationdata.f90 index 779e825..dd70f71 100644 --- a/src/clima_saturationdata.f90 +++ b/src/clima_saturationdata.f90 @@ -20,11 +20,15 @@ module clima_saturationdata real(dp) :: a_v, b_v !> Latent heat fit parameters for sublimation (no units) real(dp) :: a_s, b_s + !> Non-physical latent heat fit parameters for super-critical gas (no units) + real(dp) :: a_c, b_c contains + procedure :: latent_heat_crit procedure :: latent_heat_vap procedure :: latent_heat_sub procedure :: latent_heat + procedure :: sat_pressure_crit procedure :: sat_pressure_vap procedure :: sat_pressure_sub procedure :: sat_pressure @@ -35,6 +39,16 @@ module clima_saturationdata contains + !> Latent heat of condensation above the critical point. + !> This is non-physical. Its to approximate transitions in atmospheres + !> between super-critical and sub-critical gases. + function latent_heat_crit(self, T) result(L) + class(SaturationData), intent(inout) :: self + real(dp), intent(in) :: T !! K + real(dp) :: L !! erg/g + L = self%a_c + self%b_c*T + end function + !> Latent heat of vaporization function latent_heat_vap(self, T) result(L) class(SaturationData), intent(inout) :: self @@ -56,13 +70,29 @@ function latent_heat(self, T) result(L) class(SaturationData), intent(inout) :: self real(dp), intent(in) :: T !! K real(dp) :: L !! erg/g - if (T > self%T_triple) then + if (T >= self%T_critical) then + L = self%latent_heat_crit(T) ! non-physical approximation + elseif (T > self%T_triple .and. T < self%T_critical) then L = self%latent_heat_vap(T) else ! (T <= T_triple) then L = self%latent_heat_sub(T) endif end function + !> Saturation pressure of a super-critical gas. + !> This is non-physical. Its to approximate transitions in atmospheres + !> between super-critical and sub-critical gases. + function sat_pressure_crit(self, T) result(p_sat) + use clima_const, only: Rgas + class(SaturationData), intent(inout) :: self + real(dp), intent(in) :: T !! K + real(dp) :: p_sat !! dynes/cm2 + real(dp) :: tmp + tmp = (integral_fcn(self%a_v, self%b_v, self%T_critical) - integral_fcn(self%a_v, self%b_v, self%T_ref)) + & + (integral_fcn(self%a_c, self%b_c, T) - integral_fcn(self%a_c, self%b_c, self%T_critical)) + p_sat = self%P_ref*exp((self%mu/Rgas)*(tmp)) + end function + !> Saturation pressure over liquid function sat_pressure_vap(self, T) result(p_sat) use clima_const, only: Rgas @@ -91,7 +121,9 @@ function sat_pressure(self, T) result(p_sat) class(SaturationData), intent(inout) :: self real(dp), intent(in) :: T !! K real(dp) :: p_sat !! dynes/cm2 - if (T > self%T_triple) then + if (T >= self%T_critical) then + p_sat = self%sat_pressure_crit(T) ! non-physical + elseif (T > self%T_triple .and. T < self%T_critical) then p_sat = self%sat_pressure_vap(T) else ! (T <= T_triple) then p_sat = self%sat_pressure_sub(T) @@ -189,6 +221,17 @@ function create_SaturationData(s, name, filename, err) result(sat) sat%b_s = tmpdict%get_real('b',error = io_err) if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif + ! non-physical parameters to approximate atmospheric transitions from super-critical + ! to sub-critical states. + tmpdict => s%get_dictionary("super-critical",.true.,error = io_err) + if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif + + sat%a_c = tmpdict%get_real('a',error = io_err) + if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif + + sat%b_c = tmpdict%get_real('b',error = io_err) + if (allocated(io_err)) then; err = trim(filename)//trim(io_err%message); return; endif + ! test evaluations block real(dp) :: L(3), P(3) diff --git a/templates/EarlyMars/species.yaml b/templates/EarlyMars/species.yaml index 4eaa305..8786944 100644 --- a/templates/EarlyMars/species.yaml +++ b/templates/EarlyMars/species.yaml @@ -22,6 +22,7 @@ species: T-critical: 647.0 vaporization: {a: 2.841421e+10, b: -1.399732e+07} sublimation: {a: 2.746884e+10, b: 4.181527e+06} + super-critical: {a: 1.793161e+12, b: 0.0} note: From the NIST database - name: CO2 composition: {C: 1, O: 2} @@ -40,6 +41,7 @@ species: T-critical: 304.13 vaporization: {a: 4.656475e+09, b: -3.393595e+06} sublimation: {a: 6.564668e+09, b: -3.892217e+06} + super-critical: {a: 1.635908e+11, b: 0.0} note: From the NIST database - name: N2 composition: {N: 2} diff --git a/templates/runaway_greenhouse/species.yaml b/templates/runaway_greenhouse/species.yaml index 711e8f6..22dcc42 100644 --- a/templates/runaway_greenhouse/species.yaml +++ b/templates/runaway_greenhouse/species.yaml @@ -26,6 +26,7 @@ species: T-critical: 647.0 vaporization: {a: 2.841421e+10, b: -1.399732e+07} sublimation: {a: 2.746884e+10, b: 4.181527e+06} + super-critical: {a: 1.793161e+12, b: 0.0} note: From the NIST database - name: CO2 composition: {C: 1, O: 2} @@ -44,6 +45,7 @@ species: T-critical: 304.13 vaporization: {a: 4.656475e+09, b: -3.393595e+06} sublimation: {a: 6.564668e+09, b: -3.892217e+06} + super-critical: {a: 1.635908e+11, b: 0.0} note: From the NIST database - name: N2 composition: {N: 2}