ray/calculate1Dscalar.hで用いている計算式

(Formula used in ray/calculate1Dscalar.h)
Last Update: 2021/12/6

1. 以下の説明で用いる記号 (Notation used in the description below)

\(z\)軸を鉛直下向きに取り、深さ\(z\)における速度を\(v(z)\)、 波線が\(z\)軸となす角を\(\theta(z)\)とする。 速度の定義点の深さを\(z=z_1,\cdots,z_N\)、 それらの深さでの速度を\(v_1,\cdots,v_N\)とする。 このとき \[\begin{equation} v_n=v(z_n), n=1,\cdots,N \label{eq.vn} \end{equation}\] である。\(z_n\leq z\leq z_{n+1}\) (\(n=1,\cdots,N-1\))における速度勾配を \[\begin{equation} s_n=\frac{v_{n+1}-v_n}{z_{n+1}-z_{n}} \label{eq.sn} \end{equation}\] とする。また\(z<z_1\)における速度勾配を\(s_0(>0)\)、 \(z>z_N\)における速度勾配を\(s_N(>0)\)とする。 波線パラメータを \[\begin{equation} p=\frac{\sin\theta(z)}{v(z)} \label{eq.p} \end{equation}\] とする。
We assume the \(z\) axis to lie downward, and \(v(z)\) and \(\theta(z)\) be the velocity and the angle between the ray and the \(z\) axis, respectively. We assume the velocities \(v_1,\cdots,v_N\) to be given at depths of \(z=z_1,\cdots,z_N\), suggesting eq. (\ref{eq.vn}). Let \(s_n\) (\(n=1,\cdots,N-1\)) be the velocity slope in \(z_n\leq z\leq z_{n+1}\) given by eq. (\ref{eq.sn}), and \(s_0(>0)\) and \(s_N(>0)\) be the velocity slopes in \(z<z_1\) and \(z>z_N\), respectively. We denote the ray parameter by \(p\), which satisfies eq. (\ref{eq.p}).
以上の記号とプログラムで用いている変数との関係は 下表の通りである。
The table below represents the relation between the notation in this documentation and variables in the program code.
本マニュアルでの記号
Notation in this documentation 		プログラム中の記号
Variables in the program code
\(N\)structure.N
\(z_n\) (\(n=1,\cdots,N\))structure.depths[n-1]
\(v_n\) (\(n=1,\cdots,N\))structure.velocities[n-1]
\(s_0\)structure.slope_shallow
\(s_N\)structure.slope_deep
\(p\)p

2. 波線の最大到達深度 (The maximum depth of the ray)

波線の最大到達深度は\(\sin\theta(z)=1\)となる\(z\)として定義される。 (\ref{eq.p})式より、この条件は \[\begin{equation} pv(z)=1 \label{eq.maxdepth} \end{equation}\] と等価である。 \(s_0>0\), \(s_N>0\)という条件を付けていることで \(v(z)\)は0から\(\infty\)までの全範囲の値を取るので、 どのような波線パラメータについても (\ref{eq.maxdepth})式を満たす\(z\)が最低1つは存在することが保証される。 (\ref{eq.maxdepth})式の解が複数存在する場合は その中の一番小さい\(z\)を波線の最大到達深度と見なす。 なぜならば、ひとたび角度が90°になったら それ以上深くは潜れないからである。
The maximum depth of the ray is defined as the solution of \(\sin\theta(z)=1\), which is equivalent to eq. (\ref{eq.maxdepth}) because of eq. (\ref{eq.p}). Since \(s_0>0\) and \(s_N>0\), \(v(z)\) takes all values from 0 to \(\infty\), which ensure that at least one solution of eq. (\ref{eq.maxdepth}) exists for an arbitrary ray parameter. When two or more solutions of eq. (\ref{eq.maxdepth}) exists, the smallest \(z\) among these solutions is adopted as the maximum depth because once the angle becomes 90°, the ray does not go deeper.
波線の最大到達深度は関数ray_max_depth_1Dscalarにおいて 以下の3つのケースに場合分けして計算される。
The maximum depth of the ray is computed by function ray_max_depth_1Dscalar switching for the following three cases.

3. 速度が深さの1次関数で表される2つの深さの間の 波線に沿った水平距離・実距離・走時 (Horizontal distance, travel distance, and travel time along the ray from one depth to another between which the velocity is a linear function of depth)

2つの深さ\(z=z_s\), \(z=z_d\) (\(z_s<z_d\))があり、 \(z_s<z<z_d\)の範囲では 速度が深さの1次関数(速度勾配\(s\))で与えられているとする。 また、波線の最大到達深度は\(z\geq z_d\)の範囲にあるものとする。 この場合に\(z=z_s\)を出発した波線が最初に\(z=z_d\)に到達するまでの 波線に沿った水平距離および実距離と走時を求めることを考える。
Let \(z=z_s\) and \(z=z_d\) be two depths (\(z_s<z_d\)), between which the velocity is a linear function of depth with the slope \(s\), We assume the maximum depth of the ray to be \(z\geq z_d\). Then, let us calculate the horizontal and travel distances as well as the travel time along the ray as it travels from \(z=z_s\) and first reaches \(z=z_d\).

3.1. 水平距離 (Horizontal distance)

