地球温暖化の真実　Appendix

Published: Nov. 25, 2024, 1:04 a.m. (UTC) / Updated: Nov. 25, 2024, 1:04 a.m. (UTC) 🔖 0 Bookmarks

👍 0 👎 0

日本語

温室効果が無い場合の地表の温度は摂氏マイナス18度の証明

始めに必要な道具を準備しよう。
シュテファン=ボルツマンの法則（Stefan–Boltzmann law）は、熱輻射により黒体から放出される電磁波のエネルギーと温度の関係を表した物理法則である。
この法則によると、熱輻射により黒体から放出されるエネルギーは絶対温度の4乗に比例する。
放射発散度を $I$、絶対温度を $T$ とすれば

$$
I=\sigma T^{4}
$$

という関係が成り立つ。ここで比例係数$\sigma$はシュテファン=ボルツマン定数と呼ばれ$\sigma = 5.670374419 \times 10^{ -8} [Wm^{-2} K^{-4}]$である。

また大気圏外の人工衛星からの観測によって、太陽から地球に降り注ぐエネルギー量の単位面積当たりの流量については太陽定数と呼ばている。$E = 1,366 [Wm^{-2}]$という値を平均的に取る。細かいことを言えば、太陽活動の11年周期が有名だが、この周期に従って$E$の値も微妙に変化するが無視できる程度だ。

太陽光のエネルギーを地球の表面が反射して宇宙に放出してしまう割合（アルベド）: $R_{albedo} = 0.3$
とすると、大気圏外まで到達した太陽エネルギーのうち地表面を熱することが出来るのは$(1-R_{albedo})E $だけだ。
地表面はこの受け取ったエネルギーに相当する温度$T_{surface}$まで熱せられ、熱せられた分だけ赤外線を放出する。これを式で表すと

$$
(1-R_{albedo})E = \sigma T_{surface}^4
$$

となり、これを$T_{surface}$について解くと

$$
T_{surface} = \sqrt[4]{ \frac{(1-R_{albedo})E}{ \sigma } }
$$

結果として以下の値を得る。

$$
T_{surface} =　-18.34 ^\circ C
$$

以下のコードでは、温室効果ガスが無い場合に推定される地球の表面温度の他に、太陽系のいくつかの惑星や衛星について計算例を示してある。

import numpy as np
import math
sigma = 5.670374419e-8 ## [W/m-2 K-4]
E = 1366 ## 太陽定数　W/m2
albedo = 0.3
T_absolute = 273.15      ##　絶対温度とセ氏の換算値 

#地球
T = pow(0.25*(1-albedo)*E/sigma,0.25)-T_absolute
print("Earth: Temperature without Green House Effect: ",T, "C  ", T+T_absolute,"K")

#天王星
albedo = 0.51
Eh = E/(19.21845*19.21845)
Th = pow(0.25*(1-albedo)*Eh/sigma,0.25)-T_absolute
print("Uranus: Temperature without Green House Effect: ",Th, "C  ", Th+T_absolute,"K")

## 海王星
albedo = 0.29
Eh = E/900
Th = pow(0.25*(1-albedo)*Eh/sigma,0.25)-T_absolute
print("Neptune: Temperature without Green House Effect: ",Th, "C  ", Th+T_absolute,"K")

#トリトン
albedo = 0.71
Eh = E/900
Th = pow(0.25*(1-albedo)*Eh/sigma,0.25)-T_absolute
print("Triton: Temperature without Green House Effect: ",Th, "C  ", Th+T_absolute,"K")

#冥王星
albedo = 0.56
Eh = E/(39.44506973*39.44506973)
Th = pow(0.25*(1-albedo)*Eh/sigma,0.25)-T_absolute
print("Pluto: Temperature without Green House Effect: ",Th, "C  ", Th+T_absolute,"K")

