物理のエッセンス

# 等加速度直線運動（Uniformly Accelerated Linear Motion）

等加速度直線運動とは、**物体が一定の加速度でまっすぐ進む運動**のことです。以下に基本式とそれぞれの意味・使い方を体系的にまとめます。

---

## ✅ 基本変数の定義

| 記号 | 意味         | 単位         |
|------|--------------|--------------|
| \( x \)     | 位置（変位）       | [m]         |
| \( v \)     | 速度             | [m/s]       |
| \( a \)     | 加速度（一定）     | [m/s²]      |
| \( t \)     | 時間             | [s]         |
| \( v_0 \)   | 初速度           | [m/s]       |
| \( x_0 \)   | 初期位置         | [m]         |

---

## ✅ 公式1：速度と時間の関係式（第1式）

\[
v = v_0 + a t
\]

- 一定加速度 \( a \) において、時間 \( t \) 経過後の速度
- 速度の「変化量」は \( at \)

---

## ✅ 公式2：位置と時間の関係式（第2式）

\[
x = x_0 + v_0 t + \frac{1}{2} a t^2
\]

- 初期位置 \( x_0 \) から、時間 \( t \) の経過後の位置
- 面積の和（台形＋三角形）としても解釈可能

---

## ✅ 公式3：速度と位置の関係式（時間を含まない）

\[
v^2 = v_0^2 + 2 a (x - x_0)
\]

- 位置の変化 \( x - x_0 \) を使って速度を計算
- 時間 \( t \) が不要な「エネルギー的」な式（仕事-エネルギー原理と接続）

---

## ✅ 公式4：平均速度と移動距離

\[
x - x_0 = \frac{v + v_0}{2} \cdot t
\]

- 平均速度（初速度と終速度の平均）× 時間
- 等加速度運動における**平均速度の定義式**

---

## ✅ 活用例

### 例1：物体が静止状態から落下（自由落下）

- 初速度 \( v_0 = 0 \)
- 加速度 \( a = g \approx 9.8\ \text{m/s}^2 \)
- 公式2より：

\[
x = \frac{1}{2} g t^2
\]

---

```python
# Program Name: uniform_acceleration.py
# Creation Date: 20250615
# Overview: Simulate uniformly accelerated linear motion and visualize position and velocity over time.
# Usage: Run the script to generate motion graphs for a selected initial condition.

# ------------------------
# 必要なライブラリをインストール / Install required libraries
# ------------------------
!pip install matplotlib numpy

# ------------------------
# ライブラリのインポート / Import libraries
# ------------------------
import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
x0 = 0.0          # 初期位置 Initial position [m]
v0 = 0.0          # 初速度 Initial velocity [m/s]
a  = 9.8          # 加速度 Acceleration [m/s^2]
t_end = 5.0       # シミュレーション時間 Simulation time [s]
dt = 0.01         # 時間刻み Time step [s]

# ------------------------
# 時間配列生成 / Create time array
# ------------------------
t = np.arange(0, t_end + dt, dt)

# ------------------------
# 各式に基づく運動の計算 / Compute motion using kinematic formulas
# ------------------------
v = v0 + a * t                               # 公式1：v = v0 + at
x = x0 + v0 * t + 0.5 * a * t**2             # 公式2：x = x0 + v0*t + (1/2)*a*t^2
v_sq = v0**2 + 2 * a * (x - x0)              # 公式3：v^2 = v0^2 + 2a(x - x0)
x_mean = ((v + v0) / 2) * t                  # 公式4：x - x0 = ((v + v0)/2) * t

# ------------------------
# グラフ描画 / Plot graphs
# ------------------------
# 速度 vs 時間 / Velocity vs Time
plt.figure()
plt.plot(t, v, label='Velocity v(t)')
plt.title("Velocity vs Time")
plt.xlabel("Time [s]")
plt.ylabel("Velocity [m/s]")
plt.grid(True)
plt.legend()

# 位置 vs 時間 / Position vs Time
plt.figure()
plt.plot(t, x, label='Position x(t)')
plt.title("Position vs Time")
plt.xlabel("Time [s]")
plt.ylabel("Position [m]")
plt.grid(True)
plt.legend()

# 平均速度による位置比較 / Position from average velocity
plt.figure()
plt.plot(t, x, label='x(t) from kinematic formula')
plt.plot(t, x_mean, '--', label='x(t) from average velocity')
plt.title("Comparison of Position vs Time")
plt.xlabel("Time [s]")
plt.ylabel("Position [m]")
plt.grid(True)
plt.legend()

# ------------------------
# グラフの表示 / Show plots
# ------------------------
plt.show()

```

---

## ✅ 放物線運動の基本式まとめ

| 項目            | 物理式                                            |
| ------------- | ---------------------------------------------- |
| 水平位置 $x(t)$   | $x = v_0 \cos\theta \cdot t$                   |
| 鉛直位置 $y(t)$   | $y = v_0 \sin\theta \cdot t - \frac{1}{2}gt^2$ |
| 水平速度 $v_x(t)$ | $v_x = v_0 \cos\theta$（一定）                     |
| 鉛直速度 $v_y(t)$ | $v_y = v_0 \sin\theta - g t$                   |
| 飛行時間 $T$      | $T = \frac{2 v_0 \sin\theta}{g}$               |
| 最大高さ $H$      | $H = \frac{(v_0 \sin\theta)^2}{2g}$            |
| 水平到達距離（射程）$R$ | $R = \frac{v_0^2 \sin(2\theta)}{g}$            |

---

```python
# Program Name: parabolic_motion_plot.py
# Creation Date: 20250615
# Overview: Simulate and plot the trajectory of a projectile under parabolic motion
# Usage: Run the script to visualize motion with initial parameters

# ------------------------
# 必要なライブラリをインストール / Install required libraries
# ------------------------
!pip install numpy matplotlib

import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# 一元管理される物理定数と初期条件 / Define constants and initial parameters
# ------------------------
g = 9.8  # 重力加速度 / gravitational acceleration [m/s^2]
v0 = 20.0  # 初速度 / initial velocity [m/s]
theta_deg = 45.0  # 発射角度 / launch angle [deg]
theta = np.radians(theta_deg)  # ラジアン変換 / convert to radians

# ------------------------
# 関数定義 / Define motion functions
# ------------------------

def x_t(t):
    # 水平位置 / horizontal position
    return v0 * np.cos(theta) * t

def y_t(t):
    # 鉛直位置 / vertical position
    return v0 * np.sin(theta) * t - 0.5 * g * t**2

def flight_time():
    # 飛行時間 / total flight time
    return 2 * v0 * np.sin(theta) / g

def max_height():
    # 最大高さ / maximum height
    return (v0 * np.sin(theta))**2 / (2 * g)

def range_distance():
    # 水平到達距離（射程）/ horizontal range
    return v0**2 * np.sin(2 * theta) / g

# ------------------------
# 結果の表示 / Print computed physical quantities
# ------------------------
print(f"Launch angle θ: {theta_deg} degrees")
print(f"Initial speed v₀: {v0} m/s")
print(f"Flight time T: {flight_time():.2f} s")
print(f"Maximum height H: {max_height():.2f} m")
print(f"Horizontal range R: {range_distance():.2f} m")

# ------------------------
# シミュレーション / Simulate motion
# ------------------------
T = flight_time()
t_vals = np.linspace(0, T, 100)  # 時間の分割 / time array
x_vals = x_t(t_vals)  # x座標配列 / x values
y_vals = y_t(t_vals)  # y座標配列 / y values

# ------------------------
# プロット / Plot the trajectory
# ------------------------
plt.figure(figsize=(10, 5))
plt.plot(x_vals, y_vals)

# 軸ラベルとタイトル / Axis labels and title
plt.title("Projectile Motion Trajectory")
plt.xlabel("Horizontal Position x [m]")
plt.ylabel("Vertical Position y [m]")
plt.grid(True)
plt.ylim(bottom=0)
plt.xlim(left=0)

# 注釈 / Annotation
plt.annotate(f"Max height ≈ {max(y_vals):.2f} m", 
             xy=(x_vals[np.argmax(y_vals)], max(y_vals)), 
             xytext=(x_vals[np.argmax(y_vals)], max(y_vals) + 1),
             arrowprops=dict(arrowstyle='->'))

plt.annotate(f"Range ≈ {x_vals[-1]:.2f} m", 
             xy=(x_vals[-1], 0), 
             xytext=(x_vals[-1] - 5, 2),
             arrowprops=dict(arrowstyle='->'))

# 表示 / Show plot
plt.show()

```

# ✅ 相対速度の基本式まとめ（Summary of Relative Velocity Equations）

---

## 🔹 基本定義（Definition）

| 表現                     | 数式                                              | 意味・用途                                 |
|------------------------|--------------------------------------------------|------------------------------------------|
| Aから見たBの速度         | \(\vec{v}_{\text{B/A}} = \vec{v}_B - \vec{v}_A\) | Aを基準としたときのBの見かけの運動                |
| Bから見たAの速度         | \(\vec{v}_{\text{A/B}} = \vec{v}_A - \vec{v}_B\) | Bを基準としたときのAの運動、上記の符号反転         |
| 相対速度の大きさ         | \(|\vec{v}_{\text{rel}}| = \sqrt{(v_x)^2 + (v_y)^2}\) | 2次元相対速度ベクトルのノルム                    |
| 衝突・追いつきの時間     | \(t = \frac{d}{|\vec{v}_A - \vec{v}_B|}\)       | AがBを追いかける際の時間（1次元）               |
| 相対運動の距離           | \(x_{\text{rel}}(t) = (\vec{v}_B - \vec{v}_A) t\) | 時間に応じた相対的な位置変化                     |