一般に、深さ\(z\)において波線が微小距離\(\Delta z\)だけ深さ方向に進むときの 水平方向距離は \[\begin{equation} \Delta x =\Delta z\tan\theta(z) =\Delta z\frac{\sin\theta(z)}{\sqrt{1-\sin^2\theta(z)}} \label{eq.dx.theta} \end{equation}\] と表される。この関係は(\ref{eq.p})式を用いると \[\begin{equation} \Delta x=\Delta z\frac{pv(z)}{\sqrt{1-p^2v(z)^2}} \label{eq.dx.v} \end{equation}\] と変形できるので、これを\(z=z_s\)から\(z=z_d\)まで積分することにより、 水平距離の計算式として \[\begin{equation} X_{z_s\rightarrow z_d}=\int_{z_s}^{z_d}\frac{pv(z)}{\sqrt{1-p^2v(z)^2}}dz \label{eq.Xsd} \end{equation}\] を得る。
In general, the horizontal distance along the ray as it moves a small vertical distance \(\Delta z\) at a depth \(z\) is given by eq. (\ref{eq.dx.theta}). Inserting eq. (\ref{eq.p}) into (\ref{eq.dx.theta}) results in (\ref{eq.dx.v}), and integrating this equation from \(z=z_s\) to \(z=z_d\) gives (\ref{eq.Xsd}) as the formula for the horizontal distance.
この積分を実行するにあたり、\(z_s<z<z_d\)の範囲で \(v(z)\)が\(z\)の1次関数で書けることを利用する。 \(v(z_s)=v_s\), \(v(z_d)=v_d\)とおき、 速度勾配が\(dv/dz=s\)であることを用いると \[\begin{eqnarray} X_{z_s\rightarrow z_d} &=& \int_{v_s}^{v_d}\frac{pv}{\sqrt{1-p^2v^2}}\frac{dz}{dv}dv \nonumber \\ &=& -\frac{1}{sp}\int_{v_s}^{v_d}\frac{-p^2v}{\sqrt{1-p^2v^2}}dv \nonumber \\ &=& -\frac{1}{sp}\left[\sqrt{1-p^2v^2}\right]_{v_s}^{v_d} \nonumber \\ &=& -\frac{1}{sp}\left(\sqrt{1-p^2v_d^2}-\sqrt{1-p^2v_s^2}\right) \nonumber \\ &=& \frac{\sqrt{1-p^2v_s^2}-\sqrt{1-p^2v_d^2}}{sp} \label{eq.Xsd.integrated} \end{eqnarray}\] と積分を解析的に計算できる。
Now, remember that \(v(z)\) is a linear function of \(z\) in the \(z_s<z<z_d\) range. Using \(v(z_s)=v_s\), \(v(z_d)=v_d\), and the velocity slope \(dv/dz=s\), the integral of eq. (\ref{eq.Xsd}) can be analytically calculated as eq. (\ref{eq.Xsd.integrated}).
(\ref{eq.Xsd.integrated})式は\(v_s=v_d\)のときに0/0の不定値になってしまう。 これを避けるため、\(s=(v_d-v_s)/(z_d-z_s)\)であることを用いて変形を行うと \[\begin{eqnarray} X_{z_s\rightarrow z_d} &=& \frac{\sqrt{1-p^2v_s^2}-\sqrt{1-p^2v_d^2}}{sp} \frac{\sqrt{1-p^2v_s^2}+\sqrt{1-p^2v_d^2}} {\sqrt{1-p^2v_s^2}+\sqrt{1-p^2v_d^2}} \nonumber \\ &=& \frac{(1-p^2v_s^2)-(1-p^2v_d^2)} {sp(\sqrt{1-p^2v_s^2}+\sqrt{1-p^2v_d^2})} \nonumber \\ &=& \frac{p^2(v_d^2-v_s^2)} {\frac{v_d-v_s}{z_d-z_s}p(\sqrt{1-p^2v_s^2}+\sqrt{1-p^2v_d^2})} \nonumber \\ &=& \frac{p(v_d+v_s)(z_d-z_s)} {\sqrt{1-p^2v_s^2}+\sqrt{1-p^2v_d^2}} \label{eq.Xsd.final} \end{eqnarray}\] を得る。 この式ならば\(v_s=v_d\)の場合も適用可能であり、 (\ref{eq.Xsd})式で\(v(z)=v_s=v_d\)(定数)とおいて直接得られる式と一致する。
Eq. (\ref{eq.Xsd.integrated}) gives 0/0 when \(v_s=v_d\). This problem is avoided by arranging the equation as (\ref{eq.Xsd.final}), where \(s=(v_d-v_s)/(z_d-z_s)\) was used. Eq. (\ref{eq.Xsd.final}) is applicable to the case of \(v_s=v_d\), and is consistent with the direct result from eq. (\ref{eq.Xsd}) assuming \(v(z)=v_s=v_d\) (constant).
関数ray_horizontal_distance_1Dscalar_linearでは (\ref{eq.Xsd.final})式を用いて水平距離の計算を行う。
Function ray_horizontal_distance_1Dscalar_linear computes the horizontal distance using eq. (\ref{eq.Xsd.final}).

3.2. 実距離 (Travel distance)