#ガニメデ
albedo = 0.43
Eh = E/(5.2*5.2)
Th = pow(0.25*(1-albedo)*Eh/sigma,0.25)-T_absolute
print("Ganymede: Temperature without Green House Effect: ",Th, "C  ", Th+T_absolute,"K")

上記の結果以下が得られる。

Earth: Temperature without Green House Effect:  -18.338365975615147 C   254.81163402438483 K
Uranus: Temperature without Green House Effect:  -219.98392693659153 C   53.16607306340845 K
Neptune: Temperature without Green House Effect:  -226.46270602607254 C   46.68729397392744 K
Triton: Temperature without Green House Effect:  -235.82639238987326 C   37.323607610126714 K
Pluto: Temperature without Green House Effect:  -237.02469288327376 C   36.12530711672622 K
Ganymede: Temperature without Green House Effect:  -167.00200753276977 C   106.14799246723021 K

プランクの法則：黒体放射の分布

プランクの法則は、黒体放射のスペクトルに関する法則である。
この公式から導かれるスペクトルと温度特性は、全波長領域において、熱放射の実験結果から予想される黒体放射のスペクトルと一致する。

絶対温度$T$の黒体から放射される波長$\lambda$のスペクトルのエネルギー量$I(T, \lambda)$は次の式で与えられる。

$$
I(T, \lambda)= \frac{2hc^2}{ \lambda^5} \frac{1}{ \exp( \frac{hc}{k\lambda T} ) -1 }
$$

ここで各定数のノーテーションは以下の通り。

プランク定数: $h=6.62606896 \times 10^{-34} [Js]$
ボルツマン定数： $k= 1.3806504 \times 10^{-23} [J/K] $
光速度: $c = 2.99792458 \times 10^{8}$

この関係式を太陽の表面温度の場合と、地球の表面の平均気温の場合で描画したのが以下のコードである。

import numpy as np
import math
import matplotlib.pyplot as plt
    
#################################################################
T_sol = 5780 # 太陽の表面温度　[K]
c = 2.99792458e8 # 光の速度[m/s]
h = 6.62606896e-34 #プランク定数[Js]
k = 1.3806504e-23 #ボルツマン定数[J/K] 
b = 2.897771955e-3 # [mK] ウィーンの変移則の定数　

T_array = np.array([T_sol,287]) #黒体の温度[K]　太陽と地球の表面温度の例
wave_length = np.linspace(1e-8,1e-2,300000) #観察したい波長領域[m]
################################################################

def planck(lambd, T):  # プランクの法則の計算
    x1 = 2.0 * h * (c**2) / (lambd**5)
    x2 = 1 / (np.exp(h*c/(k*lambd*T)) - 1.0)
    return x1 * x2

fig = plt.figure(figsize=(10,5))
plt.rcParams["font.size"] = 16
ax1 = plt.subplot(111)
ax1.set_title("Blackbody radiation: solar v.s. surface of Earth")
ax2 = ax1.twinx()

for T in T_array:
    B = planck(wave_length, T) 
    x = wave_length * 1e6
    if 1000< T:
        y1 = B*1e-3*1e-6
        ln1 = ax1.plot(x, y1, 'C1',label="$T={}[^\circ K] (left)$".format(T))
    else:
        y2 = B*1e-3*1e-6
        ln2 = ax2.plot(x, y2, 'C2',label="$T={}[^\circ K] (right)$".format(T))
        max_wave_length = b/T
        print("peak at:",  max_wave_length* 1e6, "[mu m] radiation of", planck(max_wave_length, T)* 1e-6, "[Wm]")
    
h1, l1 = ax1.get_legend_handles_labels()
h2, l2 = ax2.get_legend_handles_labels()
ax1.legend(h1+h2, l1+l2, loc='upper right')

ax1.set_xlim(0, 30)

ax1.set_ylabel("Solar Radiance $[W  m^{-2} \mu m^{-1} $]")
ax1.set_xlabel("Wavelength [$\mu m$]")

ax2.set_ylabel("Surface Radiance [$W m^{-2} \mu m^{-1} $]")
ax1.grid(which="both",color="gray",linewidth=0.2)

