Return to 参考:Maxima 編(旧版アーカイブ)

Maxima で高さ h からの斜方投射の最大到達距離を求める

無次元化された解を書くと

\begin{eqnarray}
X &=& \cos\theta\cdot T \\
Y &=& H+ \sin\theta\cdot T – \frac{1}{2} T^2 \\
\end{eqnarray}

滞空時間 $T_1$ は

$$T_1 = \sin\theta + \sqrt{\sin^2\theta + 2 H}$$

この時間での水平方向の到達距離 $L$ は

$$L = X(T_1) = \cos\theta\cdot T_1$$

詳細は以下のページ

In [1]:
/* 度をラジアンに変換する関数を準備しておく */
radians(deg):= deg * %pi/180$
In [2]:
/* 角度 th は度で */
X(th, T):= cos(radians(th)) * T;
Y(th, T, H):= H + sin(radians(th)) * T - T**2/2;
Out[2]:
\[\tag{${\it \%o}_{2}$}X\left({\it th} , T\right):=\cos {\it radians}\left({\it th}\right)\,T\]
Out[2]:
\[\tag{${\it \%o}_{3}$}Y\left({\it th} , T , H\right):=H+\sin {\it radians}\left({\it th}\right)\,T+\frac{-T^2}{2}\]
In [3]:
T1(th, H):= sin(radians(th)) + sqrt(sin(radians(th))**2 + 2 * H);
L(th, H):= cos(radians(th)) * T1(th, H);
Out[3]:
\[\tag{${\it \%o}_{4}$}T_{1}\left({\it th} , H\right):=\sin {\it radians}\left({\it th}\right)+\sqrt{\sin ^2{\it radians}\left({\it th}\right)+2\,H}\]
Out[3]:
\[\tag{${\it \%o}_{5}$}L\left({\it th} , H\right):=\cos {\it radians}\left({\it th}\right)\,T_{1}\left({\it th} , H\right)\]

$H = 0$ の場合の軌道

$\theta = 45^{\circ}$ 前後の角度での斜方投射の軌道をグラフにしてみる。

In [4]:
/* θ = 45° のグラフのみ */
draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 1.1], yrange = [0, 0.5], grid = true,
  
  key = "θ = 45°", color = 3, 
  parametric(X(45,T), Y(45,T,0), T, 0, T1(45,0))
)$

In [5]:
/* 複数のグラフを一度に描く場合は... */
draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 1.1], yrange = [0, 0.5], grid = true,
  
  key = "θ = 47°", color = 1, 
  parametric(X(47,T), Y(47,T,0), T, 0, T1(47,0)),
  key = "θ = 46°", color = 2, 
  parametric(X(46,T), Y(46,T,0), T, 0, T1(46,0)),
  key = "θ = 45°", color = 3, 
  parametric(X(45,T), Y(45,T,0), T, 0, T1(45,0)),
  key = "θ = 44°", color = 4, 
  parametric(X(44,T), Y(44,T,0), T, 0, T1(44,0)),
  key = "θ = 43°", color = 5, 
  parametric(X(43,T), Y(43,T,0), T, 0, T1(43,0))
)$

In [6]:
/* 最初に makelist で作っておくといいかも */
lines: makelist([key = concat("θ = ", 48-i),
                 color = i, 
                 parametric(X(48-i,T), Y(48-i,T,0), T, 0, T1(48-i,0))], 
                i, 1, 5)$

draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 1.1], yrange = [0, 0.5], grid = true,

  lines
)$

In [7]:
/* 着地点付近を拡大表示して 45° が最大到達距離となることを確認 */
draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 1000, 
  /* 表示範囲を調整する */
  xrange = [0.996, 1.002], yrange = [0, 0.001], grid = true,

  lines
)$

$H = 0.1$ の場合の軌道

In [8]:
H: 0.1$
lines1: makelist([key = concat("θ = ", 48-i),
                 color = i, 
                 parametric(X(48-i,T), Y(48-i,T,H), T, 0, T1(48-i,H))], 
                i, 1, 10)$

draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 1.2], yrange = [0, 0.6], grid = true,

  lines1
)$

In [9]:
/* 着地点付近を拡大表示 */
draw2d(
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 1000, 
  /* 表示範囲を調整する */
  xrange = [1.08, 1.1], yrange = [0, 0.001], grid = true,

  lines1
)$

水平到達距離が最大となるのは $\theta=45^{\circ}$ のときではない!ことがわかるだろう。

水平到達距離 L

定量的に調べるために,水平到達距離 L(th, H) をグラフにしてみる。

In [10]:
H: 0.1$
draw2d(
  font = "Arial", font_size = 14, 
  grid = true, 
  xlabel = "打ち出し角度 (°)", ylabel = "規格化された水平到達距離",
  title = concat("H = ", H, "の場合の打ち出し角度と水平到達距離"),
  explicit(L(th, H), th, 38, 48)
)$

L が最大となる角度を数値微分で求める

数値微分して $\displaystyle \frac{dL(\theta, H)}{d\theta} = 0$ となる角度 $\theta$ を求めてみる。