深さ\(z\)において波線が微小距離\(\Delta z\)だけ深さ方向に進むとき、 水平方向を考慮に入れた波線に沿った微小距離は \[\begin{equation} \Delta l =\frac{\Delta z}{\cos\theta(z)} =\frac{\Delta z}{\sqrt{1-\sin^2\theta(z)}} \label{eq.dl.theta} \end{equation}\] と表される。(\ref{eq.p})式を代入すると \[\begin{equation} \Delta l=\frac{\Delta z}{\sqrt{1-p^2v(z)^2}} \label{eq.dl.v} \end{equation}\] となる。これを\(z=z_s\)から\(z=z_d\)まで積分すれば、 \(z=z_s\)を出発した波線が最初に\(z=z_d\)に到達するまでの波線に沿った実距離は \[\begin{equation} L_{z_s\rightarrow z_d}=\int_{z_s}^{z_d}\frac{1}{\sqrt{1-p^2v(z)^2}}dz \label{eq.Lsd} \end{equation}\] と書ける。
The total distance along the ray as it moves a small vertical distance \(\Delta z\) at a depth \(z\) is given by eq. (\ref{eq.dl.theta}). Inserting (\ref{eq.p}) into (\ref{eq.dl.theta}) results in (\ref{eq.dl.v}). Integrating this equation from \(z=z_s\) to \(z=z_d\) yields (\ref{eq.Lsd}) as the formula for the travel distance of the ray that propagates from the depth \(z=z_s\) to \(z=z_d\).
この積分を実行するにあたり、\(z_s<z<z_d\)の範囲で \(v(z)\)が\(z\)の1次関数で書けることを利用する。 \(v(z_s)=v_s\), \(v(z_d)=v_d\)とおき、 速度勾配が\(dv/dz=s\)であることを用いると \[\begin{eqnarray} L_{z_s\rightarrow z_d} &=& \int_{v_s}^{v_d}\frac{1}{\sqrt{1-p^2v^2}}\frac{dz}{dv}dv \nonumber \\ &=& \frac{1}{s}\int_{v_s}^{v_d}\frac{1}{\sqrt{1-p^2v^2}}dv \label{eq.Lsd.use_v_integral} \end{eqnarray}\] と変形できる。更に\(v_s=\sin{\theta}_s/p\), \(v_d=\sin{\theta}_d/p\)とおき、 \(v=\sin\theta/p\)の変数変換を施すと \[\begin{eqnarray} L_{z_s\rightarrow z_d} &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \frac{1}{\sqrt{1-\sin^2\theta}} \frac{dv}{d\theta}d\theta \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \frac{1}{\cos\theta}\frac{\cos\theta}{p}d\theta \nonumber \\ &=& \frac{1}{sp}\int_{{\theta}_s}^{{\theta}_d}d\theta \nonumber \\ &=& \frac{{\theta}_d-{\theta}_s}{sp} \label{eq.Lsd.integrated} \end{eqnarray}\] となり、積分を解析的に計算できる。
Now, remember that \(v(z)\) is a linear function of \(z\) in the \(z_s<z<z_d\) range. Using \(v(z_s)=v_s\), \(v(z_d)=v_d\), and the velocity slope \(dv/dz=s\), the integral of eq. (\ref{eq.Xsd}) can be arranged as (\ref{eq.Lsd.use_v_integral}). Further converting the integral variable from \(v\) to \(\theta\) using a relation \(v=\sin\theta/p\), with defining \(v_s=p\sin{\theta}_s\) and \(v_d=p\sin{\theta}_d\), the integral can be analytically calculated. The result is (\ref{eq.Lsd.integrated}).
(\ref{eq.Lsd.integrated})式は\(v_s=v_d\)のときに0/0の不定値になってしまう。 これを避けるため、\(s=(v_d-v_s)/(z_d-z_s)\)であることを用いて変形を行うと \[\begin{eqnarray} L_{z_s\rightarrow z_d} &=& \frac{{\theta}_d-{\theta}_s}{\frac{v_d-v_s}{z_d-z_s}p} \nonumber \\ &=& \frac{(z_d-z_s)({\theta}_d-{\theta}_s)}{(v_d-v_s)p} \nonumber \\ &=& \frac{(z_d-z_s)({\theta}_d-{\theta}_s)}{v_dp-v_sp} \nonumber \\ &=& \frac{(z_d-z_s)({\theta}_d-{\theta}_s)}{\sin{\theta}_d-\sin{\theta}_s} \nonumber \\ &=& \frac{(z_d-z_s)({\theta}_d-{\theta}_s)} {2\sin\frac{{\theta}_d-{\theta}_s}{2} \cos\frac{{\theta}_d+{\theta}_s}{2}} \nonumber \\ &=& \frac{z_d-z_s}{\cos\frac{{\theta}_d+{\theta}_s}{2}} \frac{\frac{{\theta}_d-{\theta}_s}{2}}{\sin\frac{{\theta}_d-{\theta}_s}{2}} \label{eq.Lsd.final} \end{eqnarray}\] となる。\(v_d\rightarrow v_s\)とした極限では \({\theta}_d\rightarrow {\theta}_s\)となり、このとき \[\begin{equation} \frac{\frac{{\theta}_d-{\theta}_s}{2}}{\sin\frac{{\theta}_d-{\theta}_s}{2}} \rightarrow 1 \label{eq.theta_sintheta_ratio_limit} \end{equation}\] より\(L_{z_s\rightarrow z_d}\rightarrow (z_d-z_s)/\cos{\theta}_s\)となって 波線を直線にした場合の式と一致する。
Eq. (\ref{eq.Lsd.integrated}) gives 0/0 when \(v_s=v_d\). This problem is avoided by arranging the equation as (\ref{eq.Lsd.final}). At the limit of \(v_d\rightarrow v_s\), asymptotic behaviours \({\theta}_d\rightarrow {\theta}_s\) and (\ref{eq.theta_sintheta_ratio_limit}) hold, thus giving \(L_{z_s\rightarrow z_d}\rightarrow (z_d-z_s)/\cos{\theta}_s\), which is consistent with the direct result assuming a straight ray.
関数ray_travel_distance_1Dscalar_linearでは (\ref{eq.Lsd.final})式を用いて波線に沿った実距離の計算を行う。
Function ray_travel_distance_1Dscalar_linear computes the travel distance along the ray using eq. (\ref{eq.Lsd.final}).