---

## 🔸 解説（Explanation）

### ● 1. 相対速度の意味

相対速度とは、「ある基準となる物体から見た他の物体の運動の速さと方向」を表す物理量。実生活の例：

- 電車Aに乗っている人から見た、横の線路を走る電車B → 相対速度で観察
- 船上の人から見た陸地の動き → 自分が動いているので、他が動いて見える

---

### ● 2. ベクトルとしての計算

- ベクトルの引き算によって、ある基準（観測者）から見た別の物体の速度を計算
- 相対速度は向きと大きさを持つので、**ベクトル演算が本質**

---

### ● 3. 1次元の特殊ケース

| 状況             | \(v_A\), \(v_B\)     | \(v_{B/A}\) の符号 | 解釈                   |
|----------------|----------------------|------------------|----------------------|
| 同じ方向・Aが速い   | \(v_A > v_B > 0\)    | 負                | BはAより遅いので後ろへ進むように見える |
| 同じ方向・Bが速い   | \(v_B > v_A > 0\)    | 正                | BはAを追い抜くように見える         |
| 逆向き           | \(v_A > 0, v_B < 0\) | 大きな正値または負値 | 正面衝突方向に動いているように見える   |

---

### ● 4. 応用：追いつく時間と衝突判定

- 初期距離 \(d\)、相対速度 \(v_r = v_A - v_B\)
- Aが速くて同方向なら：
  \[
  t = \frac{d}{v_A - v_B}
  \]
- \(v_A \leq v_B\) のときは**追いつけない**

---

## 🔹 図解イメージ（対応するプロット）

- ベクトルで可視化すると、直感的に「どちらが速いか」「どの方向に進んでいるか」が分かる
- 上記コードでは `matplotlib.quiver()` を用いてベクトル描画済

---


```python
# Program Name: relative_velocity_demo.py
# Creation Date: 20250615
# Overview: Visualize relative velocity vectors between two objects in 2D
# Usage: Run script to display relative velocity plot with vector arrows

# ----------------------------
# ライブラリのインストール / Install required libraries
# ----------------------------
!pip install matplotlib numpy

import numpy as np
import matplotlib.pyplot as plt

# ----------------------------
# パラメータ定義 / Define Parameters
# ----------------------------
# 物体AとBの速度ベクトル（x成分, y成分）/ Velocity vectors of A and B
v_A = np.array([2.0, 1.0])   # Object A velocity [vx, vy]
v_B = np.array([4.0, 3.0])   # Object B velocity [vx, vy]

# 相対速度計算 / Relative velocity (B with respect to A)
v_rel = v_B - v_A

# ----------------------------
# プロット設定 / Plot Setup
# ----------------------------
plt.figure(figsize=(8, 6))
origin = np.array([[0, 0], [0, 0]])  # 原点ベクトル描画用

# ベクトル描画 / Draw velocity vectors
plt.quiver(*origin, [v_A[0], v_B[0]], [v_A[1], v_B[1]],
           angles='xy', scale_units='xy', scale=1, color=['blue', 'green'], label='v_A and v_B')

# 相対速度ベクトル描画 / Draw relative velocity vector
plt.quiver(v_A[0], v_A[1], v_rel[0], v_rel[1],
           angles='xy', scale_units='xy', scale=1, color='red', label='v_rel = v_B - v_A')

# ----------------------------
# 表示設定 / Plot Decorations
# ----------------------------
plt.xlim(0, 6)
plt.ylim(0, 5)
plt.grid(True)
plt.gca().set_aspect('equal')
plt.title('Relative Velocity Visualization')
plt.xlabel('vx')
plt.ylabel('vy')
plt.legend(['v_A', 'v_B', 'v_B/A'])
plt.text(v_A[0], v_A[1], 'A', fontsize=12)
plt.text(v_B[0], v_B[1], 'B', fontsize=12)

plt.show()

# ----------------------------
# 結果の数値出力 / Print numerical results
# ----------------------------
print(f"Velocity of A: v_A = {v_A}")
print(f"Velocity of B: v_B = {v_B}")
print(f"Relative velocity of B with respect to A: v_rel = {v_rel}")
print(f"Relative speed: |v_rel| = {np.linalg.norm(v_rel):.2f}")

```


# 斜面上の力：\( mg \sin\theta \) と \( mg \cos\theta \) の物理的意味と式まとめ

物体が**傾斜角 \(\theta\)** の斜面上にあるとき、重力 \( mg \) は**斜面に沿う方向**と**斜面に垂直な方向**に分解して扱います。

---

## ✅ 基本設定

- 物体の質量：\( m \)
- 重力加速度：\( g \)
- 傾斜角：\( \theta \)（地面との角度）
- 重力：下向き \( mg \)

---

## ✅ 重力の成分分解

| 成分方向             | 式                            | 物理的意味                      |
|----------------------|-------------------------------|----------------------------------|
| 斜面に**平行**な成分   | \( F_{\parallel} = mg \sin\theta \) | 物体を**斜面に沿って滑らせる力** |
| 斜面に**垂直**な成分   | \( F_{\perp} = mg \cos\theta \)  | 斜面に**押し付ける力（垂直抗力）**  |

---

## ✅ 図解イメージ（力の分解）

```

```
↑ N（垂直抗力 = mg cosθ）
 │
 │  ／← 物体
 │／θ
↓ mg（重力）
 ↙    ↘
```

mg cosθ   mg sinθ →（斜面方向へ加速させる力）

```

---

## ✅ 運動方程式（斜面方向）

- **摩擦なし・加速度あり**の場合（斜面方向に運動）：

\[
F = ma = mg \sin\theta
\quad \Rightarrow \quad a = g \sin\theta
\]

- **静止 or 等速運動**（力の釣り合い）なら：

\[
mg \sin\theta = F_{\text{他力（摩擦など）}}
\]

---

## ✅ 垂直抗力の式（床からの反力）

- 斜面に押し付ける力は：

\[
N = mg \cos\theta
\]

摩擦力を求めるときの**基礎となる力**です。

---

## ✅ 摩擦力との関係（静止・滑り始め）

- 静止摩擦力の最大：

\[
f_{\text{max}} = \mu_s N = \mu_s mg \cos\theta
\]

- 滑り出す条件：

\[
mg \sin\theta > \mu_s mg \cos\theta
\quad \Rightarrow \quad \tan\theta > \mu_s
\]

---

## ✅ ポイントまとめ

| 観点           | キー式                     | 意味                                 |
|----------------|----------------------------|--------------------------------------|
| 斜面に沿う加速度   | \( a = g \sin\theta \)      | 摩擦がなければこの加速度で落下             |
| 斜面を押す力      | \( N = mg \cos\theta \)     | 摩擦力や圧力に関与                      |
| 静止条件の境界角度 | \( \tan\theta = \mu_s \)     | 静止 → 滑り出しの臨界条件を決める          |




```python
# Program Name: inclined_plane_forces.py
# Creation Date: 20250615
# Overview: Simulate and visualize mg*sin(θ) and mg*cos(θ) for inclined plane
# Usage: Run to plot gravitational components on an incline

# ------------------------
# 必要なライブラリのインストール / Install required libraries
# ------------------------
!pip install matplotlib numpy

# ------------------------
# 数値設定（物理定数とパラメータ） / Define parameters
# ------------------------
import numpy as np
import matplotlib.pyplot as plt

m = 1.0         # 質量 [kg] / Mass
g = 9.8         # 重力加速度 [m/s^2] / Gravitational acceleration
theta_deg = np.linspace(0, 90, 300)  # 角度範囲（度） / Angle range in degrees
theta_rad = np.radians(theta_deg)   # ラジアンに変換 / Convert to radians

# ------------------------
# 力の成分を計算 / Compute force components
# ------------------------
F_parallel = m * g * np.sin(theta_rad)    # 斜面方向の力 / Parallel component
F_perpendicular = m * g * np.cos(theta_rad)  # 垂直方向の力 / Perpendicular component

# ------------------------
# 結果をプリント / Print sample results
# ------------------------
print("Sample at θ = 30°:")
angle_sample = 30
theta_sample_rad = np.radians(angle_sample)
print(f"mg sin({angle_sample}) = {m * g * np.sin(theta_sample_rad):.2f} N")
print(f"mg cos({angle_sample}) = {m * g * np.cos(theta_sample_rad):.2f} N")

# ------------------------
# グラフ描画 / Plot results
# ------------------------
plt.figure(figsize=(10, 6))
plt.plot(theta_deg, F_parallel, label="mg sin(θ)", linewidth=2)
plt.plot(theta_deg, F_perpendicular, label="mg cos(θ)", linewidth=2)
plt.xlabel("Incline Angle θ [degrees]")
plt.ylabel("Force [N]")
plt.title("Gravitational Force Components on Inclined Plane")
plt.grid(True)
plt.legend()
plt.show()

```
以下に、**斜め2本のロープで支持された質点の釣り合い**を\*\*線形方程式の行列形式（クラメルの公式含む）\*\*で表現・解説します。

---

# ✅ 張力の釣り合いを行列で表現（クラメルの公式で解く）

## ◆ 状況設定（2つの張力と重力）