ax1.fill_between(x, y1, 0, facecolor='C1', alpha=0.5)
ax2.fill_between(x, y2, 0, facecolor='C2', alpha=0.5)
plt.show()

この結果、以下が得られる。

より高度な統合モデル

H2Oの影響がある場合$CO_2$の強制力の上限と濃度の関係

統合モデルで紹介した次の式から始めよう。

$$
T = a_0 (F_{max} (1- \exp(-\lambda \rho_{CO_2})) +F_{other}) + T_{ghe} + T_0
$$

上記モデルの前提は$CO_2$の濃度に依存して温度変化があることをモデリングするスタンスだった。ここではこれに加えて$CO_2$の濃度に依存して$H_2O$の飽和水蒸気量が増大してより多くの水分を大気が含むことが出来るようになった結果、温室効果ガスである$H_2O$が更なる温室効果を発揮して気温上昇を加速させる現象までを含めてみる。

まず、大気の温度($t[^\circ C])$と飽和水蒸気量($\alpha(t) [g/m^3]$)の関係は近似的に次の式によって与えられる。

$$
\alpha(t) = \left(\frac{217}{t+273.15}\right) \left( 6.1078 \times 10^{ \frac{7.5t}{t+237.3}} \right)
$$

水蒸気が地表面から放射されるエネルギー量を吸収することによって発揮する最大値の放射強制力を$F_{H_2O} ^{max}$とする。
これは水蒸気の吸光帯に該当する赤外線を全て水蒸気が吸収した場合に実現する仮想的な最大値という概念だ。
そして、水蒸気の飽和水蒸気量が増加するに従って、より多くのエネルギーを吸収することが出来るようになるが$F_{H_2O} ^{max}$を超えることは無い。
すると気温tの時に実現する正味の放射強制力$f_{H_2O}$は以下の式に書き下せる。

