修士学位論文

大気大循環モデルを用いた

地球気候の太陽定数依存性に関する数値的研究

2018 年 2 月 2 日

専攻名 : 惑星学専攻学籍番号 : 163s409s 氏名 : 松田幸樹

神戸大学大学院理学研究科博士課程前期課程

要旨

惑星気候の多様性を考える上で, 太陽定数は重要なパラメータの一つである. 南北 1 次元エネルギーバランスモデルや大気大循環モデルを用いた先行研究 (例えば, Ishiwatari et al., 2002, 2007) では, 惑星気候の太陽定数依存性が調べられ, 太陽定数に応じて全球凍結解, 部分凍結解, 氷なし解, 暴走温室解や, ある太陽定数に対して多重解が得られることが示された. 上記の先行研究では, 海陸分布の影響は考慮せず, 地球とは異なる単純な大気組成を仮定した. 本研究では, 惑星大気大循環モデル DCPAM を用いて, 地球の海陸分布や大気組成を与えた場合の地球気候の太陽定数依存性を明らかにし, 大気構造の太陽定数依存性について考察した. まず, DCPAM がどの程度現在の地球の気候を表現できるかを検討した. DCPAM に地球の海表面温度分布, 海表面アルベド分布を外部データとして与えた実験では, NCEP の再解析データで見られた特徴を大まかに表現していた. DCPAM に深さ 60 m の板海を与えた実験では, DCPAM に地球の海表面温度分布, 海表面アルベド分布を外部データとして与えた実験で見られた特徴を大まかに表現していた. 次に, DCPAM に 1100 W m^{−2} から 1500 W m^{−2} までの太陽定数を与え, また複数の初期条件を与えて 18 個の実験を行った. 異なる太陽定数での年平均場の構造を調べたところ, 太陽定数が大きくなると, 年平均東西平均地表面温度の南北差は小さくなり, 年平均東西平均降水量の赤道域のピークは大きくなった. しかし, ハドレー循環の強さは太陽定数の増加に対して単調増加しなかった. このことは, Ishiwatari et al. (2002) とは異なる結果である. また, 氷線緯度と太陽定数の関係を調べると, 地球の海陸分布や大気組成を考慮した大気大循環モデルを用いた実験においても先行研究と同様に, 全球凍結解, 部分凍結解, 氷なし解が存在することが確認された. しかし, 南半球では部分凍結解, 北半球では氷なし解になるような, 氷線緯度と太陽定数の関係は北半球と南半球で異なっていた. 一方, 本研究の実験では, 同じ太陽定数において異なる氷線を持つ部分凍結解が確認された (部分凍結解の初期値依存性). そこで, この部分凍結解の初期値依存性の原因を調べるために, 南北 1 次元エネルギーバランスモデルを用いて実験したところ, 南北 1 次元エネルギーバランスモデルにおいても部分凍結解の初期値依存性が現れた. このことは, 離散モデルでは氷線が移動したとき, 氷線が隣接する格子点に到達しない限りアルベドが変化しないことに起因していると考えられる. 本研究では, サブグリッドスケールの氷の面積を考慮することで新しいアルベドの表現を考案し, それを用いることで, 南北 1 次元エネルギーバランスモデルで部分凍結解を一意に求めることができることを示した.

\underline{Mthesis 1}

目次

第 1 章はじめに 5

1.1 地球大気のエネルギーバランス . . . . . . . . . . . . . . . . . . . . 5 1.1.1 大気上端での全球平均エネルギー収支 . . . . . . . . . . . . 5 1.1.2 地球全体における単位面積あたりの全球平均エネルギー収支 6 1.1.3 エネルギー収支の緯度分布 . . . . . . . . . . . . . . . . . . . 7

1.2 数値モデルを用いた過去の研究 . . . . . . . . . . . . . . . . . . . . 7 1.2.1 南北 1 次元エネルギーバランスモデルを用いた過去の研究 . 8 1.2.2 GCM を用いた過去の研究 . . . . . . . . . . . . . . . . . . . 10

1.3 本研究の目的 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11

第 2 章大気大循環モデルの概要 12 2.1 座標系 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.2 力学過程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.3 物理過程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

第 3 章標準実験 15 3.1 実験設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

Mthesis.tex  2018/03/23(松田幸樹)

\underline{Mthesis 2}

3.1.1 系の設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.2 実験結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

3.2.1 ［気候値実験］と NCEP の再解析データとの比較 . . . . . . 20 3.2.2 ［気候値実験］と［板海実験］との比較 . . . . . . . . . . . 21

3.3 まとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

第 4 章地球気候の太陽定数増減実験 22

4.1 実験設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

4.1.1 系の設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 4.2 実験結果 – 年平均場 – . . . . . . . . . . . . . . . . . . . . . . . . . 23

4.2.1 年平均地表面温度 . . . . . . . . . . . . . . . . . . . . . . . . 24 4.2.2 年平均地表面アルベド . . . . . . . . . . . . . . . . . . . . . 25 4.2.3 年平均降水量 . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4.2.4 年平均質量流れ関数の子午面分布 . . . . . . . . . . . . . . . 26 4.2.5 年平均東西風の子午面分布 . . . . . . . . . . . . . . . . . . . 27 4.2.6 年平均東西平均温度の子午面分布 . . . . . . . . . . . . . . . 28 4.2.7 年平均東西平均雲水雲氷混合比の子午面分布 . . . . . . . . . 28 4.2.8 年平均東西平均南北熱輸送の緯度分布 . . . . . . . . . . . . 29 4.2.9 表面における年平均東西平均熱フラックスの緯度分布 …. 30

4.3 実験結果 – 氷線緯度と太陽定数の関係 – . . . . . . . . . . . . . . . 31 4.3.1 氷線緯度と太陽定数の関係 . . . . . . . . . . . . . . . . . . . 31

Mthesis.tex  2018/03/23(松田幸樹)

\underline{Mthesis 3}

4.4 考察 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 4.4.1 雲水, 雲氷の消滅時定数の変化に対する影響 . . . . . . . . . 34 4.4.2 部分凍結解の初期値依存性 . . . . . . . . . . . . . . . . . . . 35

第 5 章 EBM を用いた太陽定数増減実験 37

5.1 EBM の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 5.1.1 支配方程式系 . . . . . . . . . . . . . . . . . . . . . . . . . . 37 5.1.2 境界条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 5.1.3 支配方程式の離散化 . . . . . . . . . . . . . . . . . . . . . . 38

5.2 EBM における部分凍結解の初期値依存性 . . . . . . . . . . . . . . 39 5.2.1 実験設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 5.2.2 初期条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 5.2.3 実験結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 5.2.4 まとめと考察 . . . . . . . . . . . . . . . . . . . . . . . . . . 45

5.3 EBM におけるアルベドの表現の影響 . . . . . . . . . . . . . . . . . 46 5.3.1 実験設定 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 5.3.2 初期条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 5.3.3 実験結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 5.3.4 まとめと考察 . . . . . . . . . . . . . . . . . . . . . . . . . . 52

5.4 初期値依存性の解釈 . . . . . . . . . . . . . . . . . . . . . . . . . . 52

第 6 章結論 55

Mthesis.tex  2018/03/23(松田幸樹)

\underline{Mthesis 4} 付録 A Nakajima et al. (1992) による暴走温室状態の解釈の一部 58

A.1 暴走温室状態 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58

付録 B 全球凍結状態に関する地質学的な証拠 60 B.1 全球凍結状態 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 B.2 南極以外の世界各地での氷河堆積物の発見 . . . . . . . . . . . . . . 60 B.3 縞状鉄鉱床の存在 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

B.3.1 炭素同位体比のマイナス値 . . . . . . . . . . . . . . . . . . . 61

謝辞 62

参考文献 63

第1章はじめに

1.1 地球大気のエネルギーバランス

大気や海洋のエネルギー源は太陽からの放射エネルギーである. 地球は, その太陽からの放射エネルギーを受け, それと釣り合うように長波放射エネルギーを射出している. 本節では, 大気上端と全球平均でのエネルギー収支や放射収支の緯度分布について述べる.

1.1.1 大気上端での全球平均エネルギー収支

大気上端での全球平均エネルギー収支を考える. 図 1.1 は太陽から地球への短波放射と地球が放射する長波放射の概念図である. このように, 地球は短波放射を地球の断面積で受けるのに対し, 長波放射を全球で射出する. 地球は長期間の平均としては平衡状態にあると考えられるので, 地球の大気上端において,

(1 − α)SπR^{2} = 4πR^{2}σT_{eff}^{4}  (1.1)

となる. ここで, α は地球の全球平均アルベド, S は太陽定数, R は地球半径, σ はステファンボルツマン定数 (σ = 5.67 × 10^{−8} W m^{−2} K^{−4}), T_{eff} は有効放射温度である. 地球において, 全球平均アルベドは約 0.3, 太陽定数は約 1370 W m^{−2} であるとすると, 有効放射温度は約 255 K となる. しかし, 実際の全球平均温度は約 288 K であり, 有効放射温度より約 33 K ほど高い. これは大気の温室効果が原因である. つまり, 大気は短波放射に対してはほぼ透明であるのに対して, 長波放射に対しては透明ではなく, 地表面からの長波放射を大気が吸収し地表面と宇宙へ射出しており, 結果として有効放射温度よりも高い地表面温度で平衡になるからである.

図 1.1: 地球へ入射する太陽放射と宇宙へ射出する赤外放射の概念図 (会田, 1982, 図 1.1)

1.1.2 地球全体における単位面積あたりの全球平均エネルギー収支

全球平均・年平均した地球の単位面積あたりのエネルギー収支を図 1.2 に示す. 図 1.2 は Kiehl and Trenberth (1997) の図 7 を相対値に変更したものである (浅野, 2005a). つまり, 地球大気に入射する太陽放射エネルギーフラックスを 100% としたときの, 大気中のエネルギーフラックスを図示している. 太陽放射は, 大気上端に入射したエネルギーのうち 31% が雲や大気による散乱反射によって宇宙空間に戻る. また, 20% は大気によって吸収され直接大気を加熱する. 以上を差し引いた 49% が地表面により吸収される. 一方, 長波放射は, 大気により射出される 57% と地表面による射出が大気によって吸収されずに透過した 12% が加わり, 合計 69% の長波放射が宇宙空間へ射出されている. このように, 空気分子, 水蒸気, 雲等による散乱吸収や地表面による吸収と空気分子, 水蒸気, 雲等による射出や地表面による射出の結果, 地球全体として釣り合っている.

図 1.2: 全球平均・年平均した地球の単位面積あたりのエネルギー収支 (浅野, 2010, 図 8.3)

1.1.3 エネルギー収支の緯度分布

前小節では, 全球平均したエネルギー収支を見てきた. しかし, 地球は球形であり, 自転軸が公転面に対して傾斜しているため, 各緯度でエネルギー収支が閉じているわけではない. 図 1.3 は人工衛星観測による 1962 〜1966 年の 5 年間の時間平均および帯状平均した太陽放射の吸収量と地球放射の射出量, およびアルベドの緯度分布を表す. 放射収支を見ると, 両半球ともに赤道から緯度約 35 度までは地球放射より太陽放射の方が多く, それより高緯度側では逆に太陽放射より地球放射の方が多い. つまり, 放射過程は常に低緯度側と高緯度側の温度差を拡大する方向に作用する. しかし, 実際の地球ではこの過不足は打ち消され, エネルギー収支は釣り合っている. これは, 大気や海洋の運動により低緯度から高緯度へ熱が輸送されるからである. 図 1.4 は年平均した緯度方向のエネルギー輸送である. これを見ると, 確かに大気と海洋のどちらも赤道から極へエネルギーを輸送していることがわかる.

図 1.3: 人工衛星観測による 1962 〜 1966 年の 5 年間の時間平均および帯状平均した太陽放射の吸収量と地球放射の射出量, およびアルベドの緯度分布 (Vonder Haar and Suomi, 1971, Fig 1)

1.2 数値モデルを用いた過去の研究

前節では, 現在の地球における, 大気や海洋の循環や空気分子, 水蒸気, 雲が関わる放射過程によるエネルギー輸送の釣り合いについて述べた. 放射収支の緯度分布や大気や海洋の運動は, 入射する太陽放射量が変化すると, それに応じて変化し, 現在の地球とは異なる状態で平衡に達すると考えられる. 現在までに, 地球の気候と太陽放射の関係について数値モデルを用いて研究されてきた. 以下でその一部

図 1.4: 年平均した緯度方向のエネルギー輸送. 正味の放射強制と大気による輸送は観測から評価したものであり, 海洋の輸送はエネルギー収支の残りに応じて算出されている. (Hartmann, 1994, Fig.2.14)

を示す.

1.2.1 南北 1 次元エネルギーバランスモデルを用いた過去の研究

図 1.5: 太陽定数と氷線緯度の関係 (North et al., 1981, Fig 8 を改変)

惑星気候の太陽定数依存性は, 南北 1 次元エネルギーバランスモデル (Energy Balance Model; EBM) を用いて研究された (例えば, Budyko, 1969; Sellers, 1969).

