
%% % 卒業論文 
%% % タイトル「金星大気大循環モデルにおける静的安定度の 傾圧不安定への影響」
%% %
%% % 2023/03/23 修正
%% % 2022/12/20 作成
%% %
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%%%%%%%             Textcd  Start                 %%%%%%%%
%% \begin{document}

\setcounter{chapter}{1}

\chapter{数値モデルの概要と実験設定}
\markright{第\arabic{chapter}章 モデルの概要と実験設定} %  節の題名を書き込むこと

本研究では, 金星大気大循環モデル AFES-Venus を用いて実験を行った. 
AFES-Venus は地球シミュレータ用大気大循環モデル
 (AFES: Atmospheric general circulation model For the Earth Simulator)\cite{Ohfuchi}\cite{Enomoto} 
 を金星版に改変\cite{Sugimotoetal2014}したものである. 
以下では, AFES-Venus の概要や実験設定について記す. 
なお, 安定度に関係する実験設定は各実験によって異なるため, 3-5 章の各章に記す. 

\section{座標系と離散化}
AFES-Venus では, 水平方向に球面座標系 $(\lambda$, $\phi)=(\rm{経度}, \rm{緯度})$ を, 
鉛直方向に $\sigma$ $(\sigma\equiv p/p_s)$ 座標系を用いる. 
ここで, $p$ は気圧, $p_s$ は地表面気圧である. 
経度方向の格子点は各緯度で等間隔であるが, 緯度は Gauss 緯度を用いており, 極ほど格子点間隔は小さい. 
水平方向の離散化にはスペクトル法を, 鉛直方向の離散化には差分法を用いる. 

空間解像度は T63L120 (水平方向 192 $\times$ 96 格子, 鉛直 120 層) である. 
ここで, T は三角切断を意味しており, 63 は切断波数を表している. 
鉛直方向には, 地表から高度 120 km までを考えており, 各層の間隔は 1 km である. 
また, 出力間隔は 1 地球日で, 準定常状態に達するまで 5 年間時間積分を行った. 
5 年のうち最後の 1 年を解析に用いた. 

\section{力学過程}
AFES-Venus は全球三次元プリミティブ方程式系を支配方程式とする力学モデルで, 
雲過程や光化学反応, 地形を省略している. 
この方程式系は $\sigma$ 座標系において以下のような式で構成される. 
なお, 以下で用いられる変数の定義は表 \ref{tab:val} の通りである. \\
\noindent
\textbf{運動方程式}
\begin{align}
  \frac{\partial u}{\partial t} &= (\zeta+f)v -\frac{1}{a\cos\phi}\frac{\partial}{\partial\lambda}(\Phi+E)-\frac{R_dT_v}{a\cos\phi}\frac{\partial}{\partial\lambda}\ln p_s -\dot{\sigma} \frac{\partial u}{\partial\sigma} - \mathcal{D}(u)\\
  \frac{\partial v}{\partial t} &= -(\zeta+f)u-\frac{1}{a}\frac{\partial}{\partial\phi}(\Phi+E)-\frac{R_dT_v}{a}\frac{\partial}{\partial \phi}\ln p_s -\dot{\sigma}\frac{\partial v}{\partial \sigma} -\mathcal{D}(v)
\end{align}

\noindent
\textbf{連続の式}
\begin{align}
  \frac{\partial\ln p_s}{\partial t} = -(\bm{v}_h\cdot\nabla_{\sigma}\ln p_s+\delta)-\frac{\partial\dot{\sigma}}{\partial\sigma}
\end{align}

\noindent
\textbf{静力学方程式}
\begin{align}
  \frac{\partial\Phi}{\partial(\ln\sigma)} = -R_dT_v
\end{align}