3.3. 走時 (Travel time)

深さ\(z\)において波線が微小距離\(\Delta z\)だけ深さ方向に進むとき、 水平方向を考慮に入れた波線に沿った微小距離は(\ref{eq.dl.v})で表される。 この微小距離を進むのにかかる時間は \[\begin{equation} \Delta t =\frac{\Delta l}{v(z)} =\frac{\Delta z}{v(z)\sqrt{1-p^2v(z)^2}} \label{eq.dt.v} \end{equation}\] と書ける。これを\(z=z_s\)から\(z=z_d\)まで積分すれば、 \(z=z_s\)を出発した波線が最初に\(z=z_d\)に到達するまでの所要時間は \[\begin{equation} T_{z_s\rightarrow z_d}=\int_{z_s}^{z_d}\frac{1}{v(z)\sqrt{1-p^2v(z)^2}}dz \label{eq.Tsd} \end{equation}\] と書ける。
The total distance along the ray as it moves a small vertical distance \(\Delta z\) at a depth \(z\) is given by eq. (\ref{eq.dl.v}). The time needed for the wave to propagate by this distance is given by eq. (\ref{eq.dt.v}). Integrating this equation from \(z=z_s\) to \(z=z_d\) yields (\ref{eq.Tsd}) as the formula for the travel time of the ray that propagates from the depth \(z=z_s\) to \(z=z_d\).
この積分を実行するにあたり、\(z_s<z<z_d\)の範囲で \(v(z)\)が\(z\)の1次関数で書けることを利用する。 \(v(z_s)=v_s\), \(v(z_d)=v_d\)とおき、 速度勾配が\(dv/dz=s\)であることを用いると \[\begin{eqnarray} T_{z_s\rightarrow z_d} &=& \int_{v_s}^{v_d}\frac{1}{v\sqrt{1-p^2v^2}}\frac{dz}{dv}dv \nonumber \\ &=& \frac{1}{s}\int_{v_s}^{v_d}\frac{1}{v\sqrt{1-p^2v^2}}dv \label{eq.Tsd.use_v_integral} \end{eqnarray}\] と変形できる。更に\(v_s=\sin{\theta}_s/p\), \(v_d=\sin{\theta}_d/p\)とおき、 \(v=\sin\theta/p\)の変数変換を施すと \[\begin{eqnarray} T_{z_s\rightarrow z_d} &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \frac{1}{\frac{\sin\theta}{p}\sqrt{1-\sin^2\theta}} \frac{dv}{d\theta}d\theta \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \frac{p}{\sin\theta\cos\theta}\frac{\cos\theta}{p}d\theta \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d}\frac{d\theta}{\sin\theta} \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \frac{d\theta}{2\sin\frac{\theta}{2}\cos\frac{\theta}{2}} \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d}\frac{1}{2} \left(\frac{\cos\frac{\theta}{2}}{\sin\frac{\theta}{2}} +\frac{\sin\frac{\theta}{2}}{\cos\frac{\theta}{2}}\right)d\theta \nonumber \\ &=& \frac{1}{s}\int_{{\theta}_s}^{{\theta}_d} \left(\frac{1}{2}\frac{\cos\frac{\theta}{2}}{\sin\frac{\theta}{2}} -\frac{1}{2}\frac{-\sin\frac{\theta}{2}}{\cos\frac{\theta}{2}} \right)d\theta \nonumber \\ &=& \frac{1}{s}\left[ \ln\left(\sin\frac{\theta}{2}\right)-\ln\left(\cos\frac{\theta}{2}\right) \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\sqrt{\frac{1-\cos\theta}{2}} -\ln\sqrt{\frac{1+\cos\theta}{2}} \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\sqrt{\frac{1-\cos\theta}{1+\cos\theta}} \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\sqrt{\frac{(1-\cos\theta)^2}{(1+\cos\theta)(1-\cos\theta)}} \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\sqrt{\frac{(1-\cos\theta)^2}{1-\cos^2\theta}} \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\sqrt{\frac{(1-\cos\theta)^2}{\sin^2\theta}} \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\left(\frac{1-\cos\theta}{\sin\theta}\right) \right]_{{\theta}_s}^{{\theta}_d} \nonumber \\ &=& \frac{1}{s}\left[ \ln\left(\frac{1-\cos\theta_d}{\sin\theta_d}\right) -\ln\left(\frac{1-\cos\theta_d}{\sin\theta_d}\right) \right] \nonumber \\ &=& \frac{1}{s}\left[ \ln\left(\frac{1-\sqrt{1-p^2v_d^2}}{pv_d}\right) -\ln\left(\frac{1-\sqrt{1-p^2v_s^2}}{pv_s}\right) \right] \nonumber \\ &=& \frac{1}{s}\left[ \ln\left(\frac{1}{pv_d}-\sqrt{\frac{1}{p^2v_d^2}-1}\right) -\ln\left(\frac{1}{pv_s}-\sqrt{\frac{1}{p^2v_s^2}-1}\right) \right] \label{eq.Tsd.integrated} \end{eqnarray}\] となり、積分を解析的に計算できる。
Now, remember that \(v(z)\) is a linear function of \(z\) in the \(z_s<z<z_d\) range. Using \(v(z_s)=v_s\), \(v(z_d)=v_d\), and the velocity slope \(dv/dz=s\), the integral of eq. (\ref{eq.Xsd}) can be arranged as (\ref{eq.Tsd.use_v_integral}). Further converting the integral variable from \(v\) to \(\theta\) using a relation \(v=\sin\theta/p\), with defining \(v_s=p\sin{\theta}_s\) and \(v_d=p\sin{\theta}_d\), the integral can be analytically calculated. The result is (\ref{eq.Tsd.integrated}).
(\ref{eq.Tsd.integrated})式は\(v_d=v_s\)の場合に0/0の不定値になってしまう。 また0/0が現れない形に変形するのも容易ではない。 そこで、(\ref{eq.Tsd.integrated})式において\(v_d\rightarrow v_s\) とした極限での振る舞いを考える。 そのために関数 \[\begin{equation} f(v)\equiv\ln\left(\frac{1}{pv}-\sqrt{\frac{1}{p^2v^2}-1}\right) \label{eq.f} \end{equation}\] を定義するのが便利である。この関数\(f\)の微分形は \[\begin{eqnarray} f'(v) &=& \frac{1}{\frac{1}{pv}-\sqrt{\frac{1}{p^2v^2}-1}}\left[ -\frac{1}{pv^2} -\frac{1}{2\sqrt{\frac{1}{p^2v^2}-1}}\left(-\frac{2}{p^2v^3}\right) \right] \nonumber \\ &=& \frac{1}{\frac{1}{pv}-\sqrt{\frac{1}{p^2v^2}-1}}\left( -\frac{1}{pv^2} +\frac{1}{p^2v^3\sqrt{\frac{1}{p^2v^2}-1}} \right) \nonumber \\ &=& \frac{1}{\frac{1}{pv}-\sqrt{\frac{1}{p^2v^2}-1}} \frac{1}{pv^2\sqrt{\frac{1}{p^2v^2}-1}}\left( -\sqrt{\frac{1}{p^2v^2}-1}+\frac{1}{pv} \right) \nonumber \\ &=& \frac{1}{pv^2\sqrt{\frac{1}{p^2v^2}-1}} \nonumber \\ &=& \frac{1}{v\sqrt{1-p^2v^2}} \label{eq.fdash} \end{eqnarray}\] \[\begin{eqnarray} f''(v) &=& \left(-\frac{1}{v^2}\right)\frac{1}{\sqrt{1-p^2v^2}} +\frac{1}{v}\left[-\frac{1}{2}\left(1-p^2v^2\right)^{-3/2}\right] \left(-2p^2v\right) \nonumber \\ &=& -\frac{1}{v^2\sqrt{1-p^2v^2}}+\frac{p^2}{\left(1-p^2v^2\right)^{3/2}} \nonumber \\ &=& \frac{-\left(1-p^2v^2\right)+p^2v^2}{v^2\left(1-p^2v^2\right)^{3/2}} \nonumber \\ &=& \frac{2p^2v^2-1}{v^2\left(1-p^2v^2\right)^{3/2}} \label{eq.f2dash} \end{eqnarray}\] \[\begin{eqnarray} f^{(3)}(v) &=& \frac{4p^2v\cdot v^2\left(1-p^2v^2\right)^{3/2}-\left(2p^2v^2-1\right)} {v^4\left(1-p^2v^2\right)^3} \nonumber \\ & & \left[2v\left(1-p^2v^2\right)^{3/2} +v^2\cdot\frac{3}{2}\left(1-p^2v^2\right)^{1/2}\left(-2p^2v\right) \right] \nonumber \\ &=& \frac{4p^2v^3\left(1-p^2v^2\right)^{3/2}-\left(2p^2v^2-1\right)} {v^4\left(1-p^2v^2\right)^3} \nonumber \\ & & \left[2v\left(1-p^2v^2\right)^{3/2} -3p^2v^3\left(1-p^2v^2\right)^{1/2} \right] \nonumber \\ &=& \frac{1}{v^4\left(1-p^2v^2\right)^3} \left[ 4p^2v^3\left(1-p^2v^2\right)^{3/2} -2v\left(2p^2v^2-1\right)\left(1-p^2v^2\right)^{3/2} \right. \nonumber \\ & & \left. +3p^2v^3\left(2p^2v^2-1\right)\left(1-p^2v^2\right)^{1/2} \right] \nonumber \\ &=& \frac{4p^2v^2\left(1-p^2v^2\right) -2\left(2p^2v^2-1\right)\left(1-p^2v^2\right) +3p^2v^2\left(2p^2v^2-1\right)} {v^3\left(1-p^2v^2\right)^{5/2}} \nonumber \\ &=& \frac{1}{v^3\left(1-p^2v^2\right)^{5/2}} \left[ \left(4p^2v^2-4p^4v^4\right) -\left(4p^2v^2-4p^4v^4-2+2p^2v^2\right) \right. \nonumber \\ & & \left. +\left(6p^4v^4-3p^2v^2\right) \right] \nonumber \\ &=& \frac{2-5p^2v^2+6p^4v^4}{v^3\left(1-p^2v^2\right)^{5/2}} \label{eq.f3dash} \end{eqnarray}\] と書ける。関数\(f\)を用いると(\ref{eq.Tsd.integrated})式は \[\begin{eqnarray} T_{z_s\rightarrow z_d} &=& \frac{1}{s}\left[f(v_d)-f(v_s)\right] \nonumber \\ &=& \frac{z_d-z_s}{v_d-v_s}\left[f(v_d)-f(v_s)\right] \nonumber \\ &=& (z_d-z_s)\frac{f(v_d)-f(v_s)}{v_d-v_s} \label{eq.Tsd.use_f} \end{eqnarray}\] と書けるので、\(v_d\rightarrow v_s\)の極限では \[\begin{equation} T_{z_s\rightarrow z_d} =(z_d-z_s)f'(v_s) =\frac{z_d-z_s}{v_s\sqrt{1-p^2v_s^2}} \label{eq.Tsd.homogeneous_limit} \end{equation}\] となる。これは(\ref{eq.Tsd})式において\(v(z)\)を\(z\)によらない定数と見なして 積分した形になっている。 ここでは速度が深さの線形関数で書ける場合を考えているのであるから \(v_d\rightarrow v_s\)の極限では速度は定数となる。 したがって(\ref{eq.Tsd.integrated})式で\(v_d\rightarrow v_s\)とした極限 (\ref{eq.Tsd.homogeneous_limit})が (\ref{eq.Tsd})式で\(v(z)\)を定数と見なした場合の式と一致するのは 妥当な結果である。
Eq. (\ref{eq.Tsd.integrated}) yields 0/0 when \(v_d=v_s\), and arranging this equation to avoid 0/0 is difficult. Thus, let us consider an approximate form of eq. (\ref{eq.Tsd.integrated}) at the limit of \(v_d\rightarrow v_s\). A function \(f(v)\) defined by eq. (\ref{eq.f}) is convenient for this purpose. Differential forms of this function are given by eqs. (\ref{eq.fdash}), (\ref{eq.f2dash}), and (\ref{eq.f3dash}). Using this function, the travel time is given by (\ref{eq.Tsd.use_f}), yielding (\ref{eq.Tsd.homogeneous_limit}) at the limit of \(v_d\rightarrow v_s\). This equation is equivalent to the expression obtained by (\ref{eq.Tsd}) assuming that \(v(z)\) is not dependent on \(z\). Remember that our formulation here assumed a linear velocity change with depth. Thus the limit \(v_d\rightarrow v_s\) means a constant velocity. It is therefore reasonable that the limit of eq. (\ref{eq.Tsd.integrated}) at \(v_d\rightarrow v_s\), given by (\ref{eq.Tsd.homogeneous_limit}), is equivalent to the expression obtained by assuming a constant velocity in eq. (\ref{eq.Tsd}).
ところでプログラミングでは一般に2つの実数の等号比較は不安定の原因であり、 \(v_d=v_s\)のときとそうでないときで場合分けするのは良いコードとは言えない。 ある幅を持たせて\(v_d\)と\(v_s\)が十分に近い場合とそうでない場合で場合分けする方が より安定な計算ができる。 そのためには\(v_d\)が\(v_s\)と「完全には一致しないが十分近い」場合に成り立つ 近似式が必要になる。 そこで、\(f(v_d)\)を\(v_s\)のまわりでテイラー展開して3次の項まで残すと \[\begin{equation} f(v_d)\sim f(v_s)+f'(v_s)(v_d-v_s) +\frac{f''(v_s)}{2}(v_d-v_s)^2+\frac{f^{(3)}(v_s)}{6}(v_d-v_s)^3 \label{eq.fvd.taylor3} \end{equation}\] と書ける。このとき(\ref{eq.Tsd.use_f})式より \[\begin{eqnarray} T_{z_s\rightarrow z_d} &\sim& (z_d-z_s) \frac{f'(v_s)(v_d-v_s)+\frac{f''(v_s)}{2}(v_d-v_s)^2 +\frac{f^{(3)}(v_s)}{6}(v_d-v_s)^3}{v_d-v_s} \nonumber \\ &=& f'(v_s)(z_d-z_s)+\frac{f''(v_s)}{2}(z_d-z_s)(v_d-v_s) \nonumber \\ & & +\frac{f^{(3)}(v_s)}{6}(z_d-z_s)(v_d-v_s)^2 \label{eq.Tsd.taylor3} \end{eqnarray}\] となり、(\ref{eq.fdash})-(\ref{eq.f3dash})を代入して \[\begin{eqnarray} T_{z_s\rightarrow z_d} &\sim& \frac{z_d-z_s}{v_s\sqrt{1-p^2v_s^2}} +\frac{(2p^2v_s^2-1)(z_d-z_s)} {2v_s^2\left(1-p^2v_s^2\right)^{3/2}}(v_d-v_s) \nonumber \\ & & +\frac{(2-5p^2v_s^2+6p^4v_s^4)(z_d-z_s)} {6v_s^3\left(1-p^2v_s^2\right)^{5/2}} (v_d-v_s)^2 \label{eq.Tsd.approximate3} \end{eqnarray}\] を得る。この最後の項を無視できる条件は \[\begin{equation} \left|\frac{(2p^2v_s^2-1)(z_d-z_s)} {2v_s^2\left(1-p^2v_s^2\right)^{3/2}}(v_d-v_s)\right| \gg \left|\frac{(2-5p^2v_s^2+6p^4v_s^4)(z_d-z_s)} {6v_s^3\left(1-p^2v_s^2\right)^{5/2}}(v_d-v_s)^2\right| \label{eq.approximate2_condition.original} \end{equation}\] であり、共通因子を消去すると \[\begin{equation} \left|2p^2v_s^2-1\right| \gg \left|\frac{2-5p^2v_s^2+6p^4v_s^4} {3v_s\left(1-p^2v_s^2\right)}(v_d-v_s)\right| \label{eq.approximate2_condition.arrange1} \end{equation}\] 整理して \[\begin{eqnarray} \left|\frac{v_d-v_s}{3v_s}\right| &\ll& \left|\frac{\left(1-p^2v_s^2\right)\left(2p^2v_s^2-1\right)} {2-5p^2v_s^2+6p^4v_s^4}\right| \nonumber \\ &=& \left|\frac{-1+3p^2v_s^2-2p^4v_s^4}{2-5p^2v_s^2+6p^4v_s^4}\right| \nonumber \\ &=& \left|\frac{1-3p^2v_s^2+2p^4v_s^4}{2-5p^2v_s^2+6p^4v_s^4}\right| \label{eq.approximate2_condition.final} \end{eqnarray}\] を得る。この条件が成り立つとき、(\ref{eq.Tsd.approximate3})式の 最後の項を無視した近似式 \[\begin{equation} T_{z_s\rightarrow z_d} \sim \frac{z_d-z_s}{v_s\sqrt{1-p^2v_s^2}} +\frac{(2p^2v_s^2-1)(z_d-z_s)} {2v_s^2\left(1-p^2v_s^2\right)^{3/2}}(v_d-v_s) \label{eq.Tsd.approximate2} \end{equation}\] を利用できる。
In the computer programs, two real numbers should not be compared by an equal sign because it would result in a numerical instability. Thus branching based on whether \(v_d=v_s\) or not is not a good code. Rather, the branching should be based on whether \(v_d\) and \(v_s\) are sufficiently close or not. To implement this branching, an approximated equation for the travel time that is valid when \(v_d\) and \(v_s\) are ``sufficiently close to each other but not exactly same'' is needed. We now derive this equation. A Taylor expansion of \(f(v_d)\) with respect to \(v_s\), remaining until the 3rd-order term, results in eq. (\ref{eq.fvd.taylor3}). Using this relation, eq. (\ref{eq.Tsd.use_f}) reduces to (\ref{eq.Tsd.taylor3}), and inserting eqs. (\ref{eq.fdash})-(\ref{eq.f3dash}) into it gives (\ref{eq.Tsd.approximate3}). The final term of this equation can be neglected when (\ref{eq.approximate2_condition.original}) holds. Removing the common factors from the both hand sides of (\ref{eq.approximate2_condition.original}) results in (\ref{eq.approximate2_condition.arrange1}), which can further be arranged as (\ref{eq.approximate2_condition.final}). When this condition is met, it is valid to neglect the final term of (\ref{eq.Tsd.approximate3}), resulting in (\ref{eq.Tsd.approximate2}).
(\ref{eq.Tsd.integrated})式と(\ref{eq.Tsd.approximate2})式の比較の例を 図1に示す。 ここでは(\ref{eq.Tsd.integrated})式に基づいて計算した走時を \(T_{z_s\rightarrow z_d}^{exact}\)、 (\ref{eq.Tsd.approximate2})式に基づいて計算した走時を \(T_{z_s\rightarrow z_d}^{approx}\)、 (\ref{eq.approximate2_condition.final})式の左辺を右辺で割った値を \(M_{approx}\)とする。 理論上は\(M_{approx}\ll 1\)であれば近似が成り立ち、 \(M_{approx}\)の値が小さいほど理論上の近似精度が良い。 また\(T_{z_s\rightarrow z_d}^{exact}\)と \(T_{z_s\rightarrow z_d}^{approx}\)が近いほど 実際の近似精度が良い。 これらを比較することで、実際に\(M_{approx}\)がどの程度小さければ \(T_{z_s\rightarrow z_d}^{approx}\)がどれくらい良い近似になるのかを評価できる。 図1cは\(M_{approx}\)と近似式の誤差 \(|T_{z_s\rightarrow z_d}^{approx}-T_{z_s\rightarrow z_d}^{exact}|/ T_{z_s\rightarrow z_d}^{exact}\) との関係を示したものであり、\(M_{approx}=0.01\)であれば 近似式の誤差は\(10^{-4}\)のオーダーになることが分かる。 この誤差レベルは波線パラメータ\(p\)や層厚\(z_d-z_s\)にも依存すると思われるので、 プログラムでは安全側に余裕を見て\(M_{approx}<0.0001\)のときに 近似式を用いるようにしている。
Fig. 1 is a comparison of the travel times computed by eqs. (\ref{eq.Tsd.integrated}) and (\ref{eq.Tsd.approximate2}). Here, \(T_{z_s\rightarrow z_d}^{exact}\) and \(T_{z_s\rightarrow z_d}^{approx}\) are the travel times based on eqs. (\ref{eq.Tsd.integrated}) and (\ref{eq.Tsd.approximate2}), respectively, and \(M_{approx}\) is the ratio of left- to right-hand-side of eq. (\ref{eq.approximate2_condition.final}). The approximation is theoretically valid when \(M_{approx}\ll 1\); theoretically, the smaller the \(M_{approx}\), the better the approximation. The real precision of the approximation is evaluated by comparing \(T_{z_s\rightarrow z_d}^{exact}\) and \(T_{z_s\rightarrow z_d}^{approx}\). In Fig. 1c, a relation between the theoretical goodness of the approximation (\(M_{approx}\)) and the real goodness of the approximation (\(|T_{z_s\rightarrow z_d}^{approx}-T_{z_s\rightarrow z_d}^{exact}| /T_{z_s\rightarrow z_d}^{exact}\)) is shown, indicating that a value of \(M_{approx}=0.01\) results in a travel time error in order of \(10^{-4}\) by the approximation. This error level may depend on the ray parameter \(p\) and the layer thickness \(z_d-z_s\). To have a margin for the safety side, the program switches the exact and approximate equations at \(M_{approx}=0.0001\).