それらの研究によると, 太陽定数に応じて, 惑星気候は全球凍結解, 部分凍結解, 氷なし解という 3 つの解を取ること, および, ある太陽定数に対して複数の状態が解となること (多重解の存在) が示された. 図 1.5 は太陽定数と氷線緯度の関係であり, 各緯度帯におけるエネルギー収支の式を用いて得られたものである (North et al., 1981). 長波放射は温度の一次関数, 水平熱輸送は拡散型で表現されている. 図 1.5 から, 3 つの解が存在し, およそ 1200 W m^{−2} から 1740 W m^{−2} で多重解が存在することが確認できる. また, 傾きが負の部分は不安定であり, 特に氷線緯度が高緯度に見られる不安定を小氷冠不安定, 低緯度に見られる不安定を大氷冠不安定という.

しかし, 長波放射を温度の一次関数で表現すると, 惑星放射の上限である射出限界が表現されない. 実際には, 射出限界があるため, 大きな太陽定数の場合には長波放射を温度の一次関数で表現する形式は正しくない (Nakajima et al., 1992).

図 1.6: EBM による太陽定数と氷線緯度の関係 (Ishiwatari et al., 2007, Fig 1). 白丸は氷なし解, 黒丸は部分凍結解, ダイヤモンドは全球凍結解, バツ印は暴走温室解である. A 点から F 点は臨界点を示し, 線分 α から ζ は解の分岐を示す.

そこで, Ishiwatari et al. (2007) は, Nakajima et al. (1992) の鉛直 1 次元放射対流モデルを用いて大気上端での惑星放射を表現した EBM を利用し, 惑星気候の太陽定数依存性を調べた (図 1.6). 水平熱輸送は, North (1975) と同様に拡散型で表現されている. その結果, 長波放射として温度の一次関数を用いた EBM と同様の 3 つの解だけでなく, 暴走温室解^{*1} も表現され, 氷なし解の上限が示された. また, 長波放射として温度の一次関数を用いた EBM と対照的に, 部分凍結解と氷なし解の結合は複雑な特徴を示した.

*1Nakajima et al. (1992) における暴走温室解の表現の一部は付録参照.

Mthesis.tex  2018/03/23(松田幸樹)

1.2.2 GCM を用いた過去の研究

惑星気候の太陽定数依存性は, 大気大循環モデル (General Circulation Model; GCM) を用いても研究された (Ishiwatari et al., 2007). その研究によると, GCM を用いても EBM で得られた 3 つの解や多重解の存在が確認された. さらに, Ishiwatari et al. (2007) では暴走温室解が存在することも確認された.

Ishiwatari et al. (2007) で用いられた GCM は AGCM 5.3 (Ishiwatari et al., 1998, 2002) ^{*2} である. 力学過程にはプリミティブ方程式系が用いられている. 放射過程は, Nakajima et al. (1992) の放射スキームが用いられている. 大気は水蒸気と乾燥空気で構成され, 水蒸気のみが長波放射を吸収し射出する. 乾燥空気は透明, 水蒸気の吸収係数は一定であり, 波長に依存しない. 雲の放射と放射の散乱は考慮されていない. 地形は考慮されていない. 海表面アルベドは 263 K を上回ると 0, 263 K を下回ると 0.5 とする. 積雪量や海氷の移動は考慮せず, 海は, 熱容量がゼロで海洋による熱輸送を考慮しない沼海として取り扱われている. 用いられた解像度は, 緯度経度格子間隔が約 5.6 度であり, 鉛直方向の大気層数は太陽定数が 1450 W m^{−2} より小さい場合は 26, 1450 W m^{−2} より大きい場合は 32 である.

図 1.7: GCM による太陽定数と氷線緯度の関係 (Ishiwatari et al., 2007, Fig 3). 白丸, 黒丸, ダイヤモンド, バツ印は図 1.6 と同様である. F は太陽定数 1000 W m^{−2} で計算される全球凍結解を初期条件とした結果, R は太陽定数 1600 W m^{−2} で計算される暴走温室解を初期条件とした結果, P は太陽定数 1300 W m^{−2} で計算される部分凍結解を初期条件とし太陽定数を徐々に減少させることで得られた結果である. 記号のない点は 280 K の等温状態を初期条件とした結果である.

図 1.7 に上記のモデルによって得られた太陽定数と氷線緯度の関係を示す. 図 1.7 から, GCM を用いても EBM で得られた 3 つの解や多重解の存在が確認された. さらに, 暴走温室解が存在することも確認された. 全球凍結解は太陽定数が 1200 W m^{−2} から 1710 W m^{−2} まで, 部分凍結解は太陽定数が 1250 W m^{−2} から 1570 W m^{−2} まで, 部分凍結解は太陽定数が 1300 W m^{−2} から 2000 W m^{−2} までの間に存在する. また, 大氷冠不安定は緯度 22^{∘} より低い緯度で起こっている. 一方で, EBM において見られた小氷冠不安定は存在しない.

*2AGCM は Atmospheric Global Climate Model の略である.

Mthesis.tex  2018/03/23(松田幸樹)

1.3 本研究の目的

前節までに, 地球では現在までに様々な気候状態が存在したと考えられていることを述べた. しかし, 数値モデルにより地球気候の太陽定数依存性を調べた先行研究である Budyko (1969) や Ishiwatari et al. (2007) では, 大気組成や海陸分布の影響を考慮していなかった. そこで, 本研究では, 地球流体電脳倶楽部^{*3} が開発している惑星大気大循環モデル DCPAM (Takahashi et al., 2010) を用いて, 地球の海陸分布や大気組成を考慮した地球気候における太陽定数と氷線緯度の関係について調べる.

本論文の以降の構成を簡単に記す. 第 2 章では, 本研究で用いた惑星大気大循環モデル DCPAM5 の概要を述べる. 第 3 章では, DCPAM5 を用いて現在の地球のパラメータや境界条件を与えた実験を行い, その結果と NCEP の再解析データを比較する. 第 4 章では, 現在の地球の海陸分布や大気組成を与えた数値モデルを用いて, 太陽定数増減実験を行う. 第 5 章では, 第 4 章で見られた部分凍結解の初期値依存性について南北 1 次元エネルギーバランスモデルを用いた考察を行う. 第 6 章では本論文の結論を述べる.

*3URL : http://www.gfd-dennou.org Mthesis.tex  2018/03/23(松田幸樹)

第2章大気大循環モデルの概要

本研究では, 地球流体電脳倶楽部で開発されている惑星大気大循環モデル DC-PAM5 を用いた. DCPAM5 では, 設定を変更することで, 地球や火星など様々な惑星の気候を計算することができる. 本研究では, 惑星半径や大気組成等に地球の値を用い, 太陽定数を変化させることで, 地球の気候の太陽定数増減実験を行った. ここでは, DCPAM5 のうち本研究で行った実験で使用したモデルの要素について, その概要を述べる. 詳細およびその他の部分については DCPAM5 のドキュメント^{*1} を参照されたい.

2.1 座標系

座標系は, 水平方向には経度 λ, 緯度 φ, 鉛直方向には σ = p/p_{s} をとる. ここで, p_{s} (λ, φ, z) は地表面気圧, p (λ, φ, z) は圧力である. ここで, z は高度である.

2.2 力学過程

大気の運動を記述する方程式系には, σ 座標系におけるプリミティブ方程式系を用いる. プリミティブ方程式系は, 運動方程式, 熱力学の式, 静水圧平衡の式, 連続の式, 水蒸気混合比の時間発展の式から成り, 以下のように書ける:

*1地球流体電脳倶楽部, 「DCPAM5 ドキュメント」 , 大気大循環モデル DCPAM5, URL :

https://www.gfd-dennou.org/library/dcpam/dcpam5/dcpam5_latest/doc/ Mthesis.tex  2018/03/23(松田幸樹)

連続の式:  \frac{∂p_{s}}{∂t} = −∇ · (p_{s}u) − \frac{∂}{∂σ} (p_{s}σ̇) ,  (2.1)
運動方程式:  \frac{∂u}{∂t} = − (u · ∇) u − σ̇ \frac{∂u}{∂σ} − fk × u + \frac{R^{d}T_{v}}{p_{s}} ∇p_{s} + F ,(2.2)

静水圧平衡の式:  \frac{∂Φ}{∂σ} = −\frac{R^{d}T_{v}}{σ} ,  (2.3)
熱力学の式:  \frac{∂T}{∂t} = − (u · ∇) T − σ̇ \frac{∂T}{∂σ} + \frac{κT_{v}}{p} \frac{dp}{dt} + Q (κ = \frac{R^{d}}{C_{p}} ) , (2.4)

状態方程式: p = ρR^{d}T_{v},  (2.5)
物質の連続の式:  \frac{∂q_{x}}{∂t} = − (u · ∇) q_{x} − σ̇ \frac{∂q_{x}}{∂σ} + S_{q_{x}} (x = v, w, i). (2.6)

ここで, u は水平速度ベクトル, σ̇ = Dσ/Dt は σ 座標系における鉛直速度, f はコリオリパラメータ, F は外力, Φ はジオポテンシャル, R^{d} は乾燥空気の気体定数, T は温度, T_{v} は仮温度, Q は非断熱加熱, ρ は密度, q は質量混合比, S_{q_{x}} は物質の移流以外による時間変化の項である. 下付き添字 v は水蒸気, w は雲水, i は雲氷を表す. また, 鉛直方向の境界条件は,

σ̇ = 0, (σ = 0, 1)  (2.7)

である.

空間の離散化には, 水平方向はスペクトル変換法 (Bourke, 1988) を用い, 鉛直方向は Arakawa and Suarez (1983) の差分法を用いる. 時間の離散化には複数の方法を組み合わせて用いる. 力学過程において, 水平拡散およびスポンジ層における減衰項には後方差分, その他の項にはセミインプリシット法を用いる. 物質移流ではセミ・ラグランジュ法を用いる.

2.3 物理過程

力学過程の方程式系における粘性項や非断熱加熱項は, 様々な物理過程を考慮して計算する. 以下でその概要を述べる.

放射過程には, 短波では δ-Eddington 近似した放射伝達方程式を Toon et al. (1989) の方法を用いて解く. H_{2}O, O_{3} の吸収係数, 雲による吸収・散乱, レイリー散乱の扱いはそれぞれ Chou and Lee (1996), Chou et al. (1998), Chou and Lee (1996) に従う. 長波では散乱を無視した放射伝達方程式により計算し, H_{2}O, CH_{4}, N_{2}O, CO_{2}, O_{3}, 雲の吸収を考慮する (Chou et al., 2001, Chou and Kouvaris, 1991).

乱流混合過程には, Mellor and Yamada (1982) Level 2.5 の方法に基づく乱流拡散係数を用いる. 惑星表面におけるフラックスの評価には Beljaars and Holtslag (1991) の方法を用いる.

凝結過程には, 積雲対流は Relaxed Arakawa-Schubert スキーム (Moorthi and Suarez, 1992) を氷晶を含むように拡張したものを用い, 大規模凝結は Le Treut and Li (1991) を氷晶を含むように拡張したものを用いる. 雲水, 雲氷の時間発展には, 移流に加えて, 乱流による混合, 凝結過程と定数時定数による消滅を考慮する.

地表面温度と土壌温度は, 地表面における熱収支と地中の熱伝導方程式を数値的に解くことによって求める. 土壌水分の決定にはバケツモデル (Manabe, 1969) を用いる. 積雪がある場合, 海氷上には雪は積もらない. 雪のアルベドは 0.6 として地面アルベドを変化させる.

海面温度の決定には次の 2 つの方法から選択して用いる. 1 つは海表面温度や海氷面密度を気候値などの外部データから与える方法であり, もう 1 つは海洋を板海としてモデル内部で海表面温度を計算する, つまり, 深さ 60 m の海洋の熱収支を計算することにより, 海表面温度を得る.

第3章標準実験

本章では, DCPAM5 がどの程度地球の気候を表現できるかを検討する. そのため, DCPAM5 に現在の地球のパラメータや境界条件を与えた実験を行い, その結果と National Centers for Environmental Prediction (NCEP) の再解析データとを比較する. ここで, DCPAM5 の実験においては, 海洋を 2 種類の方法で取り扱う. 1 つ目は海表面温度分布と海表面アルベド分布を外部データとして与える方法, 2 つ目は深さ 60 m の板海を考える方法である. 次章では, 板海を用いて海洋の温度の時間変化を考慮し, 氷アルベドフィードバックの効果を取り入れることで, 地球気候の多様性について調べる. そのため, ここでは板海を用いた場合の振る舞いを確認する. 以下では, 海表面温度分布と海表面アルベド分布を外部データとして与えて行う実験を［気候値実験］ , 深さ 60 m の板海を考える実験を［板海実験］と呼ぶ.

3.1 実験設定

実験の詳細な設定について以下に述べる.

3.1.1 系の設定

地形データは ETOPO1 (Amante and Eakins, 2009) を実験の解像度に合わせて平滑化したものを用いる (図 3.1). 陸面アルベド分布は Matthews (1983, 1984, 1985) によるものを使用する. ここで, 氷床におけるアルベドは極地ツンドラの値に変更している. これは, 次章で現在の地球とは異なる太陽定数を与えて実験する際に, 現在の地球で氷床が存在する領域に必ずしも氷床が存在するとは限らないためである. オゾン分布は CMIP5 ^{*1} 設定の気候値の東西平均を使用する. 海の取り扱いは実験ごとに異なる (表 3.1). [気候値実験] では, 海表面温度分布 (SST) は AMIPII ^{*2} 境界値から気候値を生成して与える. また, 海表面アルベド分布は AMIPII 境界値の海氷密度分布 (SIC) から気候値を生成して与える. 一方 [板海実験] では, 深さ 60 m の板海を考える. また, 海表面アルベド分布は, 海表面温度が 271.15 K を下回ると 0.6, 上回ると 0.1 を与える. 惑星パラメータは現在の地球のパラメータを用いる (表 3.2). 本研究で用いた解像度は T21L26 であり, これは緯度方向, 経度方向の格子点間隔が約 5.6 度, 鉛直方向の大気層数が 26 層に対応する. 与えた太陽定数は 1366 W m^{−2} であり, 時間ステップ, 積分時間は実験ごとに安定に積分できる値を選択する (表 3.1). 雲水の消滅時定数は 3240 秒, 雲氷の消滅時定数は 8400 秒とする.