\noindent
\textbf{熱力学方程式}
\begin{align}
  \frac{\partial T'}{\partial t} = -\frac{1}{a(1-\mu^2)}\frac{\partial}{\partial\lambda}(UT')-\frac{1}{a}\frac{\partial}{\partial\mu}(VT')+T'\delta-\dot{\sigma}\frac{\partial T}{\partial\sigma} + Q
\end{align}
ここで, 
\begin{align}
  \nabla_{\sigma} &= \mathbf{i} \frac{1}{a\cos\phi}\left(\frac{\partial}{\partial\lambda}\right)_{\sigma} + \mathbf{j} \frac{1}{a}\left(\frac{\partial}{\partial\phi}\right)_{\sigma}\\
  \zeta &= \frac{1}{a\cos\phi}\left[\frac{\partial v}{\partial \lambda}-\frac{\partial}{\partial\phi}(u\cos\phi)\right] \\
  f &= 2\Omega\sin\phi \\
  E &= \frac{u^2+v^2}{2} \\
  \delta &= \frac{1}{a\cos\phi}\left[\frac{\partial u}{\partial \lambda}-\frac{\partial}{\partial\phi}(v\cos\phi)\right] \\
  \bm{v}_h &= (u,v)\\
  T(\lambda,\mu, \sigma, t) &=  T_0(\sigma) + T'(\lambda, \mu, \sigma,t)
\end{align}
であり, 下付き添え字の $\sigma$ は $\sigma$ 座標における微分の意味である.  
\begin{table}[htpb]\caption[short]{変数の定義}
  \begin{center}
 \begin{tabular}{@{\vrule width 1pt\ }c|c||c|c@{\ \vrule width 1pt}}
  \noalign{\hrule height 1pt}
   文字 & 意味 & 文字 & 意味 \\
  \noalign{\hrule height 1pt}
  $t$ & 時間 & $a$ & 金星半径 \\
  \noalign{\hrule height 1pt}
  $\lambda$ & 経度 & $\mathcal{D}$ & 水平拡散項 \\
  \noalign{\hrule height 1pt}
  $\phi$ & 緯度 & $p_s$ & 地表面気圧 \\
  \noalign{\hrule height 1pt}
  $\sigma$ & $\sigma$ 高度 & $R_d$ & 乾燥空気の気体定数 \\
  \noalign{\hrule height 1pt}
  $\dot{\sigma}$ & 鉛直 $\sigma$ 速度 & $\mu$ & $\sin\phi$ \\
  \noalign{\hrule height 1pt}
  $u$ & 東西風速 & $U$ & $u\cos\phi$ \\
  \noalign{\hrule height 1pt}
  $v$ & 南北風速 & $V$ & $v\cos\phi$ \\
  \noalign{\hrule height 1pt}
  $T$ & 温度 & $Q$ & 加熱 \\
  \noalign{\hrule height 1pt}
  $T_v$ & 仮温度 & $\zeta$ & 相対渦度の鉛直成分 \\
  \noalign{\hrule height 1pt}
  $T_0$ & 標準温度 & $f$ & コリオリパラメータ \\
  \noalign{\hrule height 1pt}
  $T'$ & 温度の変動成分 & $E$ & 運動エネルギー \\  
  \noalign{\hrule height 1pt}
  $\Phi$ & ジオポテンシャル & $\delta$ & 水平発散 \\
  \noalign{\hrule height 1pt}
  $\Omega$ & 自転角速度 & $\bm{v}_h$ & 水平風 \\
  \noalign{\hrule height 0.4pt}
 \end{tabular}
 \end{center}
 \label{tab:val}
\end{table}
また, 鉛直方向の境界条件は, 
\begin{align}
  \dot{\sigma} = 0 \qquad \rm{for}\qquad \sigma=0,1 
\end{align}
である. 



\section{物理過程}

放射過程は太陽加熱関数\cite{Tomasko}\cite{Crisp1986} (太陽加熱の日変化あり) と, ニュートン冷却\cite{Crisp}
\begin{equation}
  \frac{\mathrm{d}T}{\mathrm{d}t}= \frac{T_{nc}-T}{\tau}
\end{equation}
によって簡略化されている. 
ここで, $T$ は温度, $\tau$ はニュートン冷却の緩和時定数, $T_{nc}$ はニュートン冷却の基準温度である. 
太陽加熱は高度 0-80 km の範囲に与えており, 60-80 kmで大きい (図 \ref{fig:QSOL}). 
ニュートン冷却は, 基準温度場との差に応じた冷却を与えることで, 温度場を基準温度場に近づけるはたらきがある. 
ニュートン冷却の緩和時定数の鉛直分布は図 \ref{fig:relaxation_time_of_Newtoncooling} であり, 
これは Sugimoto et al. (2014)\cite{Sugimotoetal2014} で用いたものと同じである. 
ニュートン冷却の基準温度場は東西方向に一様で, 南北方向に勾配がある. 
また、大気上端における波の反射を防ぐためのスポンジ層を高度 80 km 以上に設定している. 
地表面摩擦はレイリー摩擦によって表現している. 
渦粘性係数は 0.15 $\rm{m^2/s}$ である. 
さらに, 大気が静的不安定になった際に中立に戻すために乾燥対流調節が導入されている. 
\begin{figure}[htpb]
  \begin{tabular}{cc}
    \begin{minipage}[b]{0.5\linewidth}
      \centering
        \includegraphics[width=5.5cm]{figs/result1/case1/QSOLzh.png}\subcaption{}
        \label{fig:QSOL}
    \end{minipage}
    \begin{minipage}[b]{0.5\linewidth}
      \centering
        \includegraphics[width=6cm]{figs/relaxationtime.jpg}\subcaption{}
        \label{fig:relaxation_time_of_Newtoncooling}
    \end{minipage}

  \end{tabular}
  \caption{(a) 実験 1 の標準実験 (STD0) の, 東西時間平均した, 赤道における太陽加熱の高度分布. 横軸の単位は [K/s]. 平均期間は 1 年である. 各実験の太陽加熱の分布は付録 I の図集を見よ.    
   (b) ニュートン冷却の緩和時定数の鉛直分布. (Sugimoto et al. (2014) Figure 1a. より引用\cite{Sugimotoetal2014}. )}
\label{fig:heating}
\end{figure} 

\section{初期条件}
実験 1 では, 初期の静的安定度の鉛直分布として図 \ref{fig:init_SS1} を用いた. 
雲層 (55-60 km) では, ほぼ中立 (0.1 K/km) であり, 雲層の下端 (45 km) は安定で 4.0 K/km である. 
また, 下層 (0-35 km) は弱く成層しており, 0.5 K/km という値になっている. 
%雲層より上 (70 km 以上) は非常に安定で ** K/km である.  
この安定度分布は観測結果を模した分布であり, 
Kashimura et al. (2019)\cite{Kashimuraetal2019}, 
Sugimoto et al. (2014)\cite{Sugimotoetal2014} で用いられた安定度分布と同じである. 
実験 2, 3 では, 下層の強制安定度を 4.0 K/km という大きな値に設定したが, 図 \ref{fig:relaxation_time_of_Newtoncooling} から分かるように, 下層はニュートン冷却の緩和時定数が大きく, 
実験 1 と同じ初期安定度分布を使うと, 下層で実現される安定度が大きくならなかった. 
そのため, 実験 2, 3 では, 下層 (0-35 km) の安定度を 4.0 K/km に変更した, 図 \ref{fig:init_SS2} を初期の安定度分布として用いた. 

初期東西風速場には理想化された剛体回転のスーパーローテーションを与えた. 
図 \ref{fig:init_U1}, \ref{fig:init_U2} はそれぞれ実験 1, 実験 2,3 で用いた東西方向に一様な初期東西風の子午面分布である. 
この風速分布は
\begin{equation}
  u(\phi, z)=\left\{\begin{array}{ll} 
   u_{max} (\cos\phi) z/z_{sc}  &\text{for} z<z_{sc} \\   
  u_{max}\cos\phi & \text{for} z\ge z_{sc}, 
\end{array}\right. 
\end{equation}
で表される. 
ここで, $u_{max}$ は東西風速の最大値, $z$ は高度, $z_{sc}$ は赤道における $u_{max}$ の高度であり, 
$u_{max}=100$ m/s, $z_{sc}=70$ km とした. 
\begin{figure}[htpb]
  \begin{tabular}{cc}
\begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=7.5cm]{figs/init/initstab1.png}\subcaption{実験 1 で用いた初期安定度分布}
    \label{fig:init_SS1}
\end{minipage}
\begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=7.5cm]{figs/init/initstab2.png}\subcaption{実験 2, 3 で用いた初期安定度分布}
    \label{fig:init_SS2}