* 質量 $m$ の物体にかかる力：

  * 左ロープ：張力 $T_1$、角度 $\theta_1$
  * 右ロープ：張力 $T_2$、角度 $\theta_2$
  * 重力：$mg$、鉛直下向き

---

## ✅ 力の釣り合い：成分ごとの式

### 水平方向（x方向）：

$$
T_1 \cos\theta_1 - T_2 \cos\theta_2 = 0
\quad \text{……(1)}
$$

### 鉛直方向（y方向）：

$$
T_1 \sin\theta_1 + T_2 \sin\theta_2 = mg
\quad \text{……(2)}
$$

---

## ✅ 行列形式で表現（2元1次連立方程式）

未知数：$T_1, T_2$

$$
\begin{bmatrix}
\cos\theta_1 & -\cos\theta_2 \\
\sin\theta_1 &  \sin\theta_2
\end{bmatrix}
\begin{bmatrix}
T_1 \\
T_2
\end{bmatrix}
=
\begin{bmatrix}
0 \\
mg
\end{bmatrix}
$$

この形で線形連立方程式と見なせます。

---

## ✅ クラメルの公式で解く（Cramer’s Rule）

行列を定義：

* 係数行列 $A$：

  $$
  A =
  \begin{bmatrix}
  \cos\theta_1 & -\cos\theta_2 \\
  \sin\theta_1 &  \sin\theta_2
  \end{bmatrix}
  $$

* 定数ベクトル $\vec{b} = [0,\ mg]^T$

---

### 行列式 $\det A$：

$$
\det A = \cos\theta_1 \cdot \sin\theta_2 + \cos\theta_2 \cdot \sin\theta_1
= \sin(\theta_1 + \theta_2)
$$

---

### 解（クラメルの公式）：

$$
T_1 = \frac{
\begin{vmatrix}
0 & -\cos\theta_2 \\
mg & \sin\theta_2
\end{vmatrix}
}{\det A}
= \frac{mg \cdot \cos\theta_2}{\sin(\theta_1 + \theta_2)}
$$

$$
T_2 = \frac{
\begin{vmatrix}
\cos\theta_1 & 0 \\
\sin\theta_1 & mg
\end{vmatrix}
}{\det A}
= \frac{mg \cdot \cos\theta_1}{\sin(\theta_1 + \theta_2)}
$$

---

## ✅ 結論（一般解）

$$
\boxed{
T_1 = \frac{mg \cos\theta_2}{\sin(\theta_1 + \theta_2)}, \quad
T_2 = \frac{mg \cos\theta_1}{\sin(\theta_1 + \theta_2)}
}
$$

---

## ✅ 物理的意味

* **張力は互いの角度によって分担される**
* 角度が90度に近いと、**垂直成分が大きくなる → 張力が小さくて済む**

---

```python
# Program Name: tension_solver_cramers.py
# Creation Date: 20250615
# Overview: Solve T1 and T2 using matrix method (Cramer's rule) for static equilibrium
# Usage: Set θ1 and θ2 in degrees, then run the script to compute tensions

# ------------------------
# 必要なライブラリのインストール / Install required libraries
# ------------------------
!pip install numpy

import numpy as np

# ------------------------
# 物理パラメータの定義 / Define physical parameters
# ------------------------
m = 1.0       # 質量 [kg]
g = 9.8       # 重力加速度 [m/s^2]
theta1_deg = 30  # 左側の角度 θ1 [degrees]
theta2_deg = 45  # 右側の角度 θ2 [degrees]

# 弧度法に変換 / Convert to radians
theta1 = np.radians(theta1_deg)
theta2 = np.radians(theta2_deg)

# ------------------------
# 行列の構築 / Construct coefficient matrix and constants
# ------------------------
# A * [T1, T2]ᵗ = b
A = np.array([
    [np.cos(theta1), -np.cos(theta2)],
    [np.sin(theta1),  np.sin(theta2)]
])

b = np.array([0, m * g])  # 右辺ベクトル / Right-hand side

# ------------------------
# クラメルの公式による解法 / Solve using Cramer's Rule
# ------------------------
detA = np.linalg.det(A)

# T1 用の分子行列 / Replace first column with b
A1 = A.copy()
A1[:, 0] = b
detA1 = np.linalg.det(A1)

# T2 用の分子行列 / Replace second column with b
A2 = A.copy()
A2[:, 1] = b
detA2 = np.linalg.det(A2)

T1 = detA1 / detA
T2 = detA2 / detA

# ------------------------
# 結果出力 / Print results
# ------------------------
print(f"θ₁ = {theta1_deg}°, θ₂ = {theta2_deg}°")
print(f"Tension T₁ = {T1:.4f} N")
print(f"Tension T₂ = {T2:.4f} N")

# 理論解と一致するかチェック（参考値）
T1_theory = m * g * np.cos(theta2) / np.sin(theta1 + theta2)
T2_theory = m * g * np.cos(theta1) / np.sin(theta1 + theta2)
print(f"[Check] T₁ (theory) = {T1_theory:.4f} N")
print(f"[Check] T₂ (theory) = {T2_theory:.4f} N")

```

# ✅ 摩擦力のまとめ：式・分類・特徴・活用

摩擦力とは、**物体の接触面に沿って、運動を妨げる方向に働く力**のこと。  
「静止しているとき」と「運動中」では、式と挙動が異なります。

---

## ✅ 1. 摩擦力の分類（Types of Friction）

| 種類           | 記号           | 状態           | 式                         |
|----------------|----------------|----------------|----------------------------|
| 静止摩擦力     | \( f_s \)       | 物体が静止している | \( f_s \leq \mu_s N \)      |
| 最大静止摩擦力 | \( f_s^{\max} \) | 動き出す直前     | \( f_s^{\max} = \mu_s N \)  |
| 動摩擦力       | \( f_k \)       | 物体が滑っている   | \( f_k = \mu_k N \)         |

- \( \mu_s \)：静止摩擦係数（static friction coefficient）  
- \( \mu_k \)：動摩擦係数（kinetic friction coefficient）  
- \( N \)：垂直抗力（法線方向に働く力）

---

## ✅ 2. 静止摩擦の特徴

- 動かない限り、外力に**ちょうど等しい力で打ち消す**  
  → 式は「不等式」  
- 物体を動かすには：
  \[
  F_{\text{外}} > \mu_s N
  \]

---

## ✅ 3. 動摩擦の特徴

- 一度**滑り始めると**、摩擦力は一定の大きさになる
- 常に：
  \[
  f_k = \mu_k N \quad (\mu_k < \mu_s)
  \]

---

## ✅ 4. 斜面上の摩擦の例（静止 or 加速）

斜面角：\( \theta \)、物体質量：\( m \)

- **滑り出す条件（静止→動き出す）**：
  \[
  \mu_s \leq \tan\theta
  \]

- **動摩擦による運動方程式**（下向き加速）：
  \[
  ma = mg \sin\theta - \mu_k mg \cos\theta
  \quad \Rightarrow \quad
  a = g (\sin\theta - \mu_k \cos\theta)
  \]

---

## ✅ 5. 摩擦が関わる物理現象

| 応用例           | 関連する式や概念                 |
|------------------|------------------------------|
| 滑り出し条件     | \( \mu_s = \tan\theta_c \)（臨界角） |
| 台上の物体と外力 | \( F_{\text{必要}} = \mu_s N \) |
| 加速中の物体     | 運動方程式に \( -\mu_k N \) を追加  |
| 転がり摩擦       | 別モデル（\( f_r = \mu_r N \)）     |

---

## ✅ ポイントまとめ

| ポイント                 | 内容                                     |
|--------------------------|------------------------------------------|
| 最大静止摩擦力は変化する | 外力に応じて変動し、限界は \( \mu_s N \)     |
| 動摩擦力は一定           | 滑り始めたら \( f_k = \mu_k N \)          |
| 静止 > 動                | 常に \( \mu_s > \mu_k \) が成り立つ         |
| 摩擦力は常に逆向き       | 接触面に沿って運動や運動傾向に**逆らう方向**に働く |

---


```python
# Program Name: friction_simulation_slope.py
# Creation Date: 20250615
# Overview: Simulate motion on an inclined plane with static and kinetic friction
# Usage: Set angle, mu_s, mu_k and run to compute motion condition and acceleration

# ------------------------
# 必要なライブラリのインストール / Install required libraries
# ------------------------
!pip install numpy matplotlib

import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# 一元管理されるパラメータ / Parameter Definitions
# ------------------------
m = 1.0           # 質量 [kg]
g = 9.8           # 重力加速度 [m/s^2]
theta_deg = 30.0  # 斜面の傾斜角度 [degrees]
mu_s = 0.5        # 静止摩擦係数
mu_k = 0.4        # 動摩擦係数

# ------------------------
# ラジアン変換 / Convert angle to radians
# ------------------------
theta = np.radians(theta_deg)

# ------------------------
# 力成分の計算 / Compute forces
# ------------------------
Fg_parallel = m * g * np.sin(theta)     # 重力の斜面方向成分
N = m * g * np.cos(theta)               # 垂直抗力
f_s_max = mu_s * N                      # 最大静止摩擦力

# ------------------------
# 静止 or 加速の判定 / Determine motion status
# ------------------------
if Fg_parallel <= f_s_max:
    # 静止中 / At rest
    print(f"物体は静止しています（Fg = {Fg_parallel:.2f} N <= f_s_max = {f_s_max:.2f} N）")
    a = 0.0
else:
    # 動き出す / Starts sliding
    f_k = mu_k * N
    a = (Fg_parallel - f_k) / m
    print(f"物体は動き始めます（Fg = {Fg_parallel:.2f} N > f_s_max = {f_s_max:.2f} N）")
    print(f"動摩擦力 f_k = {f_k:.2f} N")
    print(f"加速度 a = {a:.2f} m/s²")

# ------------------------
# グラフによる力の比較 / Plot force components
# ------------------------
plt.figure(figsize=(8, 5))
plt.bar(['Fg_parallel', 'f_s_max', 'f_k'], [Fg_parallel, f_s_max, mu_k * N],
        color=['red', 'blue', 'green'])
plt.title("Force Comparison on Inclined Plane")
plt.ylabel("Force [N]")
plt.grid(axis='y')
plt.show()

```
---