$$
f_{H_2O} = F_{H_2O} ^{ max} (1 - \exp(- \lambda_{H_2O} \alpha(t))
$$

ここで$\lambda_{H_2O}$は適当な正の定数。
以下のサンプルコードでは$F_{H_2O} ^{ max}$は大気の平均気温が$14^\circ C$の際の放射強制力の4倍と仮置きする。（この値を正確に測定するための丁度いい論文は見当たらなかっため、「地球温暖化の真実: ④地球温暖化と日光のエネルギー収支」の「温暖化ガスによる吸光」の節の図を見て凡その値で与えている。
ここまでの結果をまとめると以下が得られる。

$$
T = a_0
\left[
F_{max} \left(1- e^{-\lambda \rho_{CO_2} } \right) +
F_{H_2O} ^{ max} \left(1 - e^{ - \lambda_{H_2O} \alpha(t) } \right)
+f_{other}
\right]
+ T_{ghe} + T_0
$$

ここで$f_{other}$は$CO_2$と水蒸気以外の温暖化ガス（オゾンやメタン等）による放射強制力を表す定数。

黒体放射の変移がある場合の統合モデル

ここまでの結果は、黒体放射の分布は大気の平均気温が上昇しても不変であることを前提に進めてきた。
しかし、実際には大気の気温を上昇すると、地表面も熱せされる。熱せされた物質はその温度Tに依存して黒体放射の分布の形状を変移させる。特に地表面の温度が高くなれば、分布のピークは上昇して分布は波長の短い方向（紫外線方向）へシフトする。
すると、これまで$a_0$を定数だと仮定していたものが、実は温度が高くなるとその黒体放射のピークシフトともに増大すると考えられる。この$a_0$の温度による増分を$w$と表すと以下式が得られる。

$$
T = w \space a_0
\left[
F_{max} \left(1- e^{-\lambda \rho_{CO_2} } \right) +
F_{H_2O} ^{max} \left(1 - e^{- \lambda_{H_2O} \alpha(t) } \right)
+f_{other}
\right]
+ T_{ghe} + T_0
$$

よって黒体放射のピークシフトの影響を考えることは$w$を概算する方法を考えることに帰着する。$w$の振る舞いについての定性的な要件を次のように考える。

大気の温度が$14.454^\circ C$（1997年当時の均衡温度）の時$w = 1$で(8)式と同じとなる。この際、黒体放射のピーク値は$p_{14}$とする。
それ以外の温度Tの時は、$w=p_T/p_{14}$でピークの増減に比例して$w$が増減する。実際にはピークシフトにより紫外線方向に分布がシフトすれば、水蒸気や$CO_2$の吸光帯とのオーバーラップの仕方が変わってくる影響も考慮に入れなければならないが、その辺を簡便にモデリングするための文献は見当たらなかったため、主要な効果としてピークの増減のみを取り込むことにする。

ここで$p_T$の値はプランクの法則から直接導くことが出来る。
ウィーンの変移則の定数を$b= 2.897771955 \times 10^{-3} [mK]$とすると、温度$T$の黒体のスペクトル分布におけるピークとなる波長$\lambda_{max}$は以下の式で与えられる。

$$
\lambda_{max} = \frac{b}{T}
$$

この関係式は、地球の表面温度が変わることでその表面から放出されるエネルギー分布の形状が大きく変わるが、その変化のスケールだけをピークの放射量から概算するために用いることになる。$I(T, \lambda)$が黒体放射のスペクトル分布を与える関数であったから、以下の結果を得る。

$$
p_T = I(T, \lambda_{max})
$$

$$
w = \frac{I(T, \lambda_{max})}{I(T_{eq} , \frac{b}{T_{eq} } )}
$$
ここで$T_{eq} =- 287.604[^\circ K] (=14.454[^\circ C])$。
温室効果による気温の変化が$1 \sim 2^\circ C$と僅かであれば上式$w$は1に近似して問題無いが、そうでない場合は無視できなくなってくる。実際にこうした効果が金星を猛烈に熱い星にしてる点を思い出してほしい。

気温上昇の因果関係の整理とモデリングの際の近似の考え方

例えば水蒸気の影響を考える際に、$CO_2$濃度が上昇すると気温が$\Delta t_1$だけ上昇して、その分だけ飽和水蒸気量が増大し温室効果がより発揮され、結果として追加的に気温が$\Delta t_2$だけ上昇する。ここで注意すべき点がある。
気温が$\Delta t_2$上昇すると、更にその分だけ飽和水蒸気量が増大し温室効果が発揮され、結果として更に追加的に気温が$\Delta t_3$だけ上昇するような再帰的な因果関係は考えない。

同様に、黒体放射のピークシフトにより放射強制力が増大して結果として気温が$\Delta t_4$だけ上昇したら、その$\Delta t_4$分だけ再び黒体放射のピークシフトが起きて放射強制力が増大して・・・という再帰的な因果関係は考えない。

現実の大気においては、飽和水蒸気量の僅かな増大や黒体放射のピークシフトによる放射強制力の僅かな増大に伴う僅かな気温上昇は相互作用が起きていると考えるべきだが、その詳細をモデリングするのは分かりにくい。

よって$CO_2$濃度が$ \rho_{CO_2}$で与えられた際の、$CO_2$濃度の増大のみによる温室効果考慮後の気温$T_1$は、1997年当時均衡しいる大気の気温を$t_{eq}= 14.454[^\circ C]$とすると

$$
T_1 = a_0
\left[
F_{max} \left(1- e^{-\lambda \rho_{CO_2} } \right) +
F_{H_2O} ^{max} \left(1 - e^{- \lambda_{H_2O} \alpha(t_{eq}) } \right)
+f_{other}
\right]
+ T_{ghe} + T_0
$$

で与えられる。
ここで$t_1=T_1 - 273.15 [^\circ C]$とすると、
大気の気温が$T_1$まで上昇した際に飽和水蒸気量が増大することにより発揮される温室効果考慮後の気温$T_2$は
$$
T_2 = a_0
\left[
F_{max} \left(1- e^{-\lambda \rho_{CO_2} } \right) +
F_{H_2O} ^{max} \left(1 - e^{- \lambda_{H_2O} \alpha(t_1) } \right)
+f_{other}
\right]
+ T_{ghe} + T_0
$$

によって与えられる。
ここで更に黒体放射のピークシフトによる放射強制力全体の増大効果を考慮に入れた場合の気温$T_3$は以下の式で与えられる。

$$
T_3 = w_1 a_0
\left[
F_{max} \left(1- e^{-\lambda \rho_{CO_2} } \right) +
F_{H_2O} ^{max} \left(1 - e^{- \lambda_{H_2O} \alpha(t_1) } \right)
+f_{other}
\right]
+ T_{ghe} + T_0
$$
但し、
$$
w_1 = \frac{I(T_1, \lambda_{max})}{I(T_{eq}, \frac{b}{T_{eq}} )}, \lambda_{max} = \frac{b}{T_1}
$$
とする。

モデルの評価と課題

気温上昇幅が大きい場合のシミュレーションを行う場合には、水蒸気や黒体放射の影響は無視できない点は間違いないだろう。
しかし、以下のスクリプト生成されたモデルは実際の観測点との比較において上振れているように見える。
そして推定誤差を大きくしている原因は$T_0$を見積る際に、地球の実測値が得られないため遠くの天体を参考にしていることに大きく起因している。この点の改善が今後の最大の課題と考える。

また黒体放射の影響をより正確に考慮するためには、スペクトル分布がシフトした結果水蒸気や$CO_2$の吸光帯とのオーバーラップが変わる影響について詳細な実測データが必要になると考えられる。

##H2Oの影響がある場合

import numpy as np
import math
import matplotlib.pyplot as plt

T_absolute = 273.15 # 絶対温度 K

r = 1.49597870700e11 #地球から太陽までの距離 m

#########################################################################
T_absolute = 273.15
T0 = 3.77 #  ## 太陽が無い場合の収束温度　K　海王星の平均気温
Tghe = -18 + T_absolute  ## 温室効果が無い場合の地球表面の温度　K
Teq = 14.454 + T_absolute   ## 現在の平均気温 K  14C
Feq =105.5  ## 地球の平均放射強制力 W/m^2
#########################################################################

#########################################################################
F_rho =  (32+24)/2  # 1997年当時のCO2の放射強制力 W/m^2
Fmax = F_rho*2      # CO2の放射強制力の最大値 W/m^2

rho_1997 = 362.9          # 1997年当時のCOO2濃度　[ppm]
F_h2o_14c = (75+51)/2  # 14Cの時のH2Oの放射強制力 W/m^2
F_h2o_max = 4*F_h2o_14c # H2Oの放射強制力の最大値 W/m^2
F_other_ghg = (10+7+8+4)/2  #オゾン、メタンなどその他の温室効果ガス強制力
#########################################################################

def T_co2(rho, T_0) -> float:
    a0 = (Teq-Tghe-T_0)/Feq
    
    f = a0 * (\
            Fmax * (1- np.exp(- _lambda * rho )) + F_h2o_14c\
            + F_other_ghg  \
        ) + Tghe + T_0 - T_absolute
    return f

def f1(t) -> float: 
    return 7.5*t

def f2(t) -> float:
    return t + 237.3

def g1(t) -> float:
    return 6.1078*pow(10, f1(t)/f2(t))

def g2(t) -> float:
    return t + T_absolute

def tetens(t) -> float:
    # 温度 t ℃ における飽和水蒸気量 a(t) は次式で与えられる
    a = 217*g1(t)/g2(t) # １m3の空気中に存在できる水蒸気の質量 g
    return a

def T_co2_h2o2(rho, T_0):
    a0 = (Teq-Tghe-T_0)/Feq
    
    #t0 = T_co2(rho, T_0)  # CO2濃度の効果
    T1 = a0*(\
        Fmax * (1- np.exp(- _lambda * rho)) + \
        F_h2o_max * (1 - np.exp(- _lambda_h2o * tetens(Teq-T_absolute))) + \
        F_other_ghg  \
        ) + Tghe + T_0  # CO2濃度の効果
    
    T2 = a0*(\
        Fmax * (1- np.exp(- _lambda * rho)) + \
        F_h2o_max * (1 - np.exp(- _lambda_h2o * tetens(T1-T_absolute))) + \
        F_other_ghg  \
        ) + Tghe + T_0  #水蒸気による温室効果
        
    w = planck(b/T2, T2) / planck(b/Teq, Teq)

    T3 = w*a0*(\
        Fmax * (1- np.exp(- _lambda * rho)) + \
        F_h2o_max * (1 - np.exp(- _lambda_h2o * tetens(T1-T_absolute))) + \
        F_other_ghg  \
        ) + Tghe + T_0  # 
            
    return T3 - T_absolute

#print("Fghe: ",Fghe)

_lambda =  (1/rho_1997)*math.log(Fmax/(Fmax - F_rho)) 
# print("lambda: ",_lambda)

_lambda_h2o =  (1/tetens(14))*math.log(F_h2o_max/(F_h2o_max - F_h2o_14c))
# print("lambda H2O: ",_lambda_h2o)

rho_co2 = np.linspace(1, 2000, 100)
obsx=[280,    316.9,   362.9, 372.4, 402.9, 415.7]  ##[ppm] 1891,1960, 1997,2002,2016,2021
obsy=[13.764, 13.974, 14.454, 14.544, 14.894, 14.764] ## celsius

######### ケース別の試算 ############
T =  T_co2(rho_co2, T_0=3.77)
T1 = T_co2_h2o2(rho_co2, T_0=3.77)
T2 = T_co2_h2o2(rho_co2, T_0=14.79)

############# 描画 #############
fig = plt.figure(figsize=(10,5))
plt.rcParams["font.size"] = 16
ax = plt.subplot(111)
ax.set_title("$CO_2 and H_2O$ effect on air temperature")
ax.plot(rho_co2, T, label="$CO_2$ variable effect")
ax.plot(rho_co2, T1, linestyle=":", color = "C1",label="variable $CO_2, H_2O, radiation(T_0 = 3.77)$")
ax.plot(rho_co2, T2, linestyle="--", color = "C1",label="variable $CO_2, H_2O, radiation(T_0 = 14.79)$")

ax.plot(obsx, obsy, color="C2",label="$observation$", marker = "o")
#ax.plot(obsx2, obsy2, label="$observation2[{}^{\circ}C]$", marker = "o")
ax.set_ylabel("$Temperature [{}^{\circ} C]$")
ax.set_xlabel("Concentration of $CO_2$ [$ppm$]")
ax.grid(which="both",color="gray",linewidth=0.2)
ax.legend()

ax.fill_between(rho_co2, T1, T2, facecolor='C1', alpha=0.5)


############# いくつかの試算 ##################

p = 400
print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(p, T_co2(p,T_0=3.77)))

p = 608
t1 =T_co2_h2o2(p,T_0=3.77) 
print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(p, t1 ))
#print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(p, t2))

p = 800
print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(p, T_co2(p,T_0=3.77)))

p = 800
t1 =T_co2_h2o2(p,T_0=3.77) 
print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(p, t1 ))
#print("Earth average Temperature at CO2 level of {}[ppm]: {}[C]".format(t, t2))

この結果、以下を得る。

Earth average Temperature at CO2 level of 400[ppm]: 14.974788154154112[C]
Earth average Temperature at CO2 level of 608[ppm]: 25.005256749070895[C]
Earth average Temperature at CO2 level of 800[ppm]: 18.763382333272943[C]
Earth average Temperature at CO2 level of 800[ppm]: 30.567571184111216[C]