In [11]:
/* 数値微分: h は小さいほどよい,というわけではない */
dL(th, H, h):= (L(th+h, H) - L(th-h, H))/(2*h)$
In [12]:
/* 方程式の数値解 */
fpprintprec: 0$
H: 0.1$
/* h を変えてみる */
for h in [1.e-3, 1.e-4, 1.e-5, 1.e-6, 1.e-7] do(
  print(find_root(dL(th, H, h) = 0, th, 38, 48))
)$
\(42.39204571350186\)
\(42.39204571469371\)
\(42.39204572066672\)
\(42.39204581678321\)
\(42.39204550128903\)
In [13]:
/* 方程式の数値解 */
H: 0.1$
thmax: 
find_root(dL(th, H, 1e-5) = 0, th, 38, 48)$

printf(true,"H = ~3,1f のとき,L が最大となる角度は ~9,6f °~%", H, thmax)$
printf(true,"    最大到達距離は ~9,6f", L(thmax, H))$
H = 0.1 のとき,L が最大となる角度は 42.392046 °
    最大到達距離は  1.095445

L が最大値となる角度を解析的微分で求める

$u \equiv \tan\theta$ として水平到達距離 $L$ を $u$ で表すと

\begin{eqnarray}
L(u, H)
&=& \frac{1}{1+u^2} \left\{u + \sqrt{(1+2H) u^2 + 2H} \right\}
\end{eqnarray}

In [14]:
kill(H)$
Lu(u, H):= 1/(1+u**2)*(u+sqrt((1+2*H)*u**2 + 2*H));
Out[14]:
\[\tag{${\it \%o}_{26}$}{\it Lu}\left(u , H\right):=\frac{1}{1+u^2}\,\left(u+\sqrt{\left(1+2\,H\right)\,u^2+2\,H}\right)\]
In [15]:
/* u に関する微分を解析的に */
define(dLu(u, H), diff(Lu(u, H), u));
Out[15]:
\[\tag{${\it \%o}_{27}$}{\it dLu}\left(u , H\right):=\frac{\frac{\left(2\,H+1\right)\,u}{\sqrt{\left(2\,H+1\right)\,u^2+2\,H}}+1}{u^2+1}-\frac{2\,u\,\left(\sqrt{\left(2\,H+1\right)\,u^2+2\,H}+u\right)}{\left(u^2+1\right)^2}\]
In [16]:
/* ラジアンを度に変換する関数も準備しておく */
degrees(theta):= theta*180/float(%pi)$

/* 微分がゼロとなる u は数値的に求める */
H: 0.1$
ans: find_root(dLu(u, H) = 0, u, 0, 1.2)$
thmax: degrees(atan(ans))$
printf(true,"H = ~3,1f のとき,L が最大となる角度は ~15,12f °~%", H, thmax)$
printf(true,"    最大到達距離は ~15,12f", Lu(ans, H))$
H = 0.1 のとき,L が最大となる角度は 42.392045714773 °
    最大到達距離は  1.095445115010
In [17]:
/* 上の結果を使ってグラフを描く */

draw2d(
  font = "Arial", font_size = 14, 
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 1.2], yrange = [0, 0.5], grid = true,
  xlabel = "X", ylabel = "Y",
  title = "高さ H からの斜方投射の最大水平到達距離",
  
  /* 数値を取り込んだ凡例の設定例 */
  key = 
    printf(false,"H = ~3,1f, θ = ~5,2f°, Lmax = ~5,3f",H,thmax,Lu(ans,H)),
  parametric(X(thmax, T), Y(thmax, T, H), T, 0, T1(thmax, H))
)$

問題 1

高さ $H$ からの斜方投射の水平到達距離が最大となる打ち出し角度と,そのときの最大水平到達距離がわかるようなグラフを次の $H$ についてまとめて描け。

$H = 0.2, 0.4, 0.6, 0.8, 1.0$

In [18]:
Hl: makelist(float(i/10), i, 2, 10, 2);
Out[18]:
\[\tag{${\it \%o}_{35}$}\left[ 0.2 , 0.4 , 0.6 , 0.8 , 1.0 \right] \]

解答例

In [19]:
/* 念のため */ kill(H)$

/* 到達距離が最大となる角度を追加するリスト */
thml:[]$
ansl:[]$

for H in Hl do(
  ans: find_root(dLu(u, H) = 0, u, 0, 1.2), 
  thmax: degrees(atan(ans)), 
  thml: append(thml, [thmax]), 
  ansl: append(ansl, [ans])
)$
In [20]:
/* 先に描く線をリストにしておく */
lines: makelist(
  [color=i, 
   key=printf(false, 
     "H = ~3,1f, θ = ~6,3f°, Lmax = ~5,3f", Hl[i],thml[i], Lu(ansl[i], Hl[i])),
   parametric(X(thml[i], T), Y(thml[i], T, Hl[i]), T, 0, T1(thml[i], Hl[i]))], 
  i, length(Hl), 1, -1
)$