## ✅ ばねの基本式と定数（Hookeの法則）

ばね（スプリング）の力学を記述する基本式は**Hooke（フック）の法則**で表されます。

### 🔸 フックの法則（Hooke's Law）

$$
F = -k x
$$

| 項目  | 内容                         |
| --- | -------------------------- |
| $F$ | ばねの復元力 \[N]（マイナスは変位の方向と反対） |
| $k$ | ばね定数（弾性係数）\[N/m]           |
| $x$ | 自然長からの変位（伸び or 縮み）\[m]     |

---

### ✅ ポイント解説

* **ばね定数 $k$** は「ばねの硬さ」を表し、大きいほど変位あたりの力が大きくなる。
* 符号「マイナス」は「変位方向とは逆向きの力」を意味する（復元力）。

---

### ✅ 合成ばね定数（直列・並列）

#### ● 直列接続（伸びが加算される）

$$
\frac{1}{k_{\text{eq}}} = \frac{1}{k_1} + \frac{1}{k_2} + \cdots
$$

#### ● 並列接続（力が加算される）

$$
k_{\text{eq}} = k_1 + k_2 + \cdots
$$

---

### ✅ 力学エネルギーとしての表現

$$
U = \frac{1}{2} k x^2
$$

* ばねの**位置エネルギー（弾性ポテンシャルエネルギー）**
* 変位が大きいほど、蓄えられるエネルギーも大きくなる（2乗に比例）

---


```python
# Program Name: spring_force_simulation.py
# Creation Date: 20250615
# Overview: Simulate spring force and elastic potential energy using Hooke's Law
# Usage: Run script to compute and plot spring force and potential energy

# ------------------------
# ライブラリインストール / Install required libraries
# ------------------------
!pip install matplotlib numpy

# ------------------------
# ライブラリ読み込み / Import libraries
# ------------------------
import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
k = 50         # ばね定数 [N/m] / spring constant
x_min = -0.2   # 最小変位 [m] / minimum displacement
x_max = 0.2    # 最大変位 [m] / maximum displacement
num_points = 200  # サンプル数 / number of data points

# ------------------------
# データ生成 / Generate data
# ------------------------
x = np.linspace(x_min, x_max, num_points)       # 変位 / displacement
F = -k * x                                      # フックの法則による力 / force from Hooke's law
U = 0.5 * k * x**2                              # 弾性位置エネルギー / elastic potential energy

# ------------------------
# グラフ描画 / Plot graphs
# ------------------------
# 力のグラフ / Force vs Displacement
plt.figure()
plt.plot(x, F, label='Spring Force F = -kx')
plt.title('Spring Force vs Displacement')
plt.xlabel('Displacement x [m]')
plt.ylabel('Force F [N]')
plt.grid(True)
plt.legend()
plt.show()

# エネルギーのグラフ / Potential Energy
plt.figure()
plt.plot(x, U, label='Elastic Potential Energy U = 1/2 kx²')
plt.title('Potential Energy vs Displacement')
plt.xlabel('Displacement x [m]')
plt.ylabel('Energy U [J]')
plt.grid(True)
plt.legend()
plt.show()

# ------------------------
# 数値出力 / Print sample data
# ------------------------
print("Sample displacement x [m]:", x[:5])
print("Corresponding force F [N]:", F[:5])
print("Corresponding energy U [J]:", U[:5])

```
# ✅ 浮力の基本まとめ（Archimedes' Principle）

浮力（buoyant force）は、**液体や気体中にある物体が受ける上向きの力**であり、**排除された流体の重さ**に等しいとされます。これを**アルキメデスの原理**に基づいて定式化できます。

---

## ✅ 基本式

\[
F_b = \rho_{\text{fluid}} \cdot g \cdot V
\]

| 記号               | 意味                                 | 単位       |
|------------------|------------------------------------|----------|
| \( F_b \)         | 浮力（Buoyant Force）               | [N]      |
| \( \rho_{\text{fluid}} \) | 流体の密度（Density of fluid）        | [kg/m³]  |
| \( g \)           | 重力加速度（Gravitational acceleration） | [m/s²]   |
| \( V \)           | 物体が排除した流体の体積（Volume）       | [m³]     |

---

## ✅ 力のつり合い（浮く・沈むの判定）

- **浮く場合**：\( F_b > mg \)（浮力が重力を上回る）
- **沈む場合**：\( F_b < mg \)（浮力が重力に負ける）
- **浮かんで停止する場合**：\( F_b = mg \)



```python

```
# ✅ 力積と運動量の関係（Impulse-Momentum Theorem）

力積とは、**ある時間の間に物体が受けた力の総効果**を表す物理量で、  
それは運動量の変化量に等しいという原理です。

---

## ✅ 基本式

\[
\vec{F} \Delta t = m\vec{v}' - m\vec{v}
\]

