File:Teegarden b temperature if rotates 1and2 resonance if earth like desert planet.png

1. coding: utf-8

1. 1. EBM for teegarden b, 2:3 resonance
   2. based on Bell's energy balance model

1. 1. 27.05.2023 v 0000.0001c

1. 1. you must pip install . from directory or so on Bell EBM from github
   2. https://github.com/taylorbell57/Bell_EBM
   3. you must have Python 3
   4. copy directory to your computer, go the directoty, then "pip install ."
2. 1. https://bell-ebm.readthedocs.io/en/v1.3/

import numpy as np import scipy.ndimage import matplotlib.pyplot as plt import matplotlib as mpl import astropy.constants as const import Bell_EBM as ebm import math

def simulate_ebm(groundtype1): ## Teegarden's star b albedo=0.3 star_teff=3034. star_rad=0.120 star_mass=0.097 rad=math.pow(1.16, 0.27) mass=1.16 a=0.0259 e=0.03 # 0.01 ... 0.07 argp=0.0 simulen=1000 simustepj=1000 cp_N2 = 1.039e3 # J/(kg K) cp_H2O = 4.2e6 # J/(kg K) cp_rock = .800e6 # J/(kg K) # We'll pretend the whole atmosphere absorbs and radiates P0 = const.atm.value ## quite flat, maybe greenhouse model! planet = ebm.Planet(groundtype1, rad=const.R_earth.value*rad, mass=const.M_earth.value*mass,a=a*const.au.value, e=e, argp=argp, albedo=albedo) #planet = ebm.Planet('rock', rad=const.R_earth.value*rad, mass=const.M_earth.value*mass,a=a*const.au.value, e=e, argp=argp, albedo=albedo) #planet = ebm.Planet('gas', rad=const.R_earth.value*rad, mass=const.M_earth.value*mass,a=a*const.au.value, e=e, argp=argp,cp=cp_N2, mlDepth=P0, albedo=albedo)

#star = ebm.Star() star = ebm.Star(teff=star_teff, rad=star_rad, mass=star_mass) system = ebm.System(star, planet) #system.planet.Prot = system.planet.Porb*(2./3.) system.planet.Prot = system.planet.Porb*(1./2.) Teq = np.median(system.get_teq(np.linspace(0.,system.planet.Porb,simustepj, endpoint=False))) T0 = Teq*np.ones_like(system.planet.map.values) t0 = 0 t1 = t0+system.planet.Porb*simulen dt = system.planet.Porb/simustepj times, maps = system.run_model(T0, t0, t1, dt, verbose=False) T0 = maps[-1] t0 = times[-1] t1 = t0+system.planet.Porb dt = system.planet.Porb/simustepj print("Running please wait ...") times, maps = system.run_model(T0, t0, t1, dt, verbose=False, intermediates=True) phases = system.get_phase(times) phasePeri = system.get_phase_periastron() indexPeri = np.argmin(np.abs(phases-phasePeri)) subStellLon = system.planet.orbit.get_ssp(system.planet.orbit.t_peri)[0] #mapp1=maps[indexPeri]-273.15 #mapps1=maps[indexPeri] ## take last orbit to average mapps1=maps[-simustepj:] #print(np.shape(mapps1)) mappa1=np.mean(mapps1, axis=0) mapp1=mappa1-273.15 tmean1=round(np.mean(mapp1),2) tmin1=round(np.min(mapp1),2) tmax1=round(np.max(mapp1),2) #print(" Temperature degC mean(grid) min max " , tmean1, tmin1, tmax1) return(mapp1)

1. mapp1=simulate_ebm("water")

mapp2=simulate_ebm("rock")

1. mapp3=simulate_ebm("gas")

1. mappc1=(mapp1+mapp2+mapp3)/3

plotmapp1=mapp2

def fmt(x):

   s = f"{x:.2f}"
   return rf"{s} "

plt.title("Teegarden b temperature degC if 1:2 resonance", fontsize=16) plt.xticks(fontsize=14) plt.yticks(fontsize=14)

plotmapp2 = scipy.ndimage.zoom(plotmapp1, 4)

levels1=[-300,-200,-180,-150,-120,-110,-100,-90,-80,-70,-60,-50,-40,-20,-10,-5,0,5,10,15,20,25,30,35,40,45,50,60,70,80,90,100,110,120,130,140,150,200,250,300,400,500,600,700,800,1000,2000,3000,4000]

contourcolours1=["#00003f"]

vmin1=-150 vmax1=150

cmap1="coolwarm" cmap1="jet" cmap1="bwr" cmap1="RdBu_r" cmap1="seismic" cmap1="RdYlBu_r" cmap1="Spectral_r" cmap1="hsv_r" cmap1="turbo"

1. cmap1="rainbow"
2. cmap1="gist_rainbow_r"

plt.imshow(plotmapp2,interpolation= "bicubic", cmap=cmap1, origin="lower", vmin=vmin1, vmax=vmax1, extent=[-180, 180, -90, 90]) CS=plt.contour(plotmapp2,alpha=0.5, extent=[-180, 180, -90, 90], levels=levels1 ,origin='lower', colors=contourcolours1) plt.clabel(CS, fontsize=14)

plt.show()