\end{minipage}
\end{tabular}
\begin{tabular}{cc}
  \begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=6cm]{figs/init/Uzh.png}\subcaption{実験 1 で用いた初期東西風速場}
    \label{fig:init_U1}
\end{minipage}
\begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=6cm]{figs/result2/init/Uzh.png}\subcaption{実験 2, 3 で用いた初期東西風速場}
    \label{fig:init_U2}
\end{minipage}
\end{tabular}
\begin{tabular}{cc}
\begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=6cm]{figs/init/Tzh.png}\subcaption{実験 1 で用いた初期温度場}
    \label{fig:init_T1}
\end{minipage}
\begin{minipage}[b]{0.5\linewidth}
  \centering
    \includegraphics[width=6cm]{figs/result2/init/Tzh.png}\subcaption{実験 2, 3 で用いた初期温度場}
    \label{fig:initT2}
\end{minipage}
\end{tabular}
\caption{初期場}
\label{fig:init}
\end{figure} 

初期温度場には, 初期東西風速場に釣り合う温度場を与えた. 
図 \ref{fig:init_U1}, \ref{fig:init_U2} はそれぞれ実験 1, 実験 2,3 で用いた東西方向に一様な初期温度場の子午面分布である. 
\begin{equation}
  T(\phi,\sigma)  = T_{eq}(\sigma) + \frac{1}{R}\left(\overline{u_{eq}(\sigma)}+a\Omega\right)\frac{\partial\overline{u_{eq}(\sigma)}}{\partial\ln\sigma}\sin^2\phi
\end{equation}
に従う. 
ここで, $\phi$ は緯度, $\sigma$ は $\sigma$ 高度, $T(\phi,\sigma)$ は温度, 
$T_{eq}(\sigma)$ は赤道における温度, $\overline{u_{eq}(\sigma)}$ は赤道における東西風, 
$\Omega$ は金星の自転角速度である. 

\section{物理定数}
物理定数として, 金星半径を 6051.8 km, 重力加速度を 8.87 $\rm{m/s^2}$, 自転角速度を $3.0\times 10^{-7}$ /s, 
定圧比熱を 1000 J/K, 乾燥空気の気体定数を 191.4 $\rm{J/(K\cdot mol)}$, 
地表面気圧を $9.2\times10^6$ $\rm{Pa}$ とした. 自転の方向は地球と同じ向きである. 

%\end{document}                  
%%%%%%%%              Text End                  %%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%