| 項目               | 意味                             | 単位           |
|--------------------|----------------------------------|----------------|
| \( \vec{F} \Delta t \) | 力積（Impulse）                 | [N·s] または [kg·m/s] |
| \( m \vec{v} \)        | 初期の運動量（Before）             | [kg·m/s]       |
| \( m \vec{v}' \)       | 最終の運動量（After）              | [kg·m/s]       |
| \( m \vec{v}' - m \vec{v} \) | 運動量の変化（Change of momentum） | [kg·m/s]       |

---

## ✅ ポイント解説

- 力 \( \vec{F} \) が一定の場合、その作用時間 \( \Delta t \) にわたる積をとった量が「力積」
- 力積は、**物体の運動量の変化と等しくなる**
- 単位系では：
  \[
  [\text{N·s}] = [\text{kg·m/s}]
  \]
  （実は**次元が同じ**）

---

## ✅ 応用例

- 野球のバット → ボールを打ち返す瞬間（大きな力を短時間で）
- エアバッグ → 衝突時に **時間を長くすることで力を減らす**
- ロケットの推進 → 力積の連続的適用による速度変化

---

```python
# Program Name: impulse_momentum_example.py
# Creation Date: 20250615
# Overview: Compute impulse and momentum change from initial and final velocities
# Usage: Run to see impulse from velocity change

# ------------------------
# ライブラリのインポート / Import libraries
# ------------------------
# 標準ライブラリのみ使用（外部依存なし）

# ------------------------
# パラメータの設定 / Define parameters
# ------------------------
m = 0.15       # 質量 [kg]（例：野球のボール）
v_before = -20  # 初速度 [m/s]（マイナス：向かってくる）
v_after = 40    # 最終速度 [m/s]（プラス：打ち返す）
delta_t = 0.01  # 接触時間 [s]

# ------------------------
# 力積と平均力の計算 / Compute impulse and average force
# ------------------------
impulse = m * (v_after - v_before)  # 力積 [N·s]
F_avg = impulse / delta_t           # 平均力 [N]

# ------------------------
# 結果出力 / Print results
# ------------------------
print("=== 力積と運動量の変化 ===")
print(f"質量 m = {m} kg")
print(f"初速度 v = {v_before} m/s")
print(f"最終速度 v' = {v_after} m/s")
print(f"接触時間 Δt = {delta_t} s")
print(f"→ 運動量の変化 = 力積 = {impulse:.3f} N·s")
print(f"→ 平均力 = {F_avg:.1f} N")

```

---

# ✅ 衝突の2次元運動量保存：行列形式とクラメルの公式

## 🎯 問題概要

* 質量 $m$ の物体 $P$ が角度 $\theta$ で速度 $v_0$ をもって飛来。
* 衝突後、$P$ は速度 $v$, 角度 $\phi$、質量 $M$ の物体は速度 $V$ で直線運動。
* 運動量は **x方向とy方向で保存される**。

---

## ✅ 成分ごとの保存式（連立一次方程式）

### ● x方向：

$$
m v_0 \cos\theta = m v \cos\phi + M V \quad \text{……(1)}
$$

### ● y方向：

$$
m v_0 \sin\theta = m v \sin\phi \quad \text{……(2)}
$$

---

## ✅ 行列形式に変換（未知数：$v, V$）

$$
\begin{bmatrix}
m \cos\phi & M \\
m \sin\phi & 0
\end{bmatrix}
\begin{bmatrix}
v \\
V
\end{bmatrix}
=
\begin{bmatrix}
m v_0 \cos\theta \\
m v_0 \sin\theta
\end{bmatrix}
$$

この形： $A \vec{x} = \vec{b}$ に対して、クラメルの公式で解ける。

---

## ✅ クラメルの公式による一般解

係数行列 $A$：

$$
A =
\begin{bmatrix}
m \cos\phi & M \\
m \sin\phi & 0
\end{bmatrix}
\quad , \quad
\vec{b} =
\begin{bmatrix}
m v_0 \cos\theta \\
m v_0 \sin\theta
\end{bmatrix}
$$

行列式：

$$
\det A = - M m \sin\phi
$$

---

### 解：$v$（Pの衝突後の速度）

$$
v =
\frac{
\begin{vmatrix}
m v_0 \cos\theta & M \\
m v_0 \sin\theta & 0
\end{vmatrix}
}{
\det A
}
=
\frac{
- M m v_0 \sin\theta
}{
- M m \sin\phi
}
=
\frac{v_0 \sin\theta}{\sin\phi}
$$

---

### 解：$V$（Mの衝突後の速度）

$$
V =
\frac{
\begin{vmatrix}
m \cos\phi & m v_0 \cos\theta \\
m \sin\phi & m v_0 \sin\theta
\end{vmatrix}
}{
\det A
}
=
\frac{
m^2 v_0 (\cos\phi \sin\theta - \sin\phi \cos\theta)
}{
- M m \sin\phi
}
$$

三角恒等式を使うと：

$$
V =
\frac{-m v_0 \sin(\theta - \phi)}{M \sin\phi}
$$

---

## ✅ 結論：一般化されたベクトル解（クラメルの公式による）

$$
\boxed{
\begin{aligned}
v &= \frac{v_0 \sin\theta}{\sin\phi} \\
V &= \frac{-m v_0 \sin(\theta - \phi)}{M \sin\phi}
\end{aligned}
}
$$

---

## ✅ ポイントまとめ

| 項目      | 内容                         |
| ------- | -------------------------- |
| 行列形式    | 2変数 $v, V$ に対する 2元一次方程式で構成 |
| 解法      | クラメルの公式を使って厳密に導出可能         |
| 三角関数の扱い | $\sin, \cos$ を使うため角度処理が必要  |
| 解の構造    | 運動量保存 × 方向の分解 × ベクトル計算に対応  |

---


```python
# Program Name: collision_cramers_solver.py
# Creation Date: 20250615
# Overview: Solve a 2D momentum conservation problem using Cramer's Rule
# Usage: Set m, M, v0, theta, phi and run to compute v and V

# ------------------------
# 必要ライブラリのインストール / Install required libraries
# ------------------------
!pip install numpy

# ------------------------
# ライブラリ読み込み / Import libraries
# ------------------------
import numpy as np

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
m = 1.0             # 衝突する物体の質量 [kg]
M = 2.0             # 相手物体の質量 [kg]
v0 = 10.0           # 初速度 [m/s]
theta_deg = 30      # 入射角 [deg]
phi_deg = 60        # 反射角 [deg]

# ラジアンに変換 / Convert to radians
theta = np.radians(theta_deg)
phi = np.radians(phi_deg)

# ------------------------
# 行列定義 / Define coefficient matrix A and vector b
# ------------------------
A = np.array([
    [m * np.cos(phi), M],
    [m * np.sin(phi), 0]
])

b = np.array([
    m * v0 * np.cos(theta),
    m * v0 * np.sin(theta)
])

# ------------------------
# クラメルの公式で解く / Solve using Cramer's Rule
# ------------------------
detA = np.linalg.det(A)

A_v = A.copy()
A_v[:, 0] = b
detA_v = np.linalg.det(A_v)

A_V = A.copy()
A_V[:, 1] = b
detA_V = np.linalg.det(A_V)

v = detA_v / detA
V = detA_V / detA

# ------------------------
# 結果出力 / Print results
# ------------------------
print("=== 衝突後の速度（Cramer's Rule） ===")
print(f"入力条件：θ = {theta_deg}°, φ = {phi_deg}°, 初速度 v0 = {v0} m/s")
print(f"衝突後の物体Pの速度 v = {v:.3f} m/s")
print(f"衝突後の物体Mの速度 V = {V:.3f} m/s")

```


# ✅ 衝突後の速度公式（1次元・反発係数付き）

## ⚙️ 基本設定

- 質量 \( m \) の物体 \( P \) が速度 \( v_0 \) で右向きに運動
- 質量 \( M \) の物体 \( Q \) は静止
- 衝突後の速度をそれぞれ \( v_P, v_Q \)
- 反発係数 \( e \) を考慮する

---

## ✅ 運動量保存の式（①）

\[
m v_P + M v_Q = m v_0 \tag{1}
\]

---

## ✅ 反発係数の定義（②）

\[
v_P - v_Q = -e (v_0 - 0) = -e v_0 \tag{2}
\]

---

## ✅ ① + \( M \) × ② により \( v_Q \) を消去して \( v_P \) を求める

\[
(m + M) v_P = (m - e M) v_0
\Rightarrow
\boxed{
v_P = \frac{m - e M}{m + M} v_0 \tag{3}
}
\]

---

## ✅ ① - \( m \) × ② により \( v_P \) を消去して \( v_Q \) を求める

\[
(m + M) v_Q = (1 + e) m v_0
\Rightarrow
\boxed{
v_Q = \frac{(1 + e) m}{m + M} v_0 \tag{4}
}
\]

---

## ✅ Pがはね返る条件

- \( v_P < 0 \) となる必要がある
- 上記 (3) より：

\[
\frac{m - e M}{m + M} < 0
\Rightarrow
\boxed{
m < e M
}
\]

---

## ✅ まとめ

| 項目 | 結果式 |
|------|--------|
| 衝突後の \( P \) の速度 | \( v_P = \frac{m - e M}{m + M} v_0 \) |
| 衝突後の \( Q \) の速度 | \( v_Q = \frac{(1 + e) m}{m + M} v_0 \) |
| Pがはね返る条件         | \( m < e M \)                        |

---



```python
# Program Name: elastic_collision_1d.py
# Creation Date: 20250615
# Overview: Compute 1D collision outcome using momentum and restitution coefficient
# Usage: Set m, M, v0, e and run to compute v_P and v_Q

# ------------------------
# ライブラリのインポート / Import libraries
# ------------------------
# 標準ライブラリのみ使用（追加インストール不要）

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
m = 1.0        # Pの質量 [kg]
M = 3.0        # Qの質量 [kg]
v0 = 10.0      # Pの初速度 [m/s]
e = 0.8        # 反発係数（0 ≤ e ≤ 1）

# ------------------------
# 式に基づく速度計算 / Calculate post-collision velocities
# ------------------------
v_P = ((m - e * M) / (m + M)) * v0
v_Q = ((1 + e) * m / (m + M)) * v0

# ------------------------
# 出力 / Output results
# ------------------------
print("=== 1次元衝突（反発係数付き）の結果 ===")
print(f"Pの質量 m = {m} kg, Qの質量 M = {M} kg")
print(f"初速度 v0 = {v0} m/s, 反発係数 e = {e}")
print(f"衝突後のPの速度 v_P = {v_P:.3f} m/s")
print(f"衝突後のQの速度 v_Q = {v_Q:.3f} m/s")

# ------------------------
# はね返る条件の確認 / Check rebound condition
# ------------------------
if v_P < 0:
    print(f"→ Pははね返ります（v_P < 0）: OK（条件: m < eM）")
else:
    print(f"→ Pははね返りません（v_P ≥ 0）: NG（条件: m < eM）")

```


# ✅ 完全弾性衝突の速度公式（運動量保存＋エネルギー保存）

---

## ⚙️ 条件設定（物体 P と Q の衝突）

- 質量 \( m \) の物体 P が速度 \( v_0 \) で運動し、
- 静止していた質量 \( M \) の物体 Q に衝突
- 衝突後の速度：P の速度 \( u \)、Q の速度 \( U \)
- **完全弾性衝突**なので「運動量保存」と「力学的エネルギー保存」がともに成り立つ

---

## ✅ ① 運動量保存の式

\[
mv_0 = mu + MU \tag{1}
\]

---

## ✅ ② 力学的エネルギー保存の式

\[
\frac{1}{2}mv_0^2 = \frac{1}{2}mu^2 + \frac{1}{2}MU^2 \tag{2}
\]

---

## ✅ 式(1)より \( U \) を \( u \) で表す：

\[
U = \frac{m(v_0 - u)}{M}
\]

---

## ✅ これを式(2)に代入して整理：

\[
(m + M)u^2 - 2mv_0 u + (m - M)v_0^2 = 0
\]

---

## ✅ 2次方程式として解くと：

\[
u = \frac{m \pm M}{m + M} v_0
\]

うち、  
- \( u = v_0 \) を選ぶと \( U = 0 \) となり衝突にならない（不適）  
- よって、**適切な根**は：

\[
\boxed{u = \frac{m - M}{m + M} v_0}
\quad , \quad
\boxed{U = \frac{2m}{m + M} v_0}
\]

---

## ✅ 結論（完全弾性衝突の結果）

| 物体 | 衝突後の速度 |
|------|--------------|
| P    | \( \displaystyle u = \frac{m - M}{m + M} v_0 \) |
| Q    | \( \displaystyle U = \frac{2m}{m + M} v_0 \)     |

---

## ✅ 特殊ケース

- \( m = M \)：Pは静止、Qは \( v_0 \) で飛び出す（完全交換）
- \( M \gg m \)：Pはほぼ逆向き \( -v_0 \)、Qは動かない

```python
# Program Name: elastic_collision_energy_momentum.py
# Creation Date: 20250615
# Overview: Solve a perfectly elastic 1D collision using momentum and energy conservation
# Usage: Set m, M, v0 and run to compute post-collision velocities u (P), U (Q)

# ------------------------
# ライブラリのインポート / Import libraries
# ------------------------
# 標準ライブラリのみで実行可能

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
m = 1.0       # Pの質量 [kg]
M = 2.0       # Qの質量 [kg]
v0 = 10.0     # Pの初速度 [m/s]、Qは静止していた

# ------------------------
# 弾性衝突における衝突後の速度 / Post-collision velocities
# 運動量保存 + 力学的エネルギー保存を同時に満たす解析解
# ------------------------
u = ((m - M) / (m + M)) * v0   # Pの衝突後の速度
U = (2 * m / (m + M)) * v0     # Qの衝突後の速度

# ------------------------
# 結果表示 / Print results
# ------------------------
print("=== 完全弾性衝突（P→Q）の計算結果 ===")
print(f"Pの初速度: v0 = {v0} m/s")
print(f"Pの衝突後の速度: u = {u:.3f} m/s")
print(f"Qの衝突後の速度: U = {U:.3f} m/s")

# ------------------------
# 検算: 運動量とエネルギーの保存 / Validation
# ------------------------
momentum_before = m * v0
momentum_after = m * u + M * U
energy_before = 0.5 * m * v0**2
energy_after = 0.5 * m * u**2 + 0.5 * M * U**2

print("\n--- 検算 ---")
print(f"運動量保存: {momentum_before:.3f} ≈ {momentum_after:.3f}")
print(f"エネルギー保存: {energy_before:.3f} ≈ {energy_after:.3f}")

```
# ✅ 運動量保存と重心の運動（2物体系）

---

## ✅ 各変数の定義

- 質量 \( m \)：物体 \( P \) の質量
- 質量 \( M \)：物体 \( Q \) の質量
- 座標 \( x_P, x_Q \)：それぞれの位置
- 速度 \( v_P = \frac{dx_P}{dt},\quad v_Q = \frac{dx_Q}{dt} \)

---

## ✅ 重心の位置 \( x_G \)

\[
x_G = \frac{m x_P + M x_Q}{m + M}
\]

---

## ✅ 重心の速度 \( v_G \)

座標と速度の関係（導関数）を用いる：

\[
v_G = \frac{d x_G}{dt}
= \frac{m \frac{dx_P}{dt} + M \frac{dx_Q}{dt}}{m + M}
= \frac{m v_P + M v_Q}{m + M}
\]

---

## ✅ ポイント：重心の速度と運動量保存の対応

全運動量：

\[
\text{Total momentum} = m v_P + M v_Q
\]

したがって、

\[
v_G = \frac{\text{Total momentum}}{m + M}
\quad \text{（重心速度は全運動量の "比率"）}
\]

---

## ✅ 特別な場合：全体が初め静止していたとき

- 初期運動量 \( m v_P + M v_Q = 0 \) のとき、
- 重心の速度：

\[
v_G = \frac{0}{m + M} = 0
\]

よって、  
\[
\boxed{v_G = 0 \Rightarrow x_G = \text{一定（重心は動かない）}}
\]

---

## ✅ まとめ

| 観点               | 結論                                     |
|--------------------|--------------------------------------------|
| 重心の速度         | \( v_G = \frac{m v_P + M v_Q}{m + M} \)      |
| 重心と運動量の関係 | \( v_G = \frac{\text{全運動量}}{\text{合計質量}} \) |
| 初期静止の場合     | \( v_G = 0 \Rightarrow \) 重心は不動         |

---


```python
# Program Name: center_of_mass_motion.py
# Creation Date: 20250615
# Overview: Compute center of mass velocity from two-body motion using momentum conservation
# Usage: Set masses and velocities to compute center of mass velocity

# ------------------------
# ライブラリのインポート / Import standard libraries
# ------------------------
# 追加ライブラリ不要

# ------------------------
# パラメータ定義 / Define masses and velocities
# ------------------------
m = 2.0        # Pの質量 [kg]
M = 3.0        # Qの質量 [kg]
v_P = 5.0      # Pの速度 [m/s]
v_Q = -2.0     # Qの速度 [m/s]

# ------------------------
# 重心速度の計算 / Compute center of mass velocity
# ------------------------
v_G = (m * v_P + M * v_Q) / (m + M)

# ------------------------
# 結果表示 / Print results
# ------------------------
print("=== 重心の速度の計算 ===")
print(f"Pの質量: m = {m} kg, 速度 v_P = {v_P} m/s")
print(f"Qの質量: M = {M} kg, 速度 v_Q = {v_Q} m/s")
print(f"重心の速度 v_G = {v_G:.3f} m/s")

# ------------------------
# 初期全体が静止していたか確認 / Check if center of mass is stationary
# ------------------------
if abs(v_G) < 1e-6:
    print("→ 重心は静止している（v_G ≈ 0）")
else:
    print("→ 重心は一定速度で運動している（v_G ≠ 0）")

```

# ✅ 円運動の基本式まとめ

円運動では、物体が一定の半径 \( r \) を保ちながら中心のまわりを回転します。
このとき、**速度・加速度・周期**には以下の関係が成り立ちます。

---

## 🔸 1. 接線速度（tangential speed）

\[
v = r \omega
\]

- \( r \)：円の半径 [m]  
- \( \omega \)：角速度 [rad/s]  
- 接線方向の速さ（運動の方向）

---

## 🔸 2. 向心加速度（centripetal acceleration）

\[
a = r \omega^2 = \frac{v^2}{r}
\]

- 円の中心向きに働く加速度
- 半径 \( r \) が小さいほど加速度が大きくなる
- 速度 \( v \) が大きいほど曲がりにくい（加速度増）

---

## 🔸 3. 周期（1周にかかる時間）

\[
T = \frac{2\pi}{\omega}
\]

- \( T \)：1回転の所要時間 [s]
- \( \omega \)：角速度 [rad/s]
- 角速度と周期は反比例

---

## ✅ 補足：角速度の定義

\[
\omega = \frac{2\pi}{T}
\quad \Rightarrow \quad
v = \frac{2\pi r}{T}
\]

---

## ✅ 単位まとめ

| 量         | 記号 | 単位     |
|------------|------|----------|
| 半径       | \( r \) | [m]     |
| 角速度     | \( \omega \) | [rad/s] |
| 接線速度   | \( v \) | [m/s]   |
| 加速度     | \( a \) | [m/s²]  |
| 周期       | \( T \) | [s]     |


```python
# Program Name: uniform_circular_motion.py
# Creation Date: 20250615
# Overview: Calculate velocity, acceleration, and period in uniform circular motion
# Usage: Set radius and angular velocity, then run the script

# ------------------------
# ライブラリのインポート / Import libraries
# ------------------------
import math

# ------------------------
# パラメータ定義 / Define parameters
# ------------------------
r = 2.0               # 半径 [m]
omega = 4.0           # 角速度 ω [rad/s]

# ------------------------
# 円運動の基本式 / Compute fundamental quantities
# ------------------------
v = r * omega                        # 接線速度 [m/s]
a = r * omega**2                     # 向心加速度 [m/s²]
T = (2 * math.pi) / omega           # 周期 [s]

# ------------------------
# 結果出力 / Print results
# ------------------------
print("=== 等速円運動の計算結果 ===")
print(f"半径 r = {r} m")
print(f"角速度 ω = {omega} rad/s")
print(f"→ 接線速度 v = {v:.3f} m/s")
print(f"→ 向心加速度 a = {a:.3f} m/s²")
print(f"→ 周期 T = {T:.3f} s")

# ------------------------
# 補足確認 / Additional verification
# ------------------------
# v² / r でも加速度と一致するか確認
a_check = v**2 / r
print(f"(検算) v²/r = {a_check:.3f} m/s²")

```


# 🎯 単振動の一般解とその導関数（速度・加速度まで）

## ✅ 微分方程式

$$
\frac{d^2x}{dt^2} + \omega^2 x = 0
$$

* これは、バネや単振り子などに代表される\*\*単振動（Simple Harmonic Motion, SHM）\*\*を記述する基本的な方程式です。

---

## ✅ 一般解（初期条件付き）

初期条件：

* 初期変位：$x(0) = x_0$
* 初期速度：$v(0) = v_0$

そのときの一般解は以下のように表されます：

```math
x(t) = x_0 \cos(\omega t) + \frac{v_0}{\omega} \sin(\omega t)
```

ここで：

* $x_0$：初期位置（\[m]）
* $v_0$：初期速度（\[m/s]）
* $\omega$：角振動数（\[rad/s]）

---

## ✅ 速度の式（1階微分）

```math
v(t) = \frac{dx}{dt} = -x_0 \omega \sin(\omega t) + v_0 \cos(\omega t)
```

* 位置を時間微分したもの。
* $t = 0$ のとき、$v(0) = v_0$ が確認できます。

---

## ✅ 加速度の式（2階微分）

```math
a(t) = \frac{d^2x}{dt^2} = -x_0 \omega^2 \cos(\omega t) - v_0 \omega \sin(\omega t)
```

* 速度を時間微分したもの。
* $a(t) = -\omega^2 x(t)$ を満たすことが確認できます。

---



---

## ✅ 補足：物理的意味

| 量                       | 意味        | 特徴                    |
| ----------------------- | --------- | --------------------- |
| $x(t)$                  | 位置        | 振動の主成分。コサインとサインの合成で表現 |
| $v(t)$                  | 速度        | 位置の時間微分、波形は位相ずれ       |
| $a(t)$                  | 加速度       | 常に位置のマイナス方向に比例        |
| $a(t) = -\omega^2 x(t)$ | SHMの本質的特徴 | 中心方向への復元力を表す          |

---



```python
import numpy as np
import matplotlib.pyplot as plt

# ----------------------------
# Parameters / パラメータ定義
# ----------------------------
A = 2.0               # 初期変位 Amplitude [m]
v0 = 0.0              # 初期速度 Initial velocity [m/s]
omega = 2 * np.pi     # 角振動数 Angular frequency [rad/s]
t = np.linspace(0, 4, 1000)  # 時間軸 [0, 4秒]

# ----------------------------
# Motion equations / 運動方程式
# ----------------------------
x = A * np.cos(omega * t) + (v0 / omega) * np.sin(omega * t)        # 位置
v = -A * omega * np.sin(omega * t) + v0 * np.cos(omega * t)         # 速度
a = -A * omega**2 * np.cos(omega * t) - v0 * omega * np.sin(omega * t)  # 加速度

# 周期・周波数などの物理量
T = 2 * np.pi / omega              # 周期 [s]
f = 1 / T                          # 周波数 [Hz]
wavelength = None                 # 単振動には波長の概念は不要（記載しない）
angular_freq = omega              # 角振動数

# ----------------------------
# Plot / グラフ描画
# ----------------------------
plt.figure(figsize=(12, 8))

# 位置
plt.subplot(3, 1, 1)
plt.plot(t, x, label='Displacement x(t)')
plt.ylabel("x(t) [m]")
plt.grid()
plt.legend()

# 速度
plt.subplot(3, 1, 2)
plt.plot(t, v, label='Velocity v(t)', color='orange')
plt.ylabel("v(t) [m/s]")
plt.grid()
plt.legend()

# 加速度
plt.subplot(3, 1, 3)
plt.plot(t, a, label='Acceleration a(t)', color='green')
plt.xlabel("Time [s]")
plt.ylabel("a(t) [m/s²]")
plt.grid()
plt.legend()

plt.tight_layout()
plt.show()

# 周期・周波数などを表示
T, f, omega

```



# 🧭 摩擦のあるばね振り子：運動方程式とその解（減衰単振動モデル）

## 🔁 概要

* 質量 $m$
* ばね定数 $k$
* 動摩擦係数 $\mu$
* 水平ばね振り子（自然長から引っ張って離す）
* **摩擦が存在することで振動中心が移動し、振幅が減衰する**のが特徴

---

## ✅ 行きの運動（右向き）

### 🔹 合力：

$$
F = -kx + \mu mg
$$

### 🔹 座標変換：

$$
X = x - \frac{\mu mg}{k}
\quad \Rightarrow \quad
F = -kX
$$

### 🔹 解：

$$
X(t) = A \cos(\omega t)
\quad \Rightarrow \quad
x(t) = A \cos(\omega t) + \frac{\mu mg}{k}
$$

ここで、

* $\omega = \sqrt{\frac{k}{m}}$：角振動数
* $A$：初期振幅（点 A から C₁ までの距離）

---

## ✅ 戻りの運動（左向き）

### 🔹 合力：

$$
F = -kx - \mu mg
= -k \left(x + \frac{\mu mg}{k} \right)
$$

### 🔹 座標変換：

$$
X = x + \frac{\mu mg}{k}
\quad \Rightarrow \quad
F = -kX
$$

### 🔹 解：

$$
X(t) = A \cos(\omega t)
\quad \Rightarrow \quad
x(t) = A \cos(\omega t) - \frac{\mu mg}{k}
$$

---

## ⏳ 半周期ごとの振幅と中心の遷移

| ターン数 | 振動中心 $C$            | 振幅 $A$                          |
| ---- | ------------------- | ------------------------------- |
| 1    | $+\frac{\mu mg}{k}$ | $A_0$（初期振幅）                     |
| 2    | $-\frac{\mu mg}{k}$ | $A_1 = A_0 - 2\frac{\mu mg}{k}$ |
| 3    | $+\frac{\mu mg}{k}$ | $A_2 = A_1 - 2\frac{\mu mg}{k}$ |
| …    | …                   | …                               |

→ 毎回中心が左右に交互に入れ替わりながら、**振幅は $2\mu mg/k$ ずつ減衰**していく。

---

## 🛑 最終停止条件

物体 $P$ は**端で一瞬静止**し、そのときのばねの弾性力が摩擦力以下なら、それ以上動かず停止する。

### ✅ 弾性力と摩擦力の釣り合い：

$$
k A_{\text{stop}} \leq \mu mg
\quad \Rightarrow \quad
A_{\text{stop}} \leq \frac{\mu mg}{k}
$$

---


```python
import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# Parameters / パラメータ定義
# ------------------------
mu = 0.1          # 動摩擦係数 / Coefficient of kinetic friction
m = 1.0           # 質量 [kg]
g = 9.8           # 重力加速度 [m/s^2]
k = 10.0          # ばね定数 [N/m]
A0 = 0.2          # 初期振幅 [m]

# ------------------------
# 振幅の減衰量 / Amplitude decay per cycle
# ------------------------
dA = 2 * mu * m * g / k

# ------------------------
# 周期計算 / Period of one cycle
# ------------------------
omega = np.sqrt(k / m)
T = 2 * np.pi / omega

# ------------------------
# 計算ループ / Compute amplitude decay over cycles
# ------------------------
amplitudes = []
centers = []
center_sign = 1  # 正方向から開始

A = A0
while A > dA:
    amplitudes.append(A)
    centers.append(center_sign * mu * m * g / k)
    A -= dA
    center_sign *= -1

# ------------------------
# 時系列データの生成（ステップ単振動）
# ------------------------
t = []
x = []
time = 0.0
dt = 0.01

for i, A in enumerate(amplitudes):
    center = centers[i]
    t_half = np.arange(0, T/2, dt)
    x_half = A * np.cos(omega * t_half) + center
    x.extend(x_half)
    t.extend(time + t_half)
    time += T/2

# ------------------------
# プロット / Plotting
# ------------------------
plt.figure(figsize=(10, 5))
plt.plot(t, x, label="Damped SHM with friction")
plt.xlabel("Time [s]")
plt.ylabel("Displacement [m]")
plt.title("Damped Oscillation with Friction (Piecewise SHM)")
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

```


---

# 🎯 単振り子の運動方程式と周期の導出

## ✅ 1. 状況設定

* 長さ $\ell$ の糸の先におもりをつけて振らせる（振り子）
* 糸の長さは変化しない
* 質量：$m$、重力加速度：$g$
* 振れ角：$\theta$（ラジアン）

---

## ✅ 2. 円弧方向の復元力

振り子の運動は「弧方向（接線方向）」に沿った運動です。

重力の接線方向成分が：

$$
F = -mg \sin \theta
$$

これが復元力となって、振り子を元の位置に戻そうとする。

---

## ✅ 3. 小角近似（$\theta \ll 1$）

角度が小さいとき（$\theta$ が約 $10^\circ$ 未満）：

$$
\sin \theta \approx \theta \quad \text{（ラジアン）}
$$

したがって、

$$
F = -mg \theta
$$

---

## ✅ 4. 弧方向の変位 $x$ との関係

弧長（円周上の位置）は：

$$
x = \ell \theta \quad \Rightarrow \quad \theta = \frac{x}{\ell}
$$

これを上式に代入：

$$
F = -mg \cdot \frac{x}{\ell} = -\left(\frac{mg}{\ell}\right) x
$$

---

## ✅ 5. 単振動とみなせる！

これは単振動の形式：

$$
F = -Kx \quad \text{with} \quad K = \frac{mg}{\ell}
$$

したがって、**単振動**として扱える。

---

## ✅ 6. 周期の導出

単振動の周期公式：

$$
T = 2\pi \sqrt{\frac{m}{K}} = 2\pi \sqrt{\frac{m}{mg/\ell}} = 2\pi \sqrt{\frac{\ell}{g}}
$$

---

## ✅ 結論

### 🔹 周期（Period）

$$
\boxed{T = 2\pi \sqrt{\frac{\ell}{g}}}
$$

* $\ell$：糸の長さに比例
* $g$：重力加速度に反比例
* **質量には依存しない** → これが「振り子の等時性（等周期性）」

---

## ☑️ 補足：振幅によらず周期が一定

* この性質は「単振り子の等時性」と呼ばれ、**ガリレオ**によって発見されました
* ただし、**大きな振れ角では成り立たない**（非線形化する）

---



```python
import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# Parameters / パラメータ定義
# ------------------------
l = 1.0                # 長さ [m] / Length of pendulum
g = 9.8                # 重力加速度 [m/s^2] / Gravitational acceleration
theta0 = np.radians(10)  # 初期角度 [rad] (10度)
omega = np.sqrt(g / l)   # 角振動数 [rad/s]
T = 2 * np.pi * np.sqrt(l / g)  # 周期 [s]

# ------------------------
# Time settings / 時間軸設定
# ------------------------
t = np.linspace(0, 2*T, 1000)

# ------------------------
# Small angle approximation solution / 小角近似解
# ------------------------
theta = theta0 * np.cos(omega * t)
x = l * theta  # 弧長の近似位置

# ------------------------
# Plotting / グラフ描画
# ------------------------
plt.figure(figsize=(10, 5))
plt.plot(t, x, label="Pendulum displacement x(t)")
plt.axhline(0, color='gray', linestyle='--', linewidth=0.5)
plt.title("Simple Pendulum (Small Angle Approximation)")
plt.xlabel("Time [s]")
plt.ylabel("Displacement x(t) [m]")
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

# 出力: 周期 T
T

```


# ✅ 万有引力による円軌道運動の力学整理

---

## 🔹 1. 遠心力と万有引力の釣り合い（円運動の条件）

$$
\frac{mv^2}{r} = G \frac{Mm}{r^2}
\tag{1}
$$

両辺 $m$ で割り、速度 $v$ を求めると：

$$
\boxed{v = \sqrt{\frac{GM}{r}}}
\quad \text{（軌道半径 } r \text{ の円軌道速度）}
$$

---

## 🔹 2. 円軌道の周期（1周するのにかかる時間）

軌道の長さ $2\pi r$、速度 $v$ より：

$$
T = \frac{2\pi r}{v}
\quad \Rightarrow \quad
\boxed{T = 2\pi \sqrt{\frac{r^3}{GM}}}
\tag{2}
$$

---

## 🔹 3. 地表すれすれの人工衛星（第1宇宙速度）

地球の半径 $r = R$ のとき：

$$
\boxed{v_1 = \sqrt{\frac{GM}{R}}}
\quad \text{（第1宇宙速度）}
$$

---

## 🔹 4. 力学的エネルギーの合計

運動エネルギー：

$$
E_k = \frac{1}{2}mv^2 = \frac{1}{2}m \cdot \frac{GM}{r} = \frac{GMm}{2r}
$$

位置エネルギー（万有引力による）：

$$
E_p = - \frac{GMm}{r}
$$

したがって、合計エネルギー：

$$
\boxed{E = E_k + E_p = - \frac{GMm}{2r}}
\tag{3}
$$

---

## 🔹 5. 無限遠に達するためのエネルギー（脱出条件）

無限遠でのエネルギーを 0 とするためには：

$$
E + \Delta E = 0 \Rightarrow \boxed{\Delta E = \frac{GMm}{2r}}
\quad \text{（必要な追加エネルギー）}
$$

---

## ✅ まとめ表（クラメル式っぽく）

| 項目       | 数式                                | 意味           |
| -------- | --------------------------------- | ------------ |
| 軌道速度     | $v = \sqrt{\dfrac{GM}{r}}$        | 円運動の速度       |
| 周期       | $T = 2\pi \sqrt{\dfrac{r^3}{GM}}$ | 1周にかかる時間     |
| 第1宇宙速度   | $v_1 = \sqrt{\dfrac{GM}{R}}$      | 地表から落ちない最小速度 |
| 力学的エネルギー | $E = -\dfrac{GMm}{2r}$            | 結合エネルギー（束縛）  |
| 脱出エネルギー  | $\Delta E = \dfrac{GMm}{2r}$      | 解放に必要な仕事     |

---



```python
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# ------------------------
# 定数定義
# ------------------------
G = 6.67430e-11      # 万有引力定数 [m^3 kg^-1 s^-2]
M = 5.972e24         # 地球質量 [kg]
R = 6.371e6          # 地球半径 [m]
m = 1000             # 衛星の質量 [kg]

# ------------------------
# 半径リスト（地表〜高度50000km）
# ------------------------
r_list_km = np.linspace(R/1000, (R + 5e7)/1000, 100)  # km
r_list = r_list_km * 1000  # m

# ------------------------
# 各種物理量の計算
# ------------------------
v = np.sqrt(G * M / r_list)                     # 軌道速度
T = 2 * np.pi * np.sqrt(r_list**3 / (G * M))    # 軌道周期
E_k = 0.5 * m * v**2                             # 運動エネルギー
E_p = -G * M * m / r_list                        # 位置エネルギー
E_total = E_k + E_p                              # 合計エネルギー
E_escape = G * M * m / (2 * r_list)              # 脱出に必要な追加エネルギー

# ------------------------
# データフレーム作成
# ------------------------
df_orbit = pd.DataFrame({
    "r [km]": r_list_km,
    "v [m/s]": v,
    "T [s]": T,
    "E_k [J]": E_k,
    "E_p [J]": E_p,
    "E_total [J]": E_total,
    "E_escape [J]": E_escape
})

import ace_tools as tools; tools.display_dataframe_to_user(name="Orbital Mechanics Table", dataframe=df_orbit)

# ------------------------
# グラフ描画（軌道速度とエネルギー）
# ------------------------
plt.figure(figsize=(10, 5))
plt.plot(r_list_km, v, label='Orbital Speed v(r)')
plt.title("Orbital Speed vs Radius")
plt.xlabel("Orbital Radius [km]")
plt.ylabel("Speed [m/s]")
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

```


## ✅ ケプラーの法則・角運動量保存【大学物理まとめ】

| 項目               | 数式                                                                   | 意味                   | 対応保存則            |
| ---------------- | -------------------------------------------------------------------- | -------------------- | ---------------- |
| **第1法則（軌道の法則）**  | $r(\theta) = \frac{a(1 - e^2)}{1 + e \cos\theta}$                    | 惑星の軌道は太陽を焦点とする**楕円** | -                |
| **第2法則（面積速度一定）** | $\frac{dA}{dt} = \frac{1}{2} r^2 \frac{d\theta}{dt} = \text{const.}$ | 等しい時間に等しい面積を掃く       | **角運動量保存**       |
| **角運動量保存**       | $L = m r^2 \frac{d\theta}{dt} = \text{const.}$                       | 中心力運動の特徴             | **保存量**          |
| **面積速度と角運動量の関係** | $\frac{dA}{dt} = \frac{L}{2m}$                                       | 面積速度は角運動量の半分/mass    | 幾何と力学の一致         |
| **第3法則（調和の法則）**  | $T^2 = \frac{4\pi^2 a^3}{G M}$                                       | 公転周期² ∝ 軌道長半径³       | **エネルギー保存**から導出可 |

---

## ✅ 導出：ニュートン力学と角運動量保存

### ⬛ 面積速度一定（第2法則）

$$
\frac{dA}{dt} = \frac{1}{2} r^2 \frac{d\theta}{dt}
$$

中心力 ⇒ トルク 0 ⇒ 角運動量 $L = m r^2 \frac{d\theta}{dt}$ 保存 ⇒ 面積速度一定。

---

### ⬛ 角運動量保存

$$
L = m r^2 \frac{d\theta}{dt} = \text{const.}
$$

---

### ⬛ 調和の法則（第3法則の導出）

円軌道と仮定：

1. 万有引力：$\frac{G M m}{r^2} = m \frac{v^2}{r}$
2. 速度：$v = \sqrt{\frac{G M}{r}}$
3. 周期：$T = \frac{2\pi r}{v} \Rightarrow T^2 = \frac{4\pi^2 r^3}{G M}$

---

## ✅ 円軌道における角運動量と周期の関係

$$
L = mrv = m r \cdot \sqrt{\frac{G M}{r}} = m \sqrt{G M r}
$$

$$
T = \frac{2\pi r}{v} = 2\pi \sqrt{\frac{r^3}{G M}}
$$

---

## ✅ 拡張トピック（必要に応じて追加可）

| 軌道の種類 | 離心率 $e$        | 特徴             |
| ----- | -------------- | -------------- |
| 楕円軌道  | $0 \leq e < 1$ | ケプラーの法則が適用可能   |
| 放物線軌道 | $e = 1$        | 脱出速度に相当        |
| 双曲線軌道 | $e > 1$        | 恒星間飛行やフライバイに関連 |

---



```python
# Program Name: kepler_laws_plot.py
# Creation Date: 20250615
# Overview: Visualization of Kepler's Laws using orbital mechanics in Python
# Usage: Run to generate plots showing elliptical orbit and Kepler's 2nd law (equal area sweep)

# ------------------------
# 必要なライブラリをインストール / Install required libraries
# ------------------------
!pip install numpy matplotlib

# ------------------------
# 必要なライブラリをインポート / Import necessary libraries
# ------------------------
import numpy as np
import matplotlib.pyplot as plt

# ------------------------
# 一元管理パラメータ定義 / Centralized parameter definitions
# ------------------------
a = 1.0           # 長半径 / Semi-major axis (AU)
e = 0.6           # 離心率 / Eccentricity
G = 1.0           # 万有引力定数 (簡略化のため 1) / Gravitational constant
M = 1.0           # 中心天体質量 / Mass of the central body (solar mass)
num_points = 360  # 軌道上の点数 / Number of points in orbit

# ------------------------
# 楕円軌道の極座標表現による描画 / Plot elliptical orbit using polar coordinates
# ------------------------
theta = np.linspace(0, 2 * np.pi, num_points)
r = a * (1 - e**2) / (1 + e * np.cos(theta))  # 楕円軌道式 / Elliptical orbit formula

# 極座標→デカルト座標へ変換 / Convert polar to Cartesian coordinates
x = r * np.cos(theta)
y = r * np.sin(theta)

# ------------------------
# グラフ描画 / Plot the orbit
# ------------------------
plt.figure(figsize=(8, 6))
plt.plot(x, y, label="Elliptical Orbit")
plt.plot(-a * e, 0, 'yo', label="Sun (Focus)")  # 太陽の位置（焦点） / Position of the sun (focus)
plt.gca().set_aspect('equal')
plt.title("Kepler's First Law: Elliptical Orbit")
plt.xlabel("x (AU)")
plt.ylabel("y (AU)")
plt.legend()
plt.grid(True)
plt.show()

```

```python

```