
%% % 卒業論文 
%% % タイトル「未定」
%% %
%% % 20XX/XX/XX 修正
%% % 20XX/XX/XX 作成
%% %
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%             Style  Setting             %%%%%%%%
%% % フォント: 12point (最大), 片面印刷
%% \documentclass[a4j,12pt,openbib,oneside]{jreport}

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%             Package Include            %%%%%%%%
%% \usepackage{ascmac}
%% \usepackage{tabularx}
%% \usepackage[dvipdfmx]{graphicx}
%% %\usepackage[dvips]{graphicx}
%% %\usepackage{float}
%% \usepackage{amssymb}
%% \usepackage{amsmath}
%% \usepackage{Dennou6}		% 電脳スタイル ver 6
%% \usepackage{bm}
%% %\usepackage{fancybox}

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%            PageStyle Setting           %%%%%%%%
%% \pagestyle{DAmyheadings}


%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%        Title and Auther Setting        %%%%%%%%
%% %%
%% %%  [ ] はヘッダに書き出される.
%% %%  { } は表題 (\maketitle) に書き出される.

%% \Dtitle{sotsuron}         % 変更不可
%% \Dauthor{惑星太郎}        % ゼミ担当者の名前
%% \Ddate{2019/02/09}        % ゼミの日時 (毎回変更すること)
%% \Dfile{sotsuron.tex}

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%   Set Counter (chapter, section etc. ) %%%%%%%%
%% \setcounter{chapter}{1}    % 章番号
%% \setcounter{section}{0}    % 節番号
%% \setcounter{subsection}{1}
%% \setcounter{equation}{0}   % 式番号
%% \setcounter{page}{1}     % 必ず開始ページは明記する
%% \setcounter{figure}{1}     % 図番号
%% \setcounter{table}{0}      % 表番号
%% \setcounter{footnote}{0}


%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%        Counter Output Format           %%%%%%%%
%% \def\thechapter{\arabic{chapter}}
%% \def\thesection{\arabic{chapter}.\arabic{section}}
%% \def\thesubsection{\arabic{chapter}.\arabic{section}.\arabic{subsection}}
%% \def\theequation{\arabic{chapter}.\arabic{equation}}
%% \def\thepage{\arabic{page}}
%% \def\thefigure{\arabic{chapter}.\arabic{figure}}
%% \def\thetable{\arabic{chapter}.\arabic{table}}
%% \def\thefootnote{*\arabic{footnote}}


%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%        Dennou-Style Definition         %%%%%%%%

%% %% 改段落時の空行設定
%% \Dparskip      % 改段落時に一行空行を入れる
%% %\Dnoparskip    % 改段落時に一行空行を入れない