draw2d(
  font = "Arial", font_size = 14, 
  user_preamble = "set key samplen 1;",
  /* 滑らかにするためにサンプリングを多めに */
  nticks = 100, 
  /* 縦横比。*/
  proportional_axes = xy,
  /* 表示範囲 */
  xrange = [0, 2], yrange = [0, 1.5], grid = true,
  xlabel = "X", ylabel = "Y",
  title = "高さ H からの斜方投射の最大水平到達距離",
  
  lines
)$

問題 2

高さ $H$ を横軸に,そのときの水平到達距離が最大となる打ち出し角度 $\theta$ を縦軸にしたグラフを描け。

解答例

In [21]:
/* リスト Hl を横軸,thml を縦軸に */
draw2d(
  points(Hl, thml)
)$

このままだと,あまりにショボいので,$H$ を与えると水平到達距離が最大となる打ち出し角度 $\theta$ を与える関数を作成し,陽関数として曲線のグラフにしてみる。

In [22]:
/* 念のため */ kill(H)$

thetamax(H):= block(
  assume(u > 0), /* H=0 の場合を解く際に必要 */
  ans: find_root(dLu(u, H) = 0, u, 0, 1.2),
  return(degrees(atan(ans)))
)$
In [23]:
draw2d(
  font = "Arial", font_size = 14, 
  /* 表示範囲 */
  xrange = [0, 1], yrange = [25, 45], grid = true,
  xlabel = "規格化された高さ H", ylabel = "角度 (°)",
  title = "高さ H からの斜方投射の水平到達距離を最大にする打出し角度",
  
  explicit(thetamax(H), H, 0, 1)
)$

参考:水平到達距離が最大になる角度

実は別ページで書いているように,規格化された高さ $H$ が与えられたとき,水平到達距離が最大となる角度 $\theta_{\rm max}$ は,以下の式で与えられることがわかっている。

$$\theta_{\rm max} = \arctan \frac{1}{\sqrt{1+2H}}$$

あえて2つの関数をグラフにしてみると,完全に重なっていることがわかるであろう。

問題 3

高さ $H$ を横軸に,そのときの最大水平到達距離を縦軸にしたグラフを描け。

解答例

In [24]:
/* 小数点表示を4桁に */
fpprintprec: 4$
Lumaxl: makelist(Lu(ansl[i], Hl[i]), i, 1, length(Hl))$

/* リスト Lumax の要素を1つずつ表示 */
for Lumax in Lumaxl do(
  print(Lumax)
)$
\(1.183\)
\(1.342\)
\(1.483\)
\(1.612\)
\(1.732\)
In [25]:
/* リスト Hl を横軸,Lumaxl を縦軸に */
draw2d(
  points(Hl, Lumaxl)
)$

このままだと,あまりにショボいので,$H$ を与えると最大水平到達距離を与える関数を作成し,陽関数として曲線のグラフにしてみる。

In [26]:
Lmax(H):= block(
  assume(u > 0), /* H=0 の場合を解く際に必要 */
  ans: find_root(dLu(u, H) = 0, u, 0, 1.2),
  return(Lu(ans, H))
)$
In [27]:
draw2d(
  font = "Arial", font_size = 14, 
  /* 表示範囲 */
  xrange = [0, 1], yrange = [1, 2], grid = true,
  xlabel = "規格化された高さ H", ylabel = "規格化された最大水平到達距離",
  title = "高さ H からの斜方投射の最大水平到達距離",
  
  explicit(Lmax(H), H, 0, 1)
)$

参考:最大水平到達距離

実は別ページで書いているように,規格化された高さ $H$ が与えられたときの最大水平到達距離は,規格化された量 $L_{\rm m}$ で以下の式で与えられることがわかっている。

$$L_{\rm m} = \sqrt{1+2H}$$

あえて2つの関数をグラフにしてみると,完全に重なっていることがわかるであろう。

問題 4

$h = 5\,\mbox{(m)}$ の高さから速さ $v_0 = 10\,\mbox{(m/s)}$ で斜方投射したときの最大水平到達距離は何 $\mbox{m}$ か。またそのときの打出し角度は何度か。

ヒント:

$h = 5\,\mbox{(m)}$ は規格化された高さ $H$ ではいくらか,また規格化された最大到達距離 $L_{\rm max}(H)$ を次元をもった量に直すといくらか,ということ。

\begin{eqnarray}
H &=& \frac{h}{\ell} = \frac{g h}{v_0^2} \\
x &=& \ell X = \frac{v_0^2}{g} X \\
x_{\rm max} &=& \frac{v_0^2}{g} L_{\rm max}(H)
\end{eqnarray}

重力加速度 $g = 9.80665 \,(\mbox{m}/\mbox{s}^2)$

解答例

In [28]:
h: 5$
v0: 10$
g: 9.80665$

H: h * g/v0**2;
Out[28]:
\[\tag{${\it \%o}_{55}$}0.4903\]
In [29]:
printf(true, 
  "最大到達距離は ~9,6f m~%",v0**2/g*Lmax(H))$
最大到達距離は 14.351088 m
In [30]:
printf(true, 
  "到達距離が最大となる角度は ~9,6f °~%",thetamax(H))$
到達距離が最大となる角度は 35.395688 °