図1. \(z_s=0\) m, \(z_d=1000\) m, \(p=0.0001\) s m\(^{-1}\), \(v_s=2000\) m s\(^{-1}\) のときの(\ref{eq.Tsd.integrated})式と(\ref{eq.Tsd.approximate2})式の比較。 \(v_d\)を1800 m s\(^{-1}\)から2200 m s\(^{-1}\)まで 0.1 m s\(^{-1}\)刻みで動かして計算した。 (a)\(v_d\)と\(T_{z_s\rightarrow z_d}\)の関係。 赤線は(\ref{eq.Tsd.integrated})式による厳密解\(T_{z_s\rightarrow z_d}^{exact}\)、 青線は(\ref{eq.Tsd.approximate2})式による 近似解\(T_{z_s\rightarrow z_d}^{approx}\)。 (b)\(v_d\)と理論上の近似精度\(M_{approx}\)との関係。 この値が小さいほど理論上の近似精度が良い。 (c)\(M_{approx}\)と \(|T_{z_s\rightarrow z_d}^{approx}-T_{z_s\rightarrow z_d}^{exact}|/ T_{z_s\rightarrow z_d}^{exact}\) との関係。
Fig. 1. A comparison of travel times computed by eqs. (\ref{eq.Tsd.integrated}) and (\ref{eq.Tsd.approximate2}) for a case of \(z_s=0\) m, \(z_d=1000\) m, \(p=0.0001\) s m\(^{-1}\), \(v_s=2000\) m s\(^{-1}\), and moving \(v_d\) from 1800 m s\(^{-1}\) to 2200 m s\(^{-1}\) at an interval of 0.1 m s\(^{-1}\). (a) A relation between \(v_d\) and \(T_{z_s\rightarrow z_d}\). The red and blue lines represent the exact solution \(T_{z_s\rightarrow z_d}^{exact}\) from eq. (\ref{eq.Tsd.integrated}) and an approximate form \(T_{z_s\rightarrow z_d}^{approx}\) based on eq. (\ref{eq.Tsd.approximate2}), respectively. (b) A relation between \(v_d\) and \(M_{approx}\). The value of \(M_{approx}\) represents a theoretical goodness of the approximation; smaller \(M_{approx}\) values correspond to better approximations. (c) A relation between \(M_{approx}\) and \(|T_{z_s\rightarrow z_d}^{approx}-T_{z_s\rightarrow z_d}^{exact}|/ T_{z_s\rightarrow z_d}^{exact}\).