*1 「CMIP5」は「Climate Model Intercomparison Project Phase 5」の略である. CMIP5 の詳細は Taylor et al. (2012) を参照されたい.

Mthesis.tex  2018/03/23(松田幸樹)

図 3.1: 使用した地形データ

表 3.1: 実験設定 実験名時間ステップ積分時間海の取り扱い ［気候値実験］ 32 分  12 年 SST と SIC を与えて計算 ［板海実験］  32 分  84 年深さ 60 m の板海で計算

表 3.2: 惑星パラメータ表 ; \begin{matrix}惑星パラメータ & 値 & 惑星パラメータ & 値 \\ 重力加速度 & 9.8 m s^{−2} & 自転傾斜角 & 23.44 度 \\ 惑星半径 & 6371 km & 軌道離心率 & 0.016713 \end{matrix}

*2 「AMIPII」は「Atmosphere Model Intercomparison Project II」の略である. SST と SIC の値は Taylor et al. (2000) より引用している.

Mthesis.tex  2018/03/23(松田幸樹)

初期条件として, 280 K の等温静止乾燥大気を与える.

図 3.2 は全エネルギー (運動エネルギー, ポテンシャルエネルギー, 内部エネルギー, 潜熱エネルギーの和) の全球平均値の時間変化を表す. ［気候値実験］は 4380 日 (12 年), ［板海実験］は 30660 日 (84 年) の積分で全エネルギーが経年変化しなくなっており, 統計的平衡状態に至ったと考える.

(a) ［気候値実験］ (b) ［板海実験］ ; {図 3.2: 全エネルギーの全球平均値の時間変化 (単位は J m^{−2}): (a) ［気候値実験］ , (b) \atop ［板海実験］}

3.2 実験結果

NCEP の再解析データは, 1948 年から現在まで 6 時間ごと, 日平均, 月平均データを提供している. 計算結果の比較には, その内の 1988 年から 2007 年までの月平均値を用いた. なお, この期間を選んだ理由は惑星表面条件およびオゾン分布に用いた CMIP5 のデータの期間と揃えるためである.

図 3.3, 3.4 はそれぞれ年平均した地表面温度と降水量の緯度経度分布, 図 3.5, 3.6 はそれぞれ年平均東西平均した地表面温度と降水量である. 図 3.7, 3.8, 3.9 はそれぞれ年平均した質量流れ関数, 東西風, 温度の子午面分布を示す.

(a) NCEP  (b) ［気候値実験］  (c) ［板海実験］

図 3.3: 年平均地表面温度の緯度経度分布 (図内の黒線は 271.15 K を示す): (a) NCEP, (b) ［気候値実験］ , (c) ［板海実験］

(a) NCEP  (b) ［気候値実験］  (c) ［板海実験］

図 3.4: 年平均降水量の緯度経度分布: (a) NCEP, (b) ［気候値実験］ , (c) ［板海実験］

{図 3.5: 年平均東西平均地表面温度の緯度分布. 赤線は NCEP, 緑線は［気候値実験］ , 青 \atop 線は［板海実験］を示す.}

\underline{Mthesis 19}

{図 3.6: 年平均東西平均降水量の緯度分布. 赤線は NCEP, 緑線は［気候値実験］ , 青線は \atop ［板海実験］を示す.}

(a) NCEP (b) ［気候値実験］ (c) ［板海実験］

{図 3.7: 年平均質量流れ関数の子午面分布: (a) NCEP, (b) ［気候値実験］ , (c) ［板海実 \atop 験］ )}

Mthesis.tex  2018/03/23(松田幸樹)

(a) NCEP  (b) ［気候値実験］  (c) ［板海実験］

図 3.8: 年平均東西風の子午面分布: (a) NCEP, (b) ［気候値実験］ , (c) ［板海実験］

(a) NCEP  (b) ［気候値実験］  (c) ［板海実験］

図 3.9: 年平均東西平均温度の子午面分布: (a) NCEP, (b) ［気候値実験］ , (c) ［板海実験］

3.2.1 ［気候値実験］と NCEP の再解析データとの比較

まず, ［気候値実験］の結果と NCEP の再解析データを比較する. 図 3.3 や図 3.5 を見ると, 地表面温度は赤道が最も高く, 極へ向かってなだらかに低くなるような特徴や量はおおよそ一致する. 図 3.4 や図 3.6 を見ると, 降水量は赤道と中緯度でピークを持つという特徴はおおよそ一致する. しかし, 赤道と中緯度において, 量は 1 × 10^{−5} kg m^{−2} s^{−1} 程度異なり, また赤道域を見ると, ［気候値実験］では, ピークは 1 つであるのに対し NCEP の再解析データではピークは 2 つであるという違いも見える. 図 3.7 を見ると, 質量流れ関数は 3 つの循環セルがあるという特徴はおおよそ一致する. しかし, フェレル循環は［気候値実験］の方が NCEP の再解析データよりも強い. 図 3.8 を見ると, 対流圏上層には中緯度にコアをもつ西風ジェットが存在すること, 地上付近では赤道から極へ向けて東風, 西風, 東風の分布になっているという特徴はおおよそ一致している. 図 3.9 を見ると, 対流圏界面高度, 両半球とも中緯度において水平温度勾配が急であること, また対流圏では赤道から極へ向けて温度が単調減少するという特徴はおおよそ一致している. また, 水平温度勾配が中緯度において急であることは, 温度風平衡のために, 図 3.8 において対流圏上層には中緯度にコアをもつ西風ジェットが存在することに対応していることもわかる.

3.2.2 ［気候値実験］と［板海実験］との比較

次に, ［気候値実験］の結果と［板海実験］の結果を比較する. 図 3.5 を見ると, 地表面温度は赤道が最も高く, 極へ向かってなだらかに低くなる特徴はおおよそ一致する. 図 3.6 を見ると, 降水量は赤道と中緯度でピークを持つ特徴はおおよそ一致する. しかし, 赤道において, 量は 8 × 10^{−5} kg m^{−2} s^{−1} 程度［板海実験］の方が大きい. 特に, 太平洋西の熱帯において降水量により大きな差が見える. 図 3.7 を見ると, 質量流れ関数は 3 つの循環セルがある特徴はおおよそ一致する. しかし, ハドレー循環の強さは両半球ともに［板海実験］の方が［気候値実験］よりも強い. 図 3.8 を見ると, 対流圏上層には中緯度にコアをもつ西風ジェットが存在すること, 地上付近では赤道から極へ向けて東風, 西風, 東風の分布になっている特徴はおおよそ一致している. ただし, 西風ジェットの最大値は［板海実験］の方が［気候値実験］よりも約 2 倍大きい. 図 3.9 を見ると, 対流圏界面高度, 両半球とも中緯度において水平温度勾配が急であり, また対流圏では赤道から極へ向けて温度が単調減少する特徴はおおよそ一致している. しかし, 中緯度における水平温度勾配は両半球ともに［板海実験］の方が［気候値実験］よりも大きい. その結果, 上記の西風ジェットの最大値の差が現れると考えられる.

3.3 まとめ

前節より, ［板海実験］では雲水, 雲氷の消滅時定数を正しく設定した場合に比べて地表面温度は大きくなり, それに伴い大気場も変化すると考えられる. しかし, ［板海実験］を［気候値実験］と比較した結果, 特徴はおおむね再現していると考え, 地球気候の特徴を見る上では, ［板海実験］の設定でも可能であると考える. 次章では, ［板海実験］と同様の設定で太陽定数増減実験を行う.

第4章地球気候の太陽定数増減実験

本章では, 地球気候の太陽定数依存性を調べるために, 第 2 章で示した大気大循環モデル DCPAM を用いて, 様々な太陽定数と初期条件を与えて実験を行う. 以下では, 特に 3 つの太陽定数を用いた実験について実験名を付ける. 一つ目の実験では, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1100 W m^{−2} とする. この実験を [S1100] と呼ぶ. 二つ目の実験では, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1300 W m^{−2} とする. 本章では, この実験を [S1300] と呼ぶ. 三つ目の実験では, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1475 W m^{−2} とする. この実験を [S1475] と呼ぶ.

4.1 実験設定

実験の詳細な設定について以下に述べる.

4.1.1 系の設定

地形データは ETOPO1 (Amante and Eakins, 2009) を実験の解像度に合わせて平滑化したものを用いる (図 3.1). 陸面アルベド分布は Matthews (1983, 1984, 1985) によるものを使用する. ただし, 氷床におけるアルベドは極地ツンドラの値に変更している. これは, 現在の地球とは異なる太陽定数を与えた実験において現在地球で氷床が存在する領域に必ずしも氷床が存在するとは限らないためである. オゾン分布は CMIP5 設定の気候値の東西平均を使用する. 海の取り扱いは深さ 60 m の板海を考える. また, 海表面アルベドは, 海表面温度が 271.15 K を下回ると 0.6, 上回ると 0.1 を与える. 惑星パラメータは現在の地球のパラメータを用いる (表 3.2). 本研究で用いた解像度は T21L26 であり, これは緯度方向, 経度方向の格子点間隔が約 5.6 度, 鉛直方向の大気層数が 26 層に対応する. ここまでは, 前章の［板海実験］と同様の設定である. 与えた太陽定数は 1100 W m^{−2} から 1475 W m^{−2} であり, 時間ステップ, 積分時間は実験ごとに安定に積分できる値を選択する (表 4.1). [S1100] において, 雲水の消滅時定数は 2460 秒, 雲氷の消滅時定数は 7800 秒とし, [S1300] において, 雲水の消滅時定数は 1860 秒, 雲氷の消滅時定数は 7200 秒とする. [S1475] において, 雲水の消滅時定数は 2460 秒, 雲氷の消滅時定数は 7800 秒とする. 雲水, 雲氷の消滅時定数については後ほど考察する.

表 4.1: 実験設定

実験名太陽定数初期条件時間ステップ積分時間
[S1100] 1100 W m^{−2} 280 K の等温静止乾燥大気  32 分  72 年
[S1300] 1300 W m^{−2} 280 K の等温静止乾燥大気  24 分  156 年
[S1475] 1475 W m^{−2} 280 K の等温静止乾燥大気  20 分  102 年

初期条件としては, それぞれ 280 K の等温静止乾燥大気を与える.

図 4.1 は全エネルギー (運動エネルギー, ポテンシャルエネルギー, 内部エネルギー, 潜熱エネルギーの和) の全球平均値の時間変化を表す. [S1100] は 26280 日 (72 年), [S1366] は 56940 日 (156 年), [S1475] は 37230 日 (102 年) の積分で全エネルギーが経年変化しなくなっており, 統計的平衡状態に至ったと考える.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.1: 全エネルギーの全球平均値の時間変化 (単位は J m^{−2} ): (a) [S1100], (b) [S1300], (c) [S1475]

4.2 実験結果 – 年平均場 –

ここでは, まず 3 つの実験から得られた年平均場について記述する.

4.2.1 年平均地表面温度

図 4.2 に年平均東西平均地表面温度の緯度分布を示す. 太陽定数が増加するにつれて全緯度で地表面温度は増加する. そして, [S1100] と [S1475] では南北温度差が [S1300] に比べて小さい. 図 4.3 に年平均地表面温度の緯度経度分布を示す. どの実験においても, 地表面温度は赤道が最も高く極に向かって低くなる. また, 地表面温度は東西方向にほぼ一様に分布しているが, 大陸上での平均温度が海上よりも低い. [S1300] より [S1475] の方が南北温度差が小さいのは, 太陽定数が増加すると水蒸気量が増加し温室効果による加熱と潜熱フラックスによる冷却が増えるからと考えられる (図 4.13, 図??). 一方で, [S1300] より [S1100] の方が南北温度差が小さいのは, [S1100] は全球凍結解であり, 表面アルベドの値はほぼ一定であることが原因であると考えられる (図 4.4a). 詳細については各緯度における放射収支を見る必要があり、検討が必要である.

{図 4.2: 年平均東西平均地表面温度の緯度分布. 青線は [S1100], 緑線は [S1300], 赤線は \atop [S1475] の値を示す.}

(a) [S1100]  (b) [S1300]  (c) [S1475]

{図 4.3: 年平均地表面温度の緯度経度分布 (図内の黒線は 271.15 K を示す): (a) [S1100], \atop (b) [S1300], (c) [S1475]}

4.2.2 年平均地表面アルベド

図 4.4 に年平均した地表面アルベドの緯度経度分布を示す. 地表面アルベドが 0.6 以上の部分は一年中凍っているもしくは雪が積もっていることを意味し, [S1100] は全球凍結解, [S1300] は部分凍結解であることが分かる.