%% %% 改段落時のインデント設定
%% \Dparindent    % 改段落時にインデントする
%% %\Dnoparindent  % 改段落時にインデントしない

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%  macro  
%% \newcommand{\pdiff}[3]{
%%   \ifnum #1<1%
%%   \text{エラー}(最初の引数は1以上)\else{%
%%     \ifnum #1=1 \frac{\partial #2}{\partial #3}\else%
%%       \frac{\partial^{#1} #2}{\partial #3^{#1}}\fi
%%   }\fi%
%% }

%% \newcommand*{\aboxed}[2]{%
%%   \rlap{\boxed{#1#2}}%
%%   \phantom{\hskip\fboxrule\hskip\fboxsep #1}&\phantom{#2}%
%% }

%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%             Textcd  Start                 %%%%%%%%
%% \begin{document}

\setcounter{chapter}{1}

\chapter{モデルの概要と実験設定}
\markright{第\arabic{chapter}章 モデルの概要と実験設定} %  節の題名を書き込むこと
% ここから本文を書く
本章では, 本研究で使用した南北 1 次元エネルギーバランスモデル (EBM) についての説明を行う. 今回行う実験では, 格子点数や初期値を変化させていくつかの実験を行い, 全球平均太陽放射フラックスに対しての氷線緯度の依存性を計算する. このモデルは, Sellers (1969) \cite{Sellers} の定式化に従う.
\section{南北 1 次元エネルギーバランスモデル}

\subsection{支配方程式系}

このモデルの支配方程式は, 以下の四つの項からなる方程式である.
\begin{equation}
  C\frac{\partial T(x,t)}{\partial t}=F_{s}-F_{OLR}+\frac{\partial}{\partial x}\bigg{\{}D(1-x^2)\frac{\partial T(x,t)}{\partial x}\bigg{\}}
\end{equation}
ここで, $x$ は緯度の正弦, $t~[\textrm{s}]$ は時間, $T~[\textrm{K}]$ は地表面温度, $C$ は熱容量, $D$ は拡散係数である. 本研究では, Lindzen (1990) \cite{Lindzen} の値を参照し, $C=1.0~[\textrm{J}/\textrm{K}/\textrm{m}^2]$, $D=0.2~[\textrm{W}/\textrm{m}^2/\textrm{K}]$ とした. また, $F_{s}$ は太陽放射フラックス, $F_{OLR}$ は惑星放射フラックス, そして, 式中の右辺第三項は, 水平熱輸送によって吸収されるフラックスである. さらに, 太陽放射フラックス ($F_{s}$), 惑星放射フラックス ($F_{OLR}$) は, 以下のように書ける.
\begin{align}
  F_{s}&=Qs(x)(1-\alpha(x))\\
  F_{OLR}&=A+BT(x,t)
\end{align}
ここで, $Q$ は全球平均太陽放射フラックス, $A$, $B$ は定数, $\alpha$ はアルベドである. また, これらの値に関しても, Lindzen (1990) \cite{Lindzen} の値を参照し, $A=212.15~[\textrm{W}/\textrm{m}^{2}]$, $B=1.55~[\textrm{W}/\textrm{m}^{2}/\textrm{K}$] とした. $s(x)$ は,年平均太陽放射フラックスの子午面分布であり, 二次のルジャンドル近似を用いることで, 
\begin{equation}
  s(x)=1+s_{2}p_{2}(x)=1+s_{2}\frac{(3x^2-1)}{2}
\end{equation}
と表される. ここで, $s_{2}=-0.482$ であり, これは, Lindzen (1990) \cite{Lindzen} の値を参照したものである. また, $s(x)$ は, 全球に渡って面積積分を行うことで球の表面積となるように規格化され, 今回の計算領域においては, 1 となるように規格化されている.\footnote{$s(x)$ については付録 A を見よ.}

\subsection{境界条件}

このモデルでは, 物理量を東西平均し, 南北対称としている. 計算領域は赤道 ($x=0$) から, 極 ($x=1$) である. また, 計算領域の両端においては熱の出入りはゼロとするため, 境界条件は,
\begin{align}
  \bigg{(}\frac{\partial T}{\partial x}\bigg{)}=0~~(x=0),\\
  \bigg{(}\frac{\partial T}{\partial x}\bigg{)}=0~~(x=1),
\end{align}
である.

\subsection{離散化の方法}

ここでは, 離散化の方法について説明する. 本研究で行う実験での計算領域は, 赤道から極までである. その領域を $x=0$ を赤道, $x=1$ を極として, $N$ 個の領域に分割し, それぞれの領域における物理量の平均的な値をその中心に置く. そして, その中心に置いた点を計算上の格子点とし, 方程式を空間に対して離散化する.

%離散化の図の挿入
\begin{figure}[htbp]
  \begin{center}
    \includegraphics[width=12cm]{./figs/risannka_kousiten.png}
  \caption{離散モデルにおける変数配置}
  \end{center}
\end{figure}

従って, (2.1) は, 下に示す (2.7) のように書ける.
\begin{equation}
\begin{split}
  C\frac{\partial T_{i}(t)}{\partial t}&=Qs(x_{i})(1-\alpha(x_{i}))-(A+BT_{i}(t))\\
  &+\frac{D}{\Delta x^2}[(1-x^{2}_{i+\frac{1}{2}})(T_{i+1}(t)-T_{i}(t))-(1-x^{2}_{i-\frac{1}{2}})(T_{i}(t)-T_{i-1}(t))]
\end{split}
\end{equation}
ここで, $i$ は格子点番号, $\Delta x$ は格子点間隔である. 

\section{数値計算の解法}

上記で説明したエネルギーバランスモデルを用いて, 氷線緯度の太陽定数依存性を調べるために, 二種類の方法を用いる. 一つは, 地表面温度の時間発展を計算することで定常状態を求め, 安定解を求める方法である. (2.7) を時間に対しても離散化すると, 
\begin{align}
  C\frac{T^{l+1}_{i}-T^{l}_{i}}{\Delta t}&=Qs(x_{i})(1-\alpha(x_{i}))-(A+BT^{l}_{i})\notag\\
  &+\frac{D}{\Delta x^2}[(1-x^{2}_{i+\frac{1}{2}})(T^{l}_{i+1}(t)-T^{l}_{i}(t))-(1-x^{2}_{i-\frac{1}{2}})(T^{l}_{i}(t)-T^{l}_{i-1}(t))]
\end{align}
となる. この式において右辺が, $10^{-5}~[\textrm{W}/\textrm{m}^2]$ よりも小さくなるとき収束したと判断する. もう一つは, ある格子点に水の凝固点温度, $T_{k}=271.15~[\textrm{K}]$, を与え, $T_{1},\cdots,T_{k-1},T_{k+1},\cdots,T_{N},Q$ を未知数とする連立一次方程式を作成し, その連立方程式を行列を用いて計算する. 以下にその方法を具体的に示す. 空間に対して離散化した方程式 (2.7) をまず以下のように式変形し, $i$ を 1 から $N$ まで変化させた際の連立一次方程式を作成する.
\begin{equation}
  \begin{split} 
    -\frac{D}{(\Delta x)^2}(1-x^2_{i-\frac{1}{2}})T_{i-1}+\bigg{\{}B&+\frac{D}{(\Delta x)^2}(1-x^2_{i+\frac{1}{2}})+\frac{D}{(\Delta x)^2}(1-x^2_{i-\frac{1}{2}})\bigg{\}}T_{i}\\
    &-\frac{D}{(\Delta x)^2}(1-x^2_{i+\frac{1}{2}})T_{i+1}-s(x_{i})(1-\alpha(x_{i}))Q=-A
  \end{split}
\end{equation}
次に, $T_{1},\cdots,T_{k-1},T_{k+1},\cdots,T_{N},Q$ を未知数として, $A\vec{x}=\vec{b}$ の形の行列に変形する.
\begin{align}
  \begin{pmatrix}
      M'_{1,0}+M'_{1,1} & M'_{1,2} & 0 & \dots & 0 & 0 & Q'_{1}\\
      M'_{2,1} & M'_{2,2} & M'_{2,3} & \dots & 0 & 0 & Q'_{2}\\
      0 & M'_{3,2} & M'_{3,3} & \dots & 0 & 0 & Q'_{3}\\
      \vdots &\vdots & \vdots & \ddots & \vdots & \vdots & \vdots\\
      0 & 0 & 0 & \dots & M'_{N-1,N-1} & M'_{N-1,N} & Q'_{N-1}\\
      0 & 0 & 0 & \dots & M'_{N,N-1} & M'_{N,N}+M'_{N,N+1} & Q'_{N}\\
      \delta_{1} & \delta_{2} & \delta_{3} & \dots & \delta_{N-1} & \delta_{N} & 0
  \end{pmatrix}
  \begin{pmatrix}
      T_{1}\\
      T_{2}\\
      T_{3}\\
      \vdots\\
      T_{N-1}\\
      T_{N}\\
      Q
  \end{pmatrix}
  =
  \begin{pmatrix}
      -A\\
      -A\\
      -A\\
      \vdots\\
      -A\\
      -A\\
      271.15
  \end{pmatrix}
\end{align}
ここで, $T_{i}$ は, $i$ 番目の格子点の温度, $M'_{i,j}$ は, $i$ 番目の格子点の $T_{j}$ の係数, $Q'_{i}$ は, $i$ 番目の格子点の $Q$ の係数である. また, $\delta_{i}$ は, $i=k$ のとき, $\delta_{i}=1$ をとり, $i\neq k$ のとき, $\delta_{i}=0$ をとる. この行列を計算することで, 安定解である部分凍結解を含む不安定解を求める.

\section{アルベドの設定}

本研究では, アルベドをサブグリッドスケールでの氷の面積を考慮して設定する. 以下に, それぞれの緯度帯に置ける温度を表す図を示す (図 2.2). この図を見ると, 縦向きの赤い線分と横向きの水の凝固点温度となる赤い直線が交わる点よりも左側では, 地表面温度は水の凝固点温度よりも大きく, また, 右側では, 小さくなっている. 従って, 凝固点温度である緯度よりも極側では, 氷床が発達し, アルベドが大きくなり, 赤道側では, 氷床が存在せず, アルベドが低くなることが分かる. そこで, 格子点間の温度分布を線形と仮定し, 格子内の凝固点温度を下回る面積と上回る面積を推定して, 格子領域の面積に対するそれぞれの面積の比を重みとして, 平均をとることで格子内のアルベドを求める. ここで, 氷床が存在するときのアルベドを 0.6, 氷床が存在しないときのアルベドを 0.1 とすると, アルベドは,  
\begin{equation}
  \begin{split}
    \alpha(x_{i})&=\frac{1}{\Delta x}\int_{\Delta x} \alpha(x)~dx\\
    &=\frac{1}{\Delta x}\int_{\Delta x}[0.6\beta(T)+0.1(1-\beta(T))]~dx
  \end{split}
\end{equation}
である. ここで, 
\begin{align}
  \beta(T)&=
  \begin{cases}
    1~(T(x)\leq 271.15~\textrm{K})\\
    0~(T(x)> 271.15~\textrm{K})\\
  \end{cases}\\
  T(x)&=T(x_{i})+\frac{T(x_{i+1})-T(x_{i})}{x_{i+1}-x_{i}}(x-x_{i})
\end{align}
である. これにより, 各格子内のアルベドを算出する. 
\begin{figure}[t]
  \begin{center}
    \includegraphics[width=7cm]{./figs/albedodesc_wl.png}
    \caption{地表面温度の緯度分布が水の凝固点温度を横切るときのアルベドの面積比を表した図. 横軸: 緯度, 縦軸: 地表面温度. 横向きの赤線は, 水の凝固点温度, 縦の青い線は格子領域を表し, 地表面温度が水の凝固点温度を横切る格子を濃い青線で囲んでいる. 赤の矢印線は, アルベドが $0.1$, 青の矢印線は, $0.6$ の領域を表す.}
  \end{center}
\end{figure}

\section{実験設定}

今回の実験では, 四種類の実験を行う. 実験の各名称は次のとおりである. [GIC250], [NIC300], [PIC6\_300-250], [PIC11\_300-250]. 各実験では, 初期値を変化させて実験を行う. [GIC250] では, $0\leq x\leq1$ の範囲で, 初期値を $T(x)=250~\textrm{K}$ とし, [NIC300] では, $0\leq x\leq1$ の範囲で, 初期値を $T(x)=300~\textrm{K}$ としている. また, [PIC6\_300-250] では, 初期値を $T(x)=250~(x\geq0.4)$, $T=300~(x<0.4)$ とし, [PIC11\_300-250] では, 初期値を $T(x)=250~(x\geq0.72)$, $T=300~(x<0.72)$ とし, 緯度によって初期値を変化せている. (表 2.1)

今回の実験では, 標準実験として格子点数を $N=16$ として各計算を行い, 全球平均太陽放射フラックスは, $250~[\textrm{W}/\textrm{m}^2]$ - $500~[\textrm{W}/\textrm{m}^2]$ で変化させる (表 2.2).
\begin{table}[htbp]
  \caption{初期値の設定}
\begin{center}
\begin{tabular}{|c|c|}\hline
  実験名 & 実験設定\\\hline\hline
  [GIC250] & $T(x)=250~\textrm{K}~(0\leq x\leq 1)$\\\hline
  [NIC300] & $T(x)=300~\textrm{K}~(0\leq x\leq 1)$\\\hline
  [PIC6\_300-250] & 
  $\begin{array}{c}
    T(x)=
    \begin{cases}
    300~\textrm{K}~(x<0.4) \\
    250~\textrm{K}~(x\geq 0.4)
    \end{cases}
  \end{array}$ \\\hline
  [PIC11\_300-250] & 
  $\begin{array}{c}
    T(x)=
    \begin{cases}
    300~\textrm{K}~(x<0.72) \\
    250~\textrm{K}~(x\geq 0.72)
    \end{cases}
  \end{array}$ \\\hline
\end{tabular}
\end{center}
\end{table}


\begin{figure}[htbp] % H: here（ここに配置）
  \centering
  \begin{minipage}{0.45\linewidth}
      \centering
      \includegraphics[width=\linewidth]{./figs/first_number_300250.png}
      \begin{flushleft}
        {(a) [GIC250] と [NIC300] の初期温度分布. 青線が [GIC250], 紫線が [NIC300] を表す. 黒線は, $T=271.15~[K]$}
      \end{flushleft}
  \end{minipage}
  \hspace{0.05\linewidth} % 画像の間にスペースを追加
  \begin{minipage}{0.45\linewidth}
      \centering
      \includegraphics[width=\linewidth]{./figs/first_number_04072.png}
      \begin{flushleft}
      {(b) [PIC6\_300-250] と [PIC11\_300-250] の初期温度分布. 赤線が, [PIC6\_300-250], 緑線が, [PIC11\_300-250] を表す. 黒線は, $T=271.15~[K]$}
      \end{flushleft}
  \end{minipage}
  \caption{格子点数が $N=16$ のときの, 各実験における初期温度分布}
  \label{fig:two_images}
\end{figure}



% ここで本文おしまい
%\end{document}                  
%%%%%%%%              Text End                  %%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%              Sample                    %%%%%%%%

%\begin{figure}[h]
%  \begin{center}
%    \includegraphics[width=10cm]{fig6_2.PNG}
%    \caption{\footnotesize{}}
%  \end{center}
%\end{figure}

%extractbb ***.PNG