4. 速度定義点をまたぐ2つの深さの間の波線に沿った水平距離 (Horizontal distance along the ray between two depths over a velocity definition point)

前節では\(z_s<z<z_d\)において速度が深さの1次関数で表される場合を考えた。 これは\(z_s<z<z_d\)の範囲に速度定義点が存在しない状況に対応する。 \(z_s<z<z_d\)の範囲に速度定義点が含まれる場合、その地点で速度勾配が変わる。 このようなケースでは\(z_s<z<z_d\)の範囲を速度定義点で小区間に分割し、 それらの小区間ごとに(\ref{eq.Xsd.final})式を適用すれば良い。
In the previous section, the velocity was assumed to be a linear function in the \(z_s<z<z_d\) range, which corresponds to a situation where there is no velocity definition point in this range. If the velocity definition points are present in \(z_s<z<z_d\), the velocity slope changes at these points. In this case, eq. (\ref{eq.Xsd.final}) is applicable to subsections bounded by the velocity definition points.
いま、\(z_s\)よりも深部側の\(z_s\)に最も近い速度定義点を\(z=z_l\)、 \(z_d\)よりも浅部側の\(z_d\)に最も近い速度定義点を\(z=z_m\)とする。 \(z_s<z<z_d\)の間に少なくとも1つの速度定義点が含まれる状況を考えれば \(l\geq m\)となる。 このとき、水平距離は \[\begin{equation} X_{z_s\rightarrow z_d} =X_{z_s\rightarrow z_l} +\sum_{n=l}^{m-1}X_{z_n\rightarrow z_{n+1}}+X_{z_m\rightarrow z_d} \label{eq.Xsd.multilayer} \end{equation}\] となる。この式の個々の項(\(X_{z_s\rightarrow z_l}\)など)は (\ref{eq.Xsd.final})式によって計算できる。
Let \(z=z_l\) and \(z=z_m\) be the nearest velocity definition points downward from \(z_s\) and upward from \(z_d\), respectively. When at least one velocity definition point exists in \(z_s<z<z_d\), \(l\geq m\) holds. In this case, the horizontal distance is calculated by eq. (\ref{eq.Xsd.multilayer}). Individual terms in this equation (\(X_{z_s\rightarrow z_l}\), etc.) can be calculated with eq. (\ref{eq.Xsd.final}).
関数ray_horizontal_distance_1Dscalarでは (\ref{eq.Xsd.multilayer})式と(\ref{eq.Xsd.final})式を組み合わせて 水平距離の計算を行う。
Function ray_horizontal_distance_1Dscalar computes the horizontal distance using eqs (\ref{eq.Xsd.multilayer}) and (\ref{eq.Xsd.final}).