Planet_models_LoFo14.py 3.32 KB
Newer Older
1
2
# Lopez & Fortney 2014 models
import numpy as np
3
import warnings
4
5
6
import astropy.units as u
from astropy import constants as const

7
8
import warnings
warnings.filterwarnings("ignore", category=RuntimeWarning)
9
10
11
# ignoring this warining is certainly not the best way, but it is currently
# thrown when one of the "in-between" time steps inside the Runge-Kutta
# integration method results in a planet which has no atmosphere left.
12

13
14

def calculate_core_radius(M_core):
15
16
17
    """ M-R relation for rock/iron Earth-like core. (no envelope)
    (see Lopez & Fortney 2014 for details.)
    """
18
19
20
    R_core = (M_core**0.25)
    return R_core

21

22
def calculate_planet_mass(M_core, fenv):
23
24
    """ Planet mass determined by core mass and atmosphere mass
    (specified in terms of envelope mass fraction [in % !]). """
25
26
27
    M_pl = M_core/(1-(fenv/100))
    return M_pl

28

29
def calculate_R_env(M_p, fenv, F_p, age, metallicity):
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
    """ Calculate planetary envelope radius using the parametrized
    results from Lopez & Fortney 2014 for planets with rocky core
    and H/He gaseous envelope.

    Parameters:
    -----------
    M_p (float): planet mass (in Earth masses)
    f_env (float): envelope mass fraction (given in percent!)
    F_p (float): bolometric flux revcieved by the planet (in units
                 of Earth's insolation)
    age (float): age of the planet (in Myr)
    metallicity (str): chose models with solar or enhanced metallicity
                       (set to "solarZ" or "enhZ")
    (check Lopez & Fortney 2014 for further details and range of model
    validity).

    Returns:
    --------
    R_env (float): radius of the envelope (in Earth radii)
    """
50
    age_exponent = {"solarZ": -0.11, "enhZ": -0.18}
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
    R_env = 2.06 * (M_p)**(-0.21) * (fenv/5)**0.59 * (F_p)**0.044 \
        * ((age/1e3)/5)**(age_exponent[metallicity])
    return R_env  # R_earth


def calculate_planet_radius(M_core, fenv, age,
                            flux_solar, metallicity):
    """ Calculate planetary radius (core + envelope) using the
    parametrized results from Lopez & Fortney 2014 for planets with
    rocky core and H/He gaseous envelope.

    Parameters:
    -----------
    M_core (float): core mass (in Earth masses)
    f_env (float): envelope mass fraction (given in percent!)
    flux_solar (float): bolometric flux revcieved by the planet
                        (in units of Earth's flux)
    age (float): age of the planet (in Myr)
    metallicity (str): chose models with solar or enhanced metallicity
                       (set to "solarZ" or "enhZ")
    (check Lopez & Fortney 2014 for further details and range of model
    validity).

    Returns:
    --------
    R_pl (float): radius of the planet (in Earth radii)
    """
78
79
80
81
82
83
84

    R_core = calculate_core_radius(M_core)
    M_pl = calculate_planet_mass(M_core, fenv)
    R_env = calculate_R_env(M_pl, fenv, flux_solar, age, metallicity)
    R_pl = R_core + R_env
    return R_pl

85

86
def density_planet(M_p, R_p):
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
    """ Calculate the mean density in cgs units using the radius and
    mass of the planet.

    Parameters:
    -----------
    M_p (float): mass of the planet (in Earth masses)
    R_p (float): radius of the planet (in Earth radii)

    Returns:
    --------
    rho (float): density in cgs units (g/ccm)
    """

    rho = (M_p*const.M_earth.cgs) \
        / (4./3*np.pi*(R_p*const.R_earth.cgs)**3)
    return rho.value