[S1475] では南北半球で違いが見られる. まず南半球の高緯度に着目しよう. 南極大陸の年平均地表面アルベドは約 0.6 以上であり, 部分凍結解である. 一方で, 北半球の高緯度では, グリーンランドで年中雪が積もっているが, その他の経度では凍っていないもしくは雪が積もっていない. 図 4.2 を見ると, 年平均東西平均温度は凝固点である 271.15 K を上回っている. 本研究では, 年平均東西平均温度が凝固点である 271.15 K を上回るか否かで惑星の気候状態を判定するため, [S1475] の北半球は氷なし解とする.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.4: 年平均アルベドの緯度経度分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.3 年平均降水量

図 4.5 に年平均東西平均降水量の緯度分布を示す. 太陽定数が増加するにつれて, 赤道域での最大値は大きくなり, 全降水量は増加している. また, [S1100], [S1300] では中緯度高圧帯より高緯度側でも降水量は少ないが, [S1475] では両半球で緯度 60 付近においてもピークが見られる. 図 4.6 に年平均降水量の緯度経度分布, [S1100] では, 熱帯収束帯が見られないが, [S1300] と [S1475] では熱帯収束帯が見られる.

[S1100] と [S1300] では中緯度高圧帯より高緯度側でも降水量は少ないが, [S1475] では両半球で緯度 60 付近においてもピークが見られるのは, [S1100] では全球凍結解であり, [S1300] は緯度約 30 度に氷線を持つ部分凍結解であり, 中緯度高圧帯よりも高緯度で温度が低いことが原因であると考えられる.

{図 4.5: 年平均東西平均降水量. 青線は [S1100], 緑線は [S1300], 赤線は [S1475] の値を示 \atop す.}

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.6: 年平均降水量の緯度経度分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.4 年平均質量流れ関数の子午面分布

図 4.7 に年平均した質量流れ関数の子午面分布を示す. [S1300] と [S1475] では緯度方向に 3 つの循環セルが見られる. [S1100] では, 循環は弱いが, ハドレー循環とフェレル循環は見られる. 太陽定数が増加するにつれて, 赤道域の上昇域が南方向に移動している. また, 循環の強度は太陽定数の増加に対して単調増加するわけではない. 例えば, [S1100] と [S1300] を比較すると, [S1300] の方が循環は強いが, [S1300] と [S1475] を比較すると, 南半球の循環は [S1300] の方が循環の強さが強い. 北半球と南半球の循環の強さを比較すると, [S1100] と [S1300] は南半球の方が強いが, [S1475] は北半球の方が強い. また, ハドレー循環のセルの高さは,

[S1475] の方が高い

ハドレー循環の強さについて考察する. 図 4.5 を見ると, 太陽定数が増加すると, 赤道域での降水量は大きくなる. 従って, 降水に伴う潜熱加熱により鉛直速度が大きくなり, ハドレー循環が強くなると考えられる. しかし, [S1300] と [S1475] を比較すると, 太陽定数が増加したにも関わらず, ハドレー循環の強さは南半球では弱くなり, 北半球では変化しない. この原因は大気の安定度か地表面温度の南北温度差にあると考える. 潜熱加熱は大気の安定度と鉛直速度の積に関係している. 従って, [S1475] の赤道域の降水に伴う潜熱加熱は大気の安定度が大きいことにより釣り合っている可能性がある. また, 図 4.2 より, 地表面温度の南北温度差が小さいことから, ハドレー循環は赤道から極方向へあまり熱を輸送しないとも考えることができる. Ishiwatari et al. (2002) によると, 太陽定数が増加するにつれて, ハドレー循環の強さは強くなる. この傾向は, 本研究の結果と一致しない.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.7: 年平均質量流れ関数の子午面分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.5 年平均東西風の子午面分布

図 4.8 に年平均東西風の子午面分布を示す. どの太陽定数においても地上付近では赤道から極へ向けて東風, 西風, 東風の分布になっている. [S1300] と [S1475] において対流圏上層には中緯度にコアをもつ西風ジェットが存在するが, [S1100] では見られない. [S1300] と [S1475] を比較すると, 西風ジェットは [S1475] の方が風速が強く, 西風ジェットの最大値の高度は [S1475] の方が高い.

[S1100] で中緯度にコアを持つ西風ジェットが見られないのは, 図 4.8a より南北温度差が小さいことによる. [S1300] と [S1475] で西風ジェットの最大値が異なるのは対流圏界面の高度が [S1300] より [S1475] の方が高いからである. [S1300] は 2.5 × 10^{4} Pa で最大となり, [S1475] は 0.5 × 10^{4} Pa で最大となる. 2.5 × 10^{4} Pa における東西風は [S1300] の方が [S1475] より強く, 10 × 10^{4} Pa から 2.5 × 10^{4} Pa では, 南半球で緯度約 76 度, 北半球で氷線を持たない [S1475] より両半球で緯度約 30 度に氷線をもつ [S1300] の方が中緯度における南北温度差が大きいと考えられる.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.8: 年平均東西風の子午面分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.6 年平均東西平均温度の子午面分布

図 4.9 に年平均東西平均温度の子午面分布を示す. 図 4.9 を見ると, 太陽定数が増加するにつれて対流圏界面高度が高くなっている. また, どの太陽定数においても両半球とも, 地表面から対流圏界面では, 中緯度において水平温度勾配は最大であり, 赤道から極へ向けて温度が単調減少する特徴が見られる. また, 対流圏界面では, 赤道から中緯度にかけて温度は大きくなり, 中緯度から極にかけて温度は減少している.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.9: 年平均東西平均温度の子午面分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.7 年平均東西平均雲水雲氷混合比の子午面分布

図 4.10 に年平均東西平均雲水雲氷混合比の子午面分布を示す. 図 4.10 を見ると, 太陽定数が増加するにつれて雲水雲氷混合比は大きくなっており, 特に赤道域の対流圏中層と極域の対流圏下層で大きくなっている. 従って, [S1300] と [S1475] を比較すると, 特に赤道域の対流圏上層と極域の対流圏中層, 極域の対流圏下層で雲水雲氷混合比は大きくなるため, 惑星表面が吸収する正味の太陽放射フラックスは大きくならない (図 4.13, 図??).

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.10: 年平均東西平均雲水雲氷混合比の子午面分布: (a) [S1100], (b) [S1300], (c) [S1475]

4.2.8 年平均東西平均南北熱輸送の緯度分布

図 4.11 に年平均東西平均した乾燥静的エネルギー, 潜熱の南北熱輸送と湿潤静的エネルギーの緯度分布を示す. 図 4.11a によると, [S1100], [S1300] と [S14001500] では乾燥静的エネルギーの輸送は極向きになっている. 乾燥静的エネルギーの輸送量は低緯度と高緯度では [S1100] の場合が最も小さく, 中緯度では [S14001500] の場合が最も小さい. [S1300] での乾燥静的エネルギーの輸送量は, 北半球では低緯度から中緯度にかけてほぼ一定であり, 南半球では低緯度から高緯度にかけて徐々に減少している. [S14001500] では両半球とも低緯度と高緯度で 2 つのピークをもつ. 図 4.11b によると, [S1100] では潜熱の輸送はほとんどないことが分かる.

[S1300] と [S14001500] では潜熱の輸送は低緯度では赤道向きに高緯度では極向きになっている. ただし, [S14001500] における低緯度の輸送は収束域は [S1300] よりも南半球方向に寄っており, 北半球における赤道向きの潜熱の輸送の方が南半球における赤道向きの潜熱の輸送より大きい. 図 4.11c によると, [S1300] と [S1475] はおおよそ同じ値である.

(a) 乾燥静的エネルギー  (b) 潜熱  (c) 湿潤静的エネルギー

図 4.11: 年平均東西平均の乾燥静的エネルギー (a), 潜熱 (b) と湿潤静的エネルギー (c) の南北輸送の緯度分布. 青線は [S1100], 緑線は [S1300], 赤線は [S1475] の値を示す.

4.2.9 表面における年平均東西平均熱フラックスの緯度分布

図 4.13 に表面における年平均東西平均熱フラックスの緯度分布を示す. 図??は表面における年平均東西平均熱フラックスの緯度分布の [S1475] の結果から [S1300] の結果を引いたものである. 正味太陽放射フラックスは [S1100] と [S1300] を比較すると, [S1300] の方が大きく, 赤道域で約 90 W m^{−2}, 極域で約 4 W m^{−2} 異なる. 一方で, [S1300] と [S1475] を比較すると, [S1475] の方が大きく, 赤道域で約 4.4 W m^{−2} ∼ 33.4 W m^{−2}, 極域で約 20 W m^{−2} 異なる. 下向き長波放射フラックスと潜熱フラックスは太陽定数が増加するにつれて増加している.

表面における年平均東西平均熱フラックスの緯度分布から表面温度の緯度分布の南北差が [S1300] より [S1475] の方が小さい原因を考察する. 図??によると, 下向き長波放射フラックスは全緯度で約 200 ∼ 300 W m^{−2} 程度増加している. 中でも北半球極域と南半球中緯度は相対的に大きく, 赤道域と南半球極域では相対的に小さい. 潜熱フラックスは赤道域から中緯度にかけて大きく両半球ともに極へ向かうにつれて小さくなっている. 正味の太陽放射フラックスは中緯度で約 100 W m^{−2} の加熱があるが, 赤道域と極域では差はあまりないことがわかる. 顕熱フラックスは約 50 siW.m^{−2} の加熱である. 赤道域と北半球極域において, 表面温度に特に影響を与えるのは下向き長波放射フラックスによる加熱と潜熱フラックスによる冷却であると考えられる. 赤道域と北半球極域を比較すると, 下向き長波放射フラックスによる加熱の増加率は赤道域の方が北半球極域より小さく, 潜熱フラックスによる冷却の増加率は赤道域の方が北半球極域より大きい. 以上より赤道域の方が北半球極域よりも表面温度の増加率が小さいため, 太陽定数が増加すると表面温度の赤道域と北半球極域との差が小さくなると考えられ, そのことは上向き長波放射フラックスの変化率にも表現されている. 赤道域と北半球極域を比較すると, 下向き長波放射フラックスによる加熱の増加率は赤道域と北半球極域では同程度であり, 潜熱フラックスによる冷却の増加率は赤道域の方が北半球極域より大きい. 以上より赤道域の方が北半球極域よりも表面温度の増加率が小さいため, 太陽定数が増加すると表面温度の赤道域と南半球極域との差が小さくなると考えられるが, そのことは上向き長波放射フラックスの変化を見ると, 赤道域と南半球極域で差は見られない. 従って, 太陽定数が増加すると表面温度の赤道域と南半球極域との差についてはさらなる検証が必要であると考えられる.

(a) [S1100]  (b) [S1300]  (c) [S1475]

図 4.12: 表面における年平均東西平均熱フラックスの緯度分布: (a) [S1100], (b) [S1300], (c) [S1475]. 上向き正. 青線は正味長波放射フラックス, 水色線は上向き長波放射フラックス, ピンク線は下向き長波放射フラックス, 赤線は正味太陽放射フラックス, 緑線は潜熱フラックス, オレンジ線は顕熱フラックスを示す.

4.3 実験結果 – 氷線緯度と太陽定数の関係 –

次に, 氷線緯度と太陽定数の関係を調べるために, 太陽定数と初期条件を変更した実験を 11 個追加で行った. 各実験における初期条件, 時間ステップ, 雲水, 雲氷の消滅時定数を表 4.2 に示す. 雲水, 雲氷の消滅時定数については後ほど考察する.

表 4.2: 各実験における設定

実験番号太陽定数初期条件時間ステップ雲水の消滅時定数雲氷の消滅時定数
1  1100 W m^{−2} 280 K の等温静止乾燥大気  32 分  2460 秒  7800 秒
2  1200 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
3  1210 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
4  1250 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
5  1300 W m^{−2} 280 K の等温静止乾燥大気  24 分  1860 秒  7200 秒
6  1366 W m^{−2} 280 K の等温静止乾燥大気  32 分  3240 秒  8400 秒
7  1400 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
8  1425 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
9  1475 W m^{−2} 280 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
10  1300 W m^{−2} 200 K の等温静止乾燥大気  20 分  2460 秒  7800 秒
11  1366 W m^{−2}  実験番号 4  32 分  3240 秒  8400 秒

4.3.1 氷線緯度と太陽定数の関係

図 4.14 は実験から得られた氷線緯度と太陽定数の関係である. ここで, 氷線緯度は 271.15 K 前後の格子点の年平均東西平均地表面温度の線形補間で求めた. 図 4.14 を見ると, 全球凍結解, 部分凍結解, 氷なし解の 3 つの解が得られた. また, 氷線緯度と太陽定数の関係は南北半球で異なる. これは, 南北半球で海陸分布が異なることによると考えられる. 北半球では緯度約 18 度より低緯度で, 南半球では緯度

図 4.13: 表面における年平均東西平均熱フラックスの緯度分布の [S1475] の結果から [S1300] の結果を引いたもの. 上向き正. 青線は正味長波放射フラックス, 水色線は上向き長波放射フラックス, ピンク線は下向き長波放射フラックス, 赤線は正味太陽放射フラックス, 緑線は潜熱フラックス, オレンジ線は顕熱フラックスを示す.

約 15 度より低緯度で部分凍結解は得られていない. 一方, 北半球では緯度約 69 度より高緯度で, 南半球では緯度約 78 度より高緯度で部分凍結解は得られていない.

太陽定数が 1366 W m^{−2} のとき, 初期値に応じて異なる氷線を持つ部分凍結解が得られた. 280 K の等温静止乾燥大気を初期条件として用いると氷線は北半球で緯度約 44 度, 南半球で緯度約 55 度に現れる. また, 表 4.2 の実験番号 4 の部分凍結解を初期条件として用いると, 氷線は北半球で緯度約 37 度, 南半球で緯度約 35 度に現れる.

