Find Planet Surface Temperature by Semigray Approximation

(c) 2013 by Barton Paul Levenson

Language:
Basic

Dialect:
Microsoft Visual Basic Express 2010

Discussion:
This routine implements the semigray greenhouse approximation I used in my article:

Levenson BP 2011. "Planet Temperatures with Surface Cooling Parameterized." Advances in Space Research, 47, 2044–2048.

I show the code below, along with a main (console app) program that calls it to show how to use it. Atmospheric absorption of sunlight and convective cooling are parameterized. Water vapor pressure is found with an approximation to the Clausius-Clapeyron Law. The algorithm iterates until the change in temperature falls under 0.001 K. The inputs are Pc, the partial pressure of carbon dioxide; Ph, water vapor partial pressure; F, the absorbed climate flux density; and T, a first guess at the surface temperature. Tau terms refer to optical thickness, F terms to flux density. The logic assumes a planet where carbon dioxide and water vapor are the major greenhouse gases; e.g., Venus, Earth, or Mars.

Imports System.Console
Imports System.Math

Module MainMod



    ' Main program starts here.
    Sub Main()
        Dim F As Double = 239.0#
        Dim Ph As Double = 392.0#
        Dim T As Double = 255.0#

        Semi(28.3#, Ph, F, T)
        WriteLine("Ph = " & Ph.ToString("f1"))
        WriteLine("T  = " & T.ToString("f1"))

        Dim x As String = ReadLine()
    End Sub ' Main



    ' Semi does the math.  Please import System.Math in a program that
    ' uses Semi, or alter math function calls to read Math.Abs and
    ' Math.Exp.
    Sub Semi(ByVal Pc As Double, ByRef Ph As Double, _
    ByRef F As Double, ByRef T As Double)
        Const emiss As Double = 0.996#
        Const kc As Double = 0.025#
        Const kh As Double = 0.277#

        Const Q2 As Double = 0.0698#
        Const sb As Double = 0.00000005670373#
        Const T0 As Double = 288.15#

        Dim delta As Double = 99.0#

        Do
            Dim tauCO2 As Double = kc * Pc ^ 0.53#
            Dim tauH2O As Double = kh * Ph ^ 0.3#
            Dim tau As Double = tauCO2 + tauH2O
            Dim tauVis As Double = 0.354# + 0.0157# * tau

            Dim Fsi As Double = F * Exp(-tauVis)
            Dim Fco2 As Double = 0.75# * F * tauCO2
            Dim Fh2o As Double = 0.75# * F * tauH2O

            Dim Fc As Double = 0.0#
            If (Fsi * tau) >= 55.97 Then
                Fc = -22.5# + 0.402# * Fsi * tau
            End If

            Dim Fs As Double = Fsi + Fco2 + Fh2o - Fc
            Dim lastT As Double = T
            T = (Fs / (emiss * sb)) ^ 0.25#
            delta = Abs(T - lastT)

            If T > 373.15# Then
                Ph = 147896.0#
            Else
                Ph = 392.0# * Exp(Q2 * (T - T0))
            End If
        Loop Until delta < 0.001#
    End Sub ' Semi
End Module ' MainMod