(a) 北半球

(b) 南半球

図 4.14: 氷線緯度と太陽定数の関係. 塗りつぶされた四角は 280 K の等温静止乾燥大気を初期条件としたものであり, 青色は全球凍結解, 緑色は部分凍結解, 赤色は氷なし解である.青色の白抜きの四角は200 K の等温静止乾燥大気を初期条件としたもの, 緑色の白抜きの四角は tbrefep2:cloud の実験番号 4 の部分凍結解を初期条件としたものである.

4.4 考察

4.4.1 雲水, 雲氷の消滅時定数の変化に対する影響

第 1 章で述べたように地球の放射収支には雲による散乱, 反射や吸収が重要である. DCPAM では, 雲水, 雲氷の時間発展は移流に加えて, 乱流による混合, 凝結過程と定数時定数による消滅を考慮している. 実験番号 1, 5 では, DCPAM で最適と考えられている雲水, 雲氷の消滅時定数よりも小さな時定数を与えて実験を実行した. ここではその影響について考察する.

雲水, 雲氷の消滅時定数は時間ステップに応じて適切な値を定めなければならない. 時間ステップと消滅時定数の適切な関係は表 4.3 である^{*1}.

表 4.3: 時間ステップに対する雲水, 雲氷の消滅時定数の最適値 時間ステップ [分]  32  24  20  10  6 雲水の消滅時定数 [秒] 3240 2720 2460 1860 1620 雲氷の消滅時定数 [秒] 8400 8000 7800 7200 6720

雲水の消滅時定数 τ_{w} を大きくすると, 雲水はより長く存在できるようになる. 例えば, 移流を無視すると定常状態では

q = S_{w}τ_{w}  (4.1)

より, 雲水の量は消滅時定数に比例する. ここで, S_{w} は雲水の移流以外による時間変化の量である. その結果, 太陽放射を反射する効果が大きくなり, また, 温室効果も大きくなる. 両者を比較すると, 太陽放射を反射する効果の方が温室効果よりも大きくなるため (田中, 2004), 地表面温度は低くなる.

同様に, 雲氷の消滅時定数 τ_{i} を大きくすると, 雲氷はより長く存在できるようになる. そして, 平衡状態では

q = S_{i}τ_{i}  (4.2)

より, 雲氷の量は消滅時定数に比例する. ここで, S_{i} は雲氷の移流以外による時間変化の量である. その結果, 太陽放射を反射する効果が大きくなり, また, 温室効果も大きくなる. しかし, 雲氷による雲は雲水による雲より高い高度にできるため, 温室効果はより大きくなる. 従って, 太陽放射を反射する効果の方が温室効果よりも小さくなるため (田中, 2004), 地表面温度は高くなる.

*1時間ステップと消滅時定数のより詳細な関係は以下を参照されたい. URL: https://www.

gfd-dennou.org/library/dcpam/sample/2016-10-21_yot/Earth/cloud_time_scale.pdf Mthesis.tex  2018/03/23(松田幸樹)

雲水と雲氷の効果を比較すると, 雲水の効果の方が雲氷の効果より大きいと考えられている (田中, 2004). 従って, 本研究で実施した [板海実験] では, 適切な雲水, 雲氷の消滅時定数を設定したときに比べて地表面温度は高くなると考えられる.

従って, 時間ステップに対して, 表 4.3 における値より小さい雲水, 雲氷の消滅時定数を設定した実験 (実験番号 1, 5) では, 地表面温度は大きくなると考えられる. その結果, 本来の氷線緯度は図 4.14 よりも赤道寄りに位置する可能性がある.

しかし, 氷線緯度と太陽定数の関係の特徴には本質的には影響はない. 上記の影響を考慮しても, 両半球ともに, 全球凍結解, 部分凍結解の 2 つの解が得られ, さらに北半球でのみ氷なし解が得られた. また, 氷線緯度と太陽定数の関係は南北半球で異なり, 北半球では緯度約 18 度より低緯度で, 南半球では緯度約 15 度より低緯度で部分凍結解は得られていない. 一方, 北半球では緯度約 69 度より高緯度で, 南半球では緯度約 78 度より高緯度で部分凍結解は得られていない.

なお, 上記の影響により, 氷線緯度がどの程度変化するのかは詳細な検討が必要である.

4.4.2 部分凍結解の初期値依存性

図 4.14 を見ると, 太陽定数 1366 W m^{−2} のとき, 南半球において, 異なる氷線を持つ部分凍結解が得られることが分かる. そこで, その 2 つの場合について, 氷線近傍における地表面温度を比較する. 以下では, 初期条件を表 4.2 の実験番号 4 の部分凍結解, 太陽定数を 1366 W m^{−2} とした実験を [S13661250] と呼ぶ. 図 4.15 は [S1366], [S13661250] から得られた氷線付近の地表面温度分布である. 図 4.15 によると, 最も低緯度の凍結した格子点が北半球では 1 格子点, 南半球では 4 格子点ずれることが分かる.

この異なる氷線を持つ部分凍結解は, これまでに行われてきた EBM を用いた実験 (North, 1981, 図 1.5, Ishiwatari et al., 2007, 図 1.6) や GCM を用いた実験 (Ishiwatari et al., 2007, 図 1.7) では存在しない.

Ishiwatari et al. (2007) と本研究を比較すると, Ishiwatari et al. (2007) では, 海の熱容量をゼロとし, 灰色大気を持つ全球海惑星を扱っているが, 本研究では海洋の熱容量を考慮し, 地球の大気組成, 海陸分布を持つ惑星を扱っている. これらの影響により, 部分凍結解の初期値依存性が現れた可能性はある. ただし, これらの違いのうち, どの違いが影響を与えるのかはより詳細な解析が必要である.

そこで, GCM に比べて過程が単純であり, 計算時間もかからない EBM を用いて部分凍結解の初期値依存性が現れるのか否かを次章で検討する.

なお, この異なる氷線を持った部分凍結解の存在により, 氷線緯度と太陽定数の関係における部分凍結解のブランチの詳細な構造はわからない. しかし, 全球凍結解, 部分凍結解, 氷なし解の存在や南北半球の差の存在は変わらない.

図 4.15: 年平均東西平均した氷線付近の地表面温度の緯度分布. 緑線は [S1366] の結果, オレンジ線は [S13661250] の結果である. また, 各格子点を + 印で表す. 黒線は 271.15 K の凝結温度を表す.

第5章{EBM を用いた太陽定数増減 \atop 実験}

部分凍結解の初期値依存性を調べるために, 南北 1 次元エネルギーバランスモデル (EBM) を用いて実験を行う. EBM は, GCM よりも過程が単純であるため部分凍結解の初期値依存性の解釈がしやすく, 計算量も少ないため様々な実験を行うことができる. 本章でははじめに, 作成した EBM について説明する. 次に, 作成した EBM を用いて部分凍結解の初期値依存性について調べる. ここでは, 水平格子点数や初期条件を変化させ, いくつかの実験を行う.

5.1 EBM の概要

本節では, 作成した EBM について詳細を説明する.

5.1.1 支配方程式系

本研究で作成したエネルギーバランスモデルは Sellers (1969) の定式化に従う. モデルの支配方程式は, 地表面温度の時間発展の式

C \frac{∂T (x, t)}{∂t} = F_{S} − F_{OLR} + \frac{∂}{∂x} D (1 − x^{2}) \frac{∂T (x, t)}{∂x},  (5.1)

で表される. ここで, C は熱容量, T は地表面温度, x は緯度の正弦, D は拡散係数である. また, 右辺第一項 F_{S} は太陽放射, 右辺第二項 F_{OLR} は惑星放射, 右辺第三項は水平熱輸送を表す. ここで, 太陽放射 F_{S} と惑星放射 F_{OLR} はそれぞれ

F_{S} = Qs (x) (1 − α (x))  (5.2)
F_{OLR} = A + BT (x, t)  (5.3)

である. ここで, Q は全球平均太陽放射フラックス, α はアルベド, A, B は定数である. 本研究では, Lindzen (1990) の値を参照し, A = 212.05 W m^{−2}, B = 1.55 W m^{−2} K^{−1} とした. s (x) は年平均太陽放射フラックスの子午面分布を表す関数であり,

s (x) = 1 + s_{2}p_{2} (x)  (5.4)

と近似する. ここで, s_{2} = −0.482, p_{2} (x) は 2 次のルジャンドル多項式で, p_{2} (x) = (3x^{2} − 1) /2 である. s_{2} の値は Lindzen (1990) の値を参照した.

5.1.2 境界条件

ここでは年平均の平衡状態について考えるため, 計算領域を半球対称とし, 赤道 (x = 0) から極 (x = 1) までの領域を考える. そして, 計算領域の両端において熱の出入りはないため, 境界条件は,

(\frac{∂T}{∂x} ) = 0, (x = 0) ,  (5.5)
(\frac{∂T}{∂x} ) = 0, (x = 1)  (5.6)

である.

5.1.3 支配方程式の離散化

図 5.1: 離散モデルにおける変数配置

ここでは, 5.1.1 節に示した支配方程式の離散化について説明する. 図 5.1 に示すように, 赤道から極までを N 個の緯度帯に分割する. (5.1) は下のように離散化する.

C \frac{T_{k}^{l+1} − T_{k}^{l}}{Δt} = Qs (x) (1 − α (x)) − (A + BT_{k}^{l})

+ \frac{D}{Δx^{2}} [(1 − x_{i+ \frac{1}{2}}^{2} ) (T_{i+1} − T_{i}) − (1 − x_{i− \frac{1}{2}}^{2} ) (T_{i} − T_{i−1})] . (5.7)

ここで, 上付き添字 l は時間を表し, 下付き添字 k は緯度帯を表す.

5.2 EBM における部分凍結解の初期値依存性

本節では, GCM で見られるような部分凍結解の初期値依存性が EBM においても見られるか否かを確認する. そのため, 解像度や初期条件を変化させ, いくつかの実験を行う.

5.2.1 実験設定

アルベド α (x_{i}) は,

α (x_{i}) = \{ \begin{matrix}α_{i} = 0.6 (T (x) ≤ 271.15 K) , \\ α_{g} = 0.1 (T (x) > 271.15 K) . \end{matrix}  (5.8)

を与える. つまり, 凝結温度は 271.15 K とした. 熱容量 C は 1.0 J K^{−1} とし, 拡散係数 D は 0.2 W m^{−2} K^{−1} とする. モデルの水平格子点数は 16, 50, 500 とする. ここで, 水平格子点数 16 は GCM における T21 に対応する. 与えた全球平均太陽放射フラックスは, 250 Wm^{−2} ≤ Q ≤ 500 Wm^{−2} である.

各実験では次に述べる初期条件から平衡に達するまでモデルを積分することで平衡状態を求める. (5.1) の右辺の絶対値が 10^{−5} W m^{−2} より小さくなると平衡状態になったと判断する.

5.2.2 初期条件

初期条件には次の 4 つの温度分布を与える (表 5.1, 図 5.2). 一つ目は, 全緯度で 250 K である全球凍結状態の温度分布 (図 5.2a, 青線), 二つ目は, x < 0.40 において 300 K, x ≥ 0.40 において 250 K である部分凍結状態の温度分布 (図 5.2b, 緑線), 三つ目は, x < 0.72 において 300 K, x ≥ 0.72 において 250 K である部分凍結状態の温度分布 (図 5.2b, 青線), 四つ目は, 全緯度で 300 K である氷なし状態の温度分布である (図 5.2a, 赤線). 以下では, 一つ目の初期条件を用いた実験を [GIC250], 二つ目の初期条件を用いた実験を [PIC6_300-250], 三つ目の初期条件を用いた実験を [PIC11_300-250], 四つ目の初期条件を用いた実験を [NIC300] と呼ぶ (表 5.1).

表 5.1: 実験設定

実験名初期条件
[GIC250]  T (x) = 250 K (0 ≤ x ≤ 1)
[PIC6_300-250] T (x) = \begin{cases}300 K (x < 0.40) \\ 250 K (x ≥ 0.40) \end{cases}
[PIC11_300-250] T (x) = \begin{cases}300 K (x < 0.72) \\ 250 K (x ≥ 0.72) \end{cases}
[NIC300]  T (x) = 300 K (0 ≤ x ≤ 1)

(a) [GIC250] と [NIC300] に用いた初期温度分布. 赤線: [NIC300], 青線: [GIC250], 黒線: 271.15 K.

(b) [PIC6_300-250] と [PIC11_300-250] に用いた初期温度分布. 緑線: [PIC6_300-250], 青線: [PIC11_300-250], 黒線: 271.15 K.

図 5.2: 水平格子点数 16 のときの初期温度分布

5.2.3 実験結果

まず, 水平格子点数 16, Q = 300 W m^{−2} として, [GIC250], [PIC6_300-250], [PIC11_300-250], [NIC300] の実験を行った.

図 5.3 に水平格子点数 16, Q = 300 W m^{−2} のときの地表面温度の時間発展を示す. 図 5.3a–5.3d を見ると, [GIC250], [PIC6_300-250], [PIC11_300-250], [NIC300] のいずれの場合においても地表面温度はある解に収束している. [GIC250] は全球凍結状態, [PIC6_300-250], [PIC11_300-250] は部分凍結状態, [NIC300] は氷なし状態の平衡状態となっており, 多重解が得られた. さらに, [PIC6_300-250], [PIC11_300-250] はどちらも部分凍結状態だが, 温度は異なっており, 結果として氷線も異なる緯度に位置する.

(a) [GIC250]  (b) [PIC6_300-250]

(c) [PIC11_300-250]  (d) [NIC300]

図 5.3: 水平格子点数 16, Q = 300 W m^{−2} のときの地表面温度の時間発展 (紫: 0 秒, 青: 5000 秒, 水色: 10000 秒, 緑: 25000 秒, 黄: 50000 秒, オレンジ: 75000 秒, 赤: 平衡時)

図5.4に水平格子点数 16, Q = 300 W m^{−2} のときの [PIC11_300-250] と [PIC6_300-250] から得られた地表面温度の緯度分布を示す. 2 つの結果を比較すると, 最も低緯度に位置する凍結した格子点は 3 格子点ずれている.

[PIC6_300-250] と [PIC11_300-250] において最も低緯度に位置する凍結した格子点が 3 格子点ずれることを詳しく調べるために, [PIC5_300-250] と [PIC11_300-250] だけでなく, 新たに 2 つの初期値を追加し (図 5.5 の左図の赤線, オレンジ線)

図 5.4: 水平格子点数 16, Q = 300 W m^{−2} のときの平衡状態における地表面温度の緯度分布. 緑色の × 印は [PIC6_300-250] の結果, 青色の × 印は [PIC11_300-250] の結果, 黒線は凝固点 271, 15 K である.

実験を行った. 図 5.5 の右図は図 5.4 に追加実験の結果を追加したものである. 図 5.5 の右図より, 最も低緯度に位置する凍結した格子点が 1 格子点ずつずれた部分凍結解が存在することが分かる.

図 5.5: 左図: [PIC6_300-250] と [PIC11_300-250] の初期温度分布と追加した初期温度分布 (詳細は表 5.2 参照). 右図: 水平格子点数 16, Q = 300 W m^{−2} のときの平衡状態における地表面温度の緯度分布

次に氷線緯度の太陽定数依存性を調べるために, 250 W m^{−2} ≤ Q ≤ 500 W m^{−2} の範囲で Q を 1 W m^{−2} きざみで変化させて実験を行った.

図 5.6 に水平格子点数 16 のときの氷線緯度と太陽定数の関係を示す. 氷線緯度は, 凝固点の前後の格子点の温度の線形補間で求めた. 全球凍結解, 部分凍結解, 氷

表 5.2: 追加した初期温度分布の詳細 ; \begin{matrix}色 & 初期条件 \\ 赤 & T (x) = \begin{cases}300 K (x < 0.59) \\ 250 K (x ≥ 0.59) \end{cases} \\ オレンジ & T (x) = \begin{cases}300 K (x < 0.66) \\ 250 K (x ≥ 0.66) \end{cases} \end{matrix}

なし解の 3 つの解や多重解が得られた. また, 初期条件によって, 2 種類の部分凍結解が存在することが確認でき, [PIC11_300-250] より [PIC6_300-250] から得られた平衡解の氷線の方が低緯度に位置している.

図 5.6: 水平格子点数 16 のときの太陽定数と氷線緯度の関係. ピンクは [GIC250] の結果, 青は [PIC6_300-250] の結果, 緑は [PIC11_300-250] の結果, 赤は [NIC300] を示す.

最後に, 氷線緯度と太陽定数の関係の水平格子点数依存性について調べるため, 水平格子点数を 50, 100, 500 として同様の実験を行った. 図 5.7 は水平格子点数 16, 50, 100, 500 における氷線緯度と太陽定数の関係である. 図 5.7 より, 部分凍結解の初期値依存性の幅は水平格子点数が増えると小さくなることがわかる.

(a) 水平格子点数 16 (b) 水平格子点数 50 (c) 水平格子点数 500

図 5.7: 太陽定数と氷線緯度の関係: (a) 水平格子点数 16, (b) 50, (c) 500. ピンクは [GIC250] の結果, 青は [PIC6_300-250] の結果, 緑は [PIC11_300-250] の結果, 赤は [NIC300] の結果を示す. (b), (c) には比較のために (a) の結果を黒点で描画している.

5.2.4 まとめと考察

本節では, GCM で見られた部分凍結解の初期値依存性が EBM でも見られるか調べた. その結果, EBM においても部分凍結解の初期値依存性が存在することや, その部分凍結解の初期値依存性の幅が水平格子点数を増やすと小さくなることが分かった.

図 5.8: 水平格子点数 16, Q = 300 W m^{−2} のとき得られた [PIC6_300-250] における平衡状態における太陽放射フラックス, 長波放射フラックス, 水平熱輸送フラックスの単位面積当たりの加熱率. 赤線は太陽放射フラックス, 青線は長波放射フラックス, 緑線は水平熱拡散フラックスの加熱率である

従って, 水平格子点数を無限に大きくすると, 部分凍結解は一意に求まるように思われる. つまり, 離散系における表現に部分凍結解の初期値依存性の原因があると想定される. そこで, (5.1) の右辺のうちどの項が特に離散化の影響を受けるのかを調べることにする. 図 5.8 に, 水平格子点数 16, Q = 300 W m^{−2} のとき得られた [PIC6_300-250] 平衡状態における太陽放射フラックス, 長波放射フラックス, 水平熱輸送の加熱率を示す. 拡散型の水平熱輸送を用いたことで温度は連続的な緯度分布をしているため, 温度の一次関数である長波放射フラックスは連続的な分布となっている. ところが, 吸収される太陽放射フラックスと正味の水平熱輸送の加熱率は氷線近傍で不連続的な分布を示している. 水平熱輸送の加熱率は他の二項と釣り合うように決定されるため, ここでは太陽放射フラックスの加熱率に着目する. (5.2) と (5.4) より, Qs(x) は連続的な分布であるが, α は不連続な分布になっている. つまり, 氷線近傍における太陽放射フラックスの加熱率が不連続な分布になるのはアルベドの分布に起因していると考えられる. そこで, 次節では離散モデルでのアルベドの表現を見直し, それによる解への影響を調べることにする.

5.3 EBM におけるアルベドの表現の影響

5.2 節より, 部分凍結解の初期値依存性は支配方程式の離散化に原因があることがわかった. また, 太陽放射フラックスに着目すると氷線近傍において不連続となっており, これはアルベドの分布に起因していると考えた. そこで本節では, 5.2 節のアルベドの表現 (図 5.10 の黒線) からサブグリッドスケールの氷の面積を考慮し連続系における表現に近づけたアルベドの表現 (図 5.10 の赤線) に変更し, 追加実験を行う. そして, その結果と 5.2 節の実験結果を比較する.

5.3.1 実験設定

本節では, サブグリッドスケールの氷の面積を考慮し, アルベド分布を考慮する. このとき, アルベド α (x_{i}) は,

α (x_{i}) = \frac{1}{Δx} ∫_{Δx} α (x) dx
= \frac{1}{Δx} ∫_{Δx} [α_{i}β (T) + α_{g} (1 − β (T))] dx  (5.9)

とする. ここで,

β (T) = \{ \begin{matrix}1 (T (x) ≤ 271.15 K) , \\ 0 (T (x) > 271.15 K) . \end{matrix}  (5.10)

T (x) = T (x_{i}) + \frac{T (x_{i+1}) − T (x_{i})}{x_{i+1} − x_{i}} (x − x_{i})  (5.11)

である. つまり, (5.11) よりアルベドは格子点の間の温度分布を線形と仮定して, 各格子内の温度が凝固点を下回る面積を推定し (図 5.9 の青い領域), 面積の重みをつけた平均によって格子点のアルベドを求める. その結果, 例えば, 水平格子点数 16, Q = 300 W m^{−2} のときの平衡状態の温度分布 (図 5.9 の赤線) におけるアルベド分布は, 図 5.10 の赤線のようになる. その他の設定は 5.2 と同じである.

図 5.9: 水平格子点数 16, Q = 300 W m^{−2} のときの平衡状態における地表面温度分布. 青で塗られた領域はある格子内の温度が凝固点を下回る領域であり, 青線で囲まれ白で塗られた領域はある格子内の温度が凝固点を上回る領域である.

図 5.10: 水平格子点数 16, Q = 300 W m^{−2} のときの平衡状態におけるのアルベド分布. 黒線はアルベドの扱いを変更する前, 赤線はアルベドの扱いを変更した後のアルベド分布である.

5.3.2 初期条件

実験名および初期条件は 5.2 と同じである (表 5.1, 図 5.2).

5.3.3 実験結果

まず, 水平格子点数 16, Q = 300 W m^{−2} として, [GIC250], [PIC6_300-250], [PIC11_300-250], [NIC300] の実験を行った.

図 5.11 に水平格子点数 16, Q = 300 W m^{−2} のときの地表面温度の時間発展を示す. 図 5.11a–5.11d を見ると, [GIC250], [PIC6_300-250], [PIC11_300-250],

[NIC300] のいずれの場合においてもある解に収束している. [GIC250] は全球凍結状態, [PIC6_300-250], [PIC11_300-250] は部分凍結状態, [NIC300] は氷なし状態の平衡状態となっており, 多重解が得られた. さらに, [PIC6_300-250], [PIC11_300-250] はどちらも部分凍結状態であり同じ温度となっており, 氷線も同じ緯度に位置する.

(a) [GIC250]  (b) [PIC6_300-250]

(c) [PIC11_300-250]  (d) [NIC300]

図 5.11: 水平格子点数 16, Q = 300 W m^{−2} のときの地表面温度の時間発展 (黒: 0 秒, 紫: 5000 秒, 青: 10000 秒, 水色: 50000 秒, 緑: 100000 秒, 黄: 150000 秒, オレンジ: 200000 秒, 赤: 平衡時)

次に氷線緯度の太陽定数依存性を調べるために, 250 W m^{−2} ≤ Q ≤ 500 W m^{−2} の範囲で Q を 1 W m^{−2} きざみで変化させて実験を行った.

図 5.12 に水平格子点数 16 のときの氷線緯度と太陽定数の関係を示す. 全球凍結解, 部分凍結解, 氷なし解の 3 つの解や多重解が得られた. また, アルベドの表現を変更したことによって, 初期条件によらず部分凍結解は一意に求められている.

図 5.12: 水平格子点数 16 のときの太陽定数と氷線緯度の関係. ピンクは [GIC250] の結果, 青は [PIC6_300-250] の結果, 緑は [PIC11_300-250] の結果, 赤は [NIC300] の結果を示す.

次に, 氷線緯度と太陽定数の関係の水平格子点数依存性について調べるために, 水平格子点数を 50, 100, 500 として同様の実験を行った. 図 5.13 は水平格子点数 16, 50, 100, 500 における氷線緯度と太陽定数の関係である. これによると, 部分凍結解は水平格子点数によらずおおよそ一意に求まることが分かる.

図 5.14 は各解像度ごとの氷線緯度と太陽定数の関係を 5.2 節のアルベドの扱いを変更する前の実験と 5.3 節のアルベドの扱いを変更した後の実験を比較したものである. 図 5.14 より, アルベドの扱い変更後の実験において得られる解はアルベドの扱い変更前の実験において水平格子点数を大きくしたときに収束する解とおおよそ一致する.

(a) 格子点数 16 (b) 格子点数 50 (c) 格子点数 500

図 5.13: 太陽定数と氷線緯度の関係: (a) 水平格子点数 16, (b) 50, (c) 500. ピンクは [GIC250] の結果, 青は [PIC6_300-250] の結果, 緑は [PIC11_300-250] の結果, 赤は [NIC300] の結果を示す. (b), (c) には比較のために (a) の結果を黒点で描画している.

\underline{Mthesis 51}

(a) 格子点数 16 (b) 格子点数 50 (c) 格子点数 500 図 5.14: 氷線緯度と太陽定数の関係の各解像度ごとの比較: (a) 水平格子点数 16, (b) 50, (c) 500. 赤はアルベドの扱いを変更する前, 青はアルベドの扱いを変更した後の 結果である.

Mthesis.tex  2018/03/23(松田幸樹)

5.3.4 まとめと考察

本節では, サブグリッドスケールの氷の面積を考慮してアルベドの表現を変更し, 部分凍結解の初期値依存性について調べた. その結果, 部分凍結解の初期値依存性は現れなかった. これは水平格子点数によらない結果である. また, アルベドの表現を変更する前の実験において水平格子点数を大きくしたときに収束した結果とアルベドの表現を変更した実験の結果はおおよそ一致した. 以上より, EBM における部分凍結解の初期値依存性はサブグリッドスケールのアルベド分布を無視したことに起因していることがわかった. そして, 本研究で用いたアルベドの表現は上手く機能しており, それを用いることで部分凍結解の初期値依存性はなくなることがわかった.

5.2 節では, アルベドの表現を変更する前の実験において, 水平格子点数を大きくすると部分凍結解の初期値依存性の幅は小さくなることがわかった. 従って, 5.2 節と同様のアルベドの表現をしている本研究の GCM においても解像度を大きくすると, 部分凍結解の初期値依存性の幅は小さくなることが期待されるが, 容易ではない. しかし, 5.3 と同様に, GCM においてもサブグリッドスケールの氷の面積を考慮したアルベド分布を与えることで第 4 章で見られた部分凍結解の初期値依存性は小さくなると考えられる.

5.4 初期値依存性の解釈

ここでは, アルベドの表現を変更する前になぜ EBM において部分凍結解の初期値依存性が現れたのかについて改めて考察することにする. それは, モデルの振る舞いの理解に有益であるためである. 本節では, (5.1) の左辺をゼロとした平衡状態におけるエネルギーバランスの式,

F_{S} − F_{OLR} + \frac{∂}{∂x} D (1 − x^{2}) \frac{∂T (x, t)}{∂x} = 0  (5.12)

を用いて, この式の解の振る舞いを考える. 任意の全球平均太陽放射 Q に対して, 以下のように (5.12) を解く:

1. 始めに氷線を仮定する. この時仮定した氷線を x_{0} とする.

2. 氷線を仮定したことによりアルベド分布が求まり, (5.12) を解くことができる.

3. (5.12) を解いて得られた温度分布から氷線を求める. この時求めた氷線を x とする.

4. このとき, x = x_{0} ならば x は (5.12) の解となる.

ここで, 物理定数の値やアルベドの与え方は 5.2 節と同様である. ただし, 上記の方法では, 小氷冠不安定や大氷冠不安定のような不安定解も解として求まることに注意しておく.

Q = 300 W m^{−2} のときの仮定した x_{0} と x の関係を図 5.15 に示す. 赤線は水平格子点数 16 のときの結果, 青線は水平格子点数 500 のときの結果である. また, 緑線は x = x_{0} の線であり, 緑線と赤線もしくは青線との交点が解となる. まず, 水平格子点数 16 のとき, x_{0} に対応する x が階段状に分布する. これは, 仮定する氷線の位置 x_{0} を変化させても, x_{0} が隣の格子点の位置に至るまではアルベドの分布に変化がなくそれに対応する x の値も変化しないためである. このため, 緑線と赤線の交点となっている部分凍結解が多数の x で実現する. 一方, 水平格子点数 500 では一定の x となる x_{0} の範囲が小さいので, 部分凍結解はおおよそ一意に求まることがわかる. 本節の実験では水平格子点数 500 を用いたが, 水平格子点数を無限にすれば連続系の表現と考えることができ, このとき x は連続な分布となると推測される. 従って, 連続系では部分凍結解は一意に求まると考えられる. また, 5.3 節のようにアルベドの表現を変更すると, x_{0} の変化に対して連続的にアルベドの分布が変化することになる. これによって, 5.3 節では部分凍結解が一意に求められたと考えられる.

GCM では, EBM には含まれていない多数の複雑な過程が取り入れられている. しかし, 本節で述べた離散化された EBM の振る舞いと同様のことは GCM においても起こっているだろう. GCM で見られた異なる氷線を持つ部分凍結解の存在には, ここで述べた機構が重要な役割を果たしていると考えられる.

図 5.15: Q = 300 W m^{−2} のときの x_{0} と x の関係 (青点は水平格子点数 16 のときのアルベドの扱い変更前の結果, 緑点は水平格子点数 500 のときのアルベドの扱い変更前の結果, 緑点は水平格子点数 16 のときのアルベドの扱い変更後の結果, 黒線は x = x_{0} の線).

第6章結論

惑星気候の多様性を考える上で, 太陽定数は重要なパラメータの一つである. 南北 1 次元エネルギーバランスモデルを用いた先行研究 (例えば, Budyko, 1969) において, 惑星気候の太陽定数依存性が調べられた. それによると, 太陽定数に応じて全球凍結解, 部分凍結解, 氷なし解の 3 つの解が現れること, および, ある太陽定数に対して複数の解が解となることが示された. また, 大気大循環や惑星放射の射出限界を陽に計算できる大気大循環モデルを用いた先行研究においても惑星気候の太陽定数依存性が調べられ, EBM で得られた 3 つの解や多重解が現れることが確認された (例えば, Ishiwatari et al., 2007). さらに, GCM を用いた研究では, 暴走温室解が存在することが確認された. 上記で挙げた 2 つの先行研究では, 海陸分布の影響を考慮しておらず, 地球とは異なる単純な大気組成を仮定した. 本研究では, 惑星大気大循環モデル DCPAM を用いて, 地球の海陸分布や大気組成を与えた実験を行い, 地球気候の太陽定数依存性を明らかにし, 大気構造の太陽定数依存性について考察した.

本研究で用いたモデルは惑星大気大循環モデル DCPAM である. 大気大循環モデルは次に述べるような様々な過程を考慮して, 惑星全球の温度, 風速, 密度分布などを計算する. DCPAM では, 力学過程はプリミティブ方程式系を解いている. また, 物理過程では, 放射過程, 乱流混合過程, 凝結過程, 雲過程, 陸面過程を考慮している. DCPAM に地球の海表面温度分布, 海表面アルベド分布を外部データとして与えた実験では NCEP の再解析データで見られた特徴を大まかに表現している. DCPAM に深さ 60 m の板海を与えた実験と DCPAM に地球の海表面温度分布, 海表面アルベド分布を外部データとして与えた実験で見られた特徴を大まかに表現している.

まず, DCPAM で 3 つの太陽定数と初期条件を与えた実験を行った. 一つ目の実験は, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1100 W m^{−2} とし, この実験を [S1100] と呼ぶ. 二つ目の実験は, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1300 W m^{−2} とし, この実験を [S1300] と呼ぶ. 三つ目の実験は, 初期条件を 280 K の等温静止乾燥大気, 太陽定数を 1475 W m^{−2} とし, この実験を [S1475] と呼ぶ. そして, この 3 つの実験から得られた年平均場を調べた. 年平均地表面温度は太陽定数が増加するにつれて大きくなる. また, [S1300] と [S1475] を比較すると, [S1475] の方が年平均東西平均地表面温度の南北差が小さくなる. これは太陽定数が増加するにつれて水蒸気量も増加し, 惑星表面において, 温室効果による加熱と潜熱フラックスによる冷却が増えるためと考えられる. 全降水量は太陽定数が増加するにつれて大きくなった. また, 太陽定数が大きくなるにつれて赤道域での降水量の最大値は大きくなった. 一方で, ハドレー循環の強さは太陽定数の増加に対して単調増加しない. 太陽定数が大きくなると, 南半球の循環の強さは弱くなり, 北半球の循環の強さは変わらなった. このことは Ishiwatari et al. (2002) とは異なる結果である. ハドレー循環の強さは太陽定数の増加に対して単調増加しないのは, 大気の安定度や地表面温度の南北差が小さいことが影響していると考えられるが, 詳細については検討が必要である.

次に, 1100 W m^{−2} から 1475 W m^{−2} までの太陽定数を与え, また複数の初期条件を与えて 11 個の実験を行い, 氷線緯度の太陽定数依存性を調べた. その結果, 全球凍結解, 部分凍結解, 氷なし解の 3 つの解と多重解が存在することが確認できた. これは, 灰色大気を持つ全球海惑星を考えた Ishiwatari et al. (2007) と同様の結果である. しかし, 本研究では海陸分布を考慮したことなどが影響し, 氷線緯度と太陽定数の関係は北半球と南半球で異なる. 例えば, 本研究では南半球では氷なし解を見つけることができていない. また, 両半球ともに暴走温室解を得ることはできなかった. これは与えた太陽定数の範囲が狭かったことによるだろう. 暴走温室解の発見は今後の課題である.

また, 北半球では緯度約 18 度より低緯度に, 南半球では緯度約 15 度より低緯度に氷線を持つ部分凍結解は得られていない. 一方, 北半球では緯度約 69 度より高緯度に氷線を持つ部分凍結解は得られていない. 以上のことから, 両半球ともに大氷冠不安定, 北半球においては小氷冠不安定が存在するように考えられる. その不安定の存在を明らかにするには, この領域に氷線をもつ部分凍結解を初期条件として与え, その解が維持されるか確認する必要があるだろう. これは今後の課題である.

太陽定数1366 W m^{−2} の場合, 初期値によって複数の部分凍結解が得られた. 280 K の等温静止乾燥大気を初期条件として与えると氷線緯度は北半球で緯度約 44 度, 南半球で緯度約 55 度に現れる. また, 太陽定数を 1250 W m^{−2} とし 280 K の等温静止乾燥大気を初期条件として 96 年計算した部分凍結解を初期条件として与えると氷線緯度は北半球で緯度約 37 度, 南半球で緯度約 35 度となった. この部分凍結解の初期値依存性は離散系であることが原因の一つであると考えらえる. より詳細に部分凍結解の初期値依存性を調べるために南北一次元エネルギーバランスモデルを用いて実験を行った. その結果, 南北一次元エネルギーバランスモデルにおいても部分凍結解の初期値依存性が見られ, 水平格子点数を大きくするとその初期値依存性の幅は小さくなることがわかった. また, サブグリッドスケールの氷面積を考慮したアルベドを用いると, 低解像度でも部分凍結解の初期値依存性が見られなくなることがわかった. 従って, 本研究で用いた GCM においてもサブグリッドスケールの氷の面積を考慮したアルベドの表現を考案し, それを用いることで, 低解像度でも部分凍結解を一意に求めることができる可能性がある. 実際に GCM で部分凍結解を一意に求めることも今後の課題である.

付録A Nakajima et al. (1992) に よる暴走温室状態の解釈の 一部

A.1 暴走温室状態

図 A.1: 対流圏界面温度 T_{tp} と 対流圏における厚さ τ_{tp} の関係 (Nakajima et al. 1992, Fig 1). 太線は (A.15), 細線は (A.17) を表す. 放射伝達方程式

\frac{2}{3} \frac{dF^{↑}(τ)}{dτ} = F^{↑}(τ) − πB(T),  (A.1)
−\frac{2}{3} \frac{dF^{↓}(τ)}{dτ} = F^{↓}(τ) − πB(T)  (A.2)

を考える. 上 2 式の和と差を取ると,

\frac{2}{3} \frac{d}{dτ} (F^{↑} + F^{↓}) = F^{↑} − F^{↓},  (A.3)
\frac{2}{3} \frac{d}{dτ} (F^{↑} − F^{↓}) = F^{↑} + F^{↓} − 2πB(T) (A.4)

を得る. 平衡状態において大気には放射エネルギーの収束がないので,

\frac{dF^{net}}{dτ} = \frac{d}{dτ} (F^{↑} − F^{↓}) = 0,  (A.5)
F^{net} = const  (A.6)

である. また, 大気上端における熱収支は,

F^{↑}(0) + 0 = F^{net} = F^{IR}  (A.7)

となる. ここで, F^{IR} は太陽放射である. (A.3) を光学的厚さ τ について大気上端 (τ = 0) から τ まで積分すると,

F^{↑} + F^{↓}−F^{↑}(0) = \frac{2}{3}F^{net}τ,  (A.8)
F^{↑} + F^{↓} = F^{IR} (\frac{3}{2}τ + 1)  (A.9)

を得る. また, (A.4) と (A.5) より,

F^{↑} + F^{↓} − 2πB(T) = 0,  (A.10)
F^{↑} + F^{↓} = 2πB(T) = 2σT^{4}  (A.11)

となる. 以上より,

σT^{4} = \frac{F^{IR}}{2} (\frac{3}{2}τ + 1)  (A.12)

を得る.

次に, 光学的厚さ τ を考える. ここでは, 簡単のため水蒸気のみが長波放射を吸収すると考えると, τ は以下のように定義される:

dτ = κ_{v}m_{v}χ_{v} \frac{dp}{m̄g} .  (A.13)

ここで, p は全圧, κ_{v} は水蒸気の吸収係数, p^{∗}(T) は水蒸気の圧力, m_{v} は水蒸気の分子量, m̄ は大気の分子量, χ_{v} ≡ p^{∗}(T)/p は飽和状態の水蒸気のモル分率である. 大気上端 (τ = 0) から τ まで積分すると,

τ = κ_{v}m_{v}χ_{v} \frac{p}{m̄g}  (A.14)

を得る.

今, 対流圏界面を考えると, (A.12) と (A.14) は,

σT_{tp}^{4} = \frac{F^{IR}}{2} (\frac{3}{2}τ_{tp} + 1) ,  (A.15)
τ_{tp} = κ_{v}m_{v}χ_{v} \frac{p_{tp}}{m̄g}  (A.16)

である. また, 対流圏界面ではモル分率は定数であり, χ_{v} ≡ p^{∗}(T_{tp})/p_{tp} とかけるので,

τ_{tp} = κ_{v}p^{∗}(T_{tp})\frac{m_{v}}{m̄g}  (A.17)

となる. (A.15) と (A.17) の解を求めるために, 2 式を τ_{tp} − T_{tp} 平面上に描画すると図 A.1 となる.

図 A.1 より, 太陽放射 F^{IR} がある値を超えると, (A.15) と (A.17) は交わらない. このように, 惑星放射の値には上限がある (Komabayashi 1967; Ingersoll 1969; Nakajima et al. 1992). この惑星放射の上限値を超える太陽放射が入射すると平衡にならない. この状態を暴走温室状態という.

付録B {全球凍結状態に関する地質 \atop 学的な証拠}

前節では, 地球気候の太陽定数依存性について, 数値計算を用いて調べた過去の研究を紹介した. その中で, 全球凍結解, 部分凍結解, 氷なし解と暴走温室解の 4 つの解が存在し得ることが示された. 部分凍結状態は現在の地球の状態である. その他 3 つの気候状態は現在の地球とは大きく異なる状態であるが, 全球凍結状態に関して地質学的な証拠からその存在が示唆されている. 本章では, 全球凍結状態の地質学的な証拠について述べる.

B.1 全球凍結状態

地質学的なデータから, 8 億年前から 6 億年前の間に全球凍結状態が現れた期間が少なくとも 2 度あったという仮説が立てられている. この仮説はスノーボールアース仮説と呼ばれており, Kirschvink (1992) や Hoffman et al. (1998) により提唱された. 地質学的なデータの中で代表的なものとして, 南極以外の世界各地での氷河堆積物の発見, 縞状鉄鉱床の存在, 炭素同位体比のマイナス値が挙げられる.

B.2 南極以外の世界各地での氷河堆積物の発見

イギリスの地質学者 Brian Harland は氷床が地球表面を全面的に覆った可能性について問題提起した. Harland は 1950 年代にスピッツベルゲンの探検を行ったところ, スピッツベルゲンの先カンブリア時代末の地層に氷河作用の痕跡を発見した. また, 広く文献を調べたところ, 北欧地域だけでなく, アフリカ, オーストラリア, 中央アジア, 南北アメリカ大陸にも氷河堆積物があることがわかった.

1980 年代, 地層が堆積した場所の緯度を推定する古地磁気学的研究から, 南オーストラリアの氷河作用でできた地層が赤道から約 10 度しか離れていない場所に堆積したものであるというデータが提示された. このデータからカリフォルニア工科大学の Kirschvink (1992) は全球凍結状態の存在を考え出した.

B.3 縞状鉄鉱床の存在

世界各地の 7 億年前の氷河堆積物中から縞状鉄鉱床が見つかった. 縞状鉄鉱床とは, 酸化鉄鉱物, 鉄炭酸塩鉱物, 硫化鉄鉱物などの鉄分にとんだ鉱物と細かい石英粒子が数センチメートルごとに規則的に繰り返された層状の地層であり, 先カンブリア時代に形成された. Kirschvink (1992) はこの縞状鉄鉱床の全球凍結状態との関係性について考察した. もし地球が全球凍結状態に陥ったならば, 海洋は大気と分離されるため酸素が乏しい状態となり, 一時的に大量の鉄イオンが蓄積する. その後, 全球凍結状態が終了すると, 大気から海洋に酸素が供給され, 海洋中の鉄イオンと酸素が反応して大量の酸化鉄が堆積し, 縞状鉄鉱床が形成されると考えた.

B.3.1 炭素同位体比のマイナス値

先カンブリア時代の地層からは化石がほとんど産出されないので, 地層の年代の比較が難しい. そこで, 地球科学者 Hoffman et al. (1998) は炭素同位体比の変動曲線を用いて, 地層の比較を行った. 炭素には質量数 12 (C12) と 13 (C13) の安定な同位体があり, その比率は物質により異なる. また, 同じ物質でも形成された時代により値のずれが生じ, このずれのパターンを用いて地層の比較を行う. Hoffman et al. (1998) はナミビア北部で地質調査を行い, その際に採取した岩石の同位体比を調べると, 先カンブリア時代後期の氷河堆積物のところで, 石灰岩の炭素同位体比がマイナスの値にずれていることがわかった. C13 と C12 を比較すると, より軽い C12 の方が植物の光合成に使われやすい. 従って, 炭素同位体比がマイナスの値にずれることは C12 の方が C13 より大きく, 植物の光合成が行われていないことを示す^{*1}. Hoffman et al. (1998) は全球凍結状態となれば, 光合成生物が生息できず光合成が行われないと考え, スノーボールアース仮説を立てるに至った.

δ^{13}C = (\frac{R}{R_{0}} − 1) × 10^{3}  (B.1)

*1炭素同位体比の計算はで表される. ここで, R =^{13} C/^{12}C であり, R_{0} は基準値である.

Mthesis.tex  2018/03/23(松田幸樹)

謝辞

本研究を行うに当たり, 多くの方にご協力を頂きました. 指導教員である神戸大学の林祥介教授には, 実験結果の考察や方向性の指導はもちろん, 研究内容や発表資料について貴重なコメントを頂きました. 神戸大学の岩山隆寛准教授には, セミナー発表において研究内容や発表資料に対する貴重なコメントを頂きました. 神戸大学の高橋芳幸准教授には, 気象学の知識, 数値モデルの使い方, 実験結果の考察等の研究内容に関することからセミナー資料や本論文の添削等非常にたくさんのことをご指導して頂きました. 北海道大学の石渡正樹准教授には, 実験結果の考察, 研究の方向性の指導, セミナー資料等の添削までテレビ会議システムやメールを通じてご指導して頂きました. 神戸大学の樫村博基特命助教には, セミナー発表において研究内容や発表資料に対する貴重なコメントをして頂き, 実験の考察等多くの助言を頂きました. 神戸大学博士後期課程の河合佑太氏には, 気象学の知識や実験結果の考察等研究内容に関して非常に多くの助言を頂きました. 北海道大学の荻原弘尭氏には, テレビ会議システムを通じて, 気象学の知識や発表資料に対して多くの助言を頂きました. 流体地球物理学教育研究分野の皆様には研究室活動でお世話になりました. この場で改めて, 皆様に感謝申し上げます.

参考文献

[1] Amante, C., Eakins, B. W., 2009: ETOPO1 1 Arc-Minute Global Relief Model: Procedures, Data Sources and Analysis., NOAA Technical Memorandum NESDIS, NGDC-24, 19 pp

[2] Arakawa, A., Suarez, M. J., 1983: Vertical differencing of the primitive equations in sigma coordinates., Mon. Wea. Rev., 111, 34-35.

[3] Beljaars, A. C. M., and Holtslag, A. A. M., 1991: Flux Parameterization over Land Surfaces for Atmospheric Models, J. Appl. Meteor., 30, 327-341.

[4] Bourke, W.P., 1988: Spectral methods in global climate and weather prediction models. Physically-Based Modelling and Simulation of Climates and Climatic Change. Part I., M.E. Schlesinger (ed.), Kluwer Academic Publishers, Dordrecht, 169-220.

[5] Budyko, M. I., 1969: The effect of solar radiation variations on the climate of the earth, Tellus, 21, 611-619.

[6] Chou, M. D., and Lee, K. T., 1996: Parameterizations for the absorption of solar radiation by water vapor and ozone, J. Atmos. Sci., 53, 1203-1208.

[7] Chou, M.-D., M. J. Suarez, C.-H. Ho, M. M.-H. Yan, and K.-T. Lee, 1998: Parameterizations for cloud overlapping and shortwave single-scattering properties for use in general circulation and cloud ensemble models, J. Climate, 11, 202-214.

[8] Chou, M. D., Suarez, M. J., Liang, X. Z., and Yan, M. M. H., 2001: A thermal infrared radiation parameterization for atmospheric studies, NASA Technical Report Series on Global Modeling and Data Assimilation, 19, NASA/TM-2001-104606.

[9] Hartmann, D. L., 1994: Global Physical Climatology, Academic Press, 411 pp.

[10] Hoffman, P. F., Kaufman, A. J., Halverson, G. P., and Schrag, D. P., 1998: A Neoproterozoic Snowball Earth, Science, 281, 1342-1346.

[11] Ishiwatari, M., Nakajima, K., Takehiro, S., and Hayashi, Y.-Y., 1998, The dependency of the strucure of the three-dimensional gray atmosphere on the solar constant and the runaway greenhouse states, Nagare Multimedia, 17, pap. 3

[12] Ishiwatari, M., Takehiro, S., Nakajima, K., and Hayashi, Y. Y., 2002: A numerical study on appearance of the runaway greenhouse state of a threedimensional gray atmosphere, J. Atmos, Sci., 59, 3223-3238

[13] Ishiwatari, M., Nakajima, K., Takehiro, S., and Hayashi, Y. Y., 2007: Dependence of climate states of gray atmosphere on solar constant: From the runaway greenhouse to the snowball states, J. Geophys., Res., 112, D13120, doi:10.1029/2006JD007368.

[14] Kiehl, J. T., and Trenberth, K. E., 1997: Earth’s annual global mean energy budget, Bull. Am. Meteor. Soc., 78, 199-208

[15] Kirschvink, J. L., 1992, Late Proterozoic low-latitude global glaciation: The Snowball Earth, Section 2.3 in: J.W. Schopf, C. Klein, D. Des Maris (eds), The Proterozoic Biosphere: A Multidisciplinary Study., Cambridge University Press, pp. 51-52.

[16] Le Treut, H., and Li, Z.-X., 1991: Sensitivity of an atmospheric general circulation model to prescribed SST changes: feedback effects associated with the simulation of cloud optical properties, Clim. Dyn., 5, 175-187.

[17] Lindzen, R. S., 1990: Dynamics in Atmospheric Physics, CAMBRIDGE, 310 pp.

[18] Manabe, S., 1969: Climate and the ocean circulation I. The atmospheric circulation and the hydrology of the Earth ’ s surface, Mon. Wea. Rev., 97, 739-774.

[19] Matthews, E., 1983: Global vegetation and land cover: New high-resolution data bases for climate studies, J. Clim. Appl. Meteor., 22, 474-487.

[20] Matthews, E., 1984: Prescription of Land-surface Boundary Conditions in GISS GCM II: A Simple Method Based on High-resolution Vegetation Data Sets. NASA TM-86096. National Aeronautics and Space Administration. Washington, D.C.

[21] Mellor, G. L., and Yamada, T., 1982: Development of a Turbulent Closure Model for Geophysical Fluid Problems, Rev. Geophys. Space Phys., 20, 851-875.

[22] Moorthi, S., and Suarez, M. J., 1992: Relaxed Arakawa-Schubert: A parameterization of moist convection for general circulation models, Mon. Wea. Rev., 120, 978-1002.

[23] Nakajima, S., Hayashi, Y. Y., and Abe, Y., 1992: A Study on the “Runaway Greenhouse Effect” with a One-Dimensional Radiative-Convective Equilibrium Model, J. Atmos. Sci., 49, 2256-2266.

[24] North, G. R., 1975: Analytical Solution to a Simple Climate Model with Diffusive Heat Transport, J. Atmos. Sci., 32, 1301-1307.

[25] North, G. R., Cahalan, R. F., and Coakley, J. A., 1981: Energy Balance Climate Models, Geophys. Res. Lett., 39, 91-121.

[26] Sellers, W. D., 1969: A climate model based on the energy balance of the earth-atmosphere system, J. Appl. Meteor., 8, 392-400.

[27] Takahashi, Y. O., M., Ishiwatari, S., Noda, M., Odaka, T., Horinouchi, Y., Morikawa, Y.-Y., Hayashi, DCPAM Devlopment Group, 2010: DCPAM: planetary atmosphere model, http://www.gfddennou.org/library/dcpam/, GFD Dennou Club

[28] Taylor, K. E., Williamson, D., and Zwiers, F., 2000: The sea surface temperature and sea ice concentration boundary conditions for AMIP II simulations, PCMDI Report 60, 29 pp

[29] Taylor, K. E., Stouffer, R. J., and Meehl, G. A., 2012: An Overview of CMIP5 and the Experiment Design, Bull. Amer. Meteoro. Soc., 93, 485-498.

[30] Toon, O. B., McKay, C. P., and Ackerman, T. P., 1989: Rapid calculation of radiative heating rates and photodissociation rates in inhomogeneous multiple scattering atmospheres, J. Geophys. Res., 94, 16287-16301.

[31] Vonder, Haar, T. H., and Suomi, V. E., 1971: Measurements of Earth’s Radiation Budget from Satellites During a Five-Year Period. Part I: Extended Time and Space Means, J. Atmos. Sci., 28, 305-314.

[32] 会田勝, 1954: 気象学のプロムナード 8 大気と放射過程 - 大気の熱源と放射収支を探る -, 東京堂出版, 280 pp

[33] 浅野正二, 2010: 大気放射学の基礎, 朝倉書店 267 pp

[34] 浅野正二, 2005a: 大気の放射過程. 気象ハンドブック第 3 版, 朝倉書店, 18-29 pp

[35] 井谷優花, 2012: 大気大循環モデルによる地球気候の表現とその自転軸傾斜角依存性, 神戸大学理学研究科修士論文

[36] 小倉義光, 1984: 一般気象学, 東京大学出版, 301 pp [37] 川上紳一, 東條文治, 2009: 図解入門地球史がよくわかる本, 秀和システム, 331

pp

[38] 岸保勘三郎, 田中正之, 時岡達志, 1982: 大気科学講座 4 大気の大循環, 東京大学出版会, 260 pp

[39] 田中浩, 2004: 大気科学とその周辺, 現代図書, 641 pp

[40] 地球流体電脳倶楽部, http://www.gfd-dennou.org/

[41] 地球流体電脳倶楽部, DCPAM5 支配方程式系の導出に関する参考

資料, http://www.gfd-dennou.org/library/dcpam/dcpam5/dcpam5_
latest/doc/derivation/pub/derivation.pdf

[42] 地球流体電脳倶楽部, DCPAM5 支配方程式とその離散化, http:

//www.gfd-dennou.org/library/dcpam/dcpam5/dcpam5_latest/
doc/basic_equations/pub/basic_equations.pdf

[43] 地球流体電脳倶楽部, DCPAM5 ドキュメント, http://www.gfd-dennou.org/

library/dcpam/dcpam5/dcpam5_latest/doc/

[44] 地球流体電脳倶楽部, 時間ステップ依存性雲の消滅時定数まとめ,

https://www.gfd-dennou.org/library/dcpam/sample/2016-10-21_
yot/Earth/cloud_time_scale.pdf