# 重力場及電場
> 作者：王一哲
> 日期：2018/5/7
<br />
若在空間中有一個質量為 $M$ 且質量均勻分布的球體，若以 $M$ 的球心為原點，在空間中位置向量為 $M$ 處的重力場為
$$ \vec g = - \frac{GM}{r^2} \hat r $$
上式中的負號帶表重力場方向指向球心。若有多個球體，第 $i$ 個球體質量為 $M_i$、球心位置為 $r_i$，則重力場為
$$ \vec g = \sum \vec g_i = - G \sum \frac{M_i}{(\vec r - \vec r_i)^2} \frac{\vec r - \vec r_i}{|\vec r - \vec r_i|} $$
如果要在黑板上畫出各個位置的重力場強度及方向，這幾乎是不可能的任務，下圖是我畫出來的地球重力場示意圖
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/sCJld1d.png">
<div style="text-align:center">地球重力場示意圖</div>
<br />
在以下的課程中，我們想要在空間每隔一段距離取一個點，以箭頭的長度及方向來表示該點的重力場。因此我們在程式 18-1 中先練習使用 for 迴圈，在空間中每隔一段距離畫一個箭頭。在程式 18-2 中除了畫出箭頭還要計算此處的重力場，再更新箭頭的長度及方向。如果能夠成功畫出一個球體的重力場，在程式 18-3 中則更進一步地畫出兩個球體，甚至是多個球體建立的重力場。
<br />
在畫完重力場之後，我們可以用相同的方法畫出帶電球體在空間中建立的電場，其數學型式為
$$ \vec E = \sum \vec E_i = k \sum \frac{Q_i}{(\vec r - \vec r_i)^2} \frac{\vec r - \vec r_i}{|\vec r - \vec r_i|} $$
我們只要稍微修改程式 18-3 就能畫出兩個帶電球體在空間中建立的電場。
<br />
## 程式 18-1. 用 for 迴圈產生箭頭 （[取得程式碼](https://github.com/YiZheWangTw/VPythonTutorial/blob/master/18.%E9%87%8D%E5%8A%9B%E5%A0%B4%E5%8F%8A%E9%9B%BB%E5%A0%B4/18-1_arrows.py)） （[GlowScript 網站動畫連結](http://www.glowscript.org/#/user/yizhe/folder/Public/program/18-1arrows)）
```python=
"""
VPython教學: 18-1. 用 for 迴圈產生箭頭
日期: 2018/3/2
作者: 王一哲
"""
from vpython import *
N, L = 4, 20
scene = canvas(title="Arrows", width=600, height=600, x=0, y=0, background=color.black)
fields = []
for i in range(N):
for j in range(N):
for k in range(N):
fields.append(arrow(pos=vec(L/N*i - L/2, L/N*j - L/2, L/N*k - L/2),
axis=vec(4, 0, 0), radius=1, color=color.green))
```
<br />
這個程式與之前的程式相比十分簡短，只要將箭頭畫完就可以了，因此連產生動畫的部分都沒有。比較特別之處在於使用了 3 層的 for 迴圈，第 1 層迴圈的變數是 i ，用來產生位置的 x 坐標；第 2 層迴圈的變數是 j ，用來產生位置的 y 坐標；第 3 層迴圈的變數是 k ，用來產生位置的 z 坐標。由於使用了 range(N)、N = 4 來產生數值，也就是只會產生 0、1、2、3，所以畫箭頭的位置為 -L/2、-L/4、0、L/4。程式執行的成果如下圖。
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://i.imgur.com/1NnZkKr.png">
<div style="text-align:center">程式18-1畫面截圖</div>
<br />
## 程式 18-2. 重力場, 可以改成畫不同的星球 （[取得程式碼](https://github.com/YiZheWangTw/VPythonTutorial/blob/master/18.%E9%87%8D%E5%8A%9B%E5%A0%B4%E5%8F%8A%E9%9B%BB%E5%A0%B4/18-2_gravitational_field.py)） （[GlowScript 網站動畫連結](http://www.glowscript.org/#/user/yizhe/folder/Public/program/18-2gravitationalfield)）
```python=
"""
VPython教學: 18-2. 重力場, 可以改成畫不同的星球
Ver. 1: 2018/3/2
Ver. 2: 2019/9/14 箭頭的顏色隨著量值改變
作者: 王一哲
"""
from vpython import *
"""
1. 參數設定, 設定變數及初始值, 太陽及行星半徑、質量、遠日距、遠日點速率, 資料來源
"""
radius = {"Mercury": 2439700, "Venus": 6051800, "Earth": 6371000, "Mars": 3389500, "Sun": 696392000}
mass = {"Mercury": 0.33011E24, "Venus": 4.8675E24, "Earth": 5.9723E24, "Mars": 0.64171E24, "Sun": 1988500E24}
G = 6.67408E-11 # 重力常數
N = 5 # 將顯示的空間每邊切成 N 等份
"""
2. 產生行星類別, 回傳行星產生的重力場
"""
class Planet(sphere):
def g(self, pos):
return -G*self.m / mag2(pos-self.pos) * norm(pos-self.pos)
"""
3. 畫面設定
"""
# 產生動畫視窗
L = radius["Earth"]*4
scene = canvas(title="Gravitatinoal Field", width=600, height=600, x=0, y=0, background=color.black, range=L)
# 產生地球
earth = Planet(pos=vec(0, 0, 0), radius=radius["Earth"], m=mass["Earth"], texture=textures.earth)
# 計算畫箭頭的位置, 如果不在地球內則加到串列 locations 當中
locations = []
for i in range(N+1):
for j in range(N+1):
for k in range(N+1):
location = vec(L/N*i - L/2, L/N*j - L/2, L/N*k - L/2)
if mag(location - earth.pos) > earth.radius:
locations.append(location)
# 依序讀取串列 locations 的元素, 在對應的位置產生箭頭
fields = []
for location in locations:
fields.append(arrow(pos=location, axis=vec(0, 0, 0), color=color.green))
# 更新箭頭的長度及方向, 長度乘以 1E6 才能看見, 記錄重力場強度最大值, 量值接近最大值偏紅色, 量值接近 0 偏綠色
fmax = 0
for field in fields:
field.axis = earth.g(field.pos)*1E6
if field.axis.mag >= fmax: fmax = field.axis.mag
for field in fields:
field.color = vec(field.axis.mag/fmax, 1 - field.axis.mag/fmax, 0)
```
<br />
由於我們在之前的課程〈[行星運動](hhttps://hackmd.io/@yizhewang/B1hUDpmGm)〉當中已經蒐集了不少行星資料，並且在程式 12-4 中學到如何利用 class 產生自訂類別的物件，因此我們可以將程式 12-4 和 18-1 合併，寫出程式 18-2。主要的步驟為：
1. 參數設定
2. 自訂類別 Planet，並在當中新增自訂的方法(method) g，用來回傳此物件在輸入的位置 pos 處產生的重力場。
3. 產生動畫視窗、地球。
4. 計算畫箭頭的位置，如果不在地球內則加到串列 locations 當中。
5. 依序讀取串列 locations 的元素，在對應的位置產生箭頭。
6. 第 49 ~ 55 行，依序讀取箭頭的位置，計算此處的重力場，更新箭頭的長度及方向。由於重力場的量值與畫面的寬度相比太短，因此需要將重力場乘以 1 × 10<sup>6</sup> 才能在畫面上看到箭頭。將重力場量值當中的最大值記錄為 fmax，量值越大箭頭的顏色越接近紅色，量值越小箭頭的顏色越接近綠色。
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/sCJld1d.png">
<div style="text-align:center">程式18-2畫面截圖</div>
<br />
## 程式 18-3. 重力場, 可以改成畫不同的星球 （[取得程式碼](https://github.com/YiZheWangTw/VPythonTutorial/blob/master/18.%E9%87%8D%E5%8A%9B%E5%A0%B4%E5%8F%8A%E9%9B%BB%E5%A0%B4/18-3_gravitational_field_2.py)） （[GlowScript 網站動畫連結](http://www.glowscript.org/#/user/yizhe/folder/Public/program/18-3gravitationalfield2)）
```python=
"""
VPython教學: 18-3. 重力場, 2個星球, 自行調整距離
Ver. 1: 2018/3/2
Ver. 2: 2019/9/14 箭頭的顏色隨著量值改變
作者: 王一哲
"""
from vpython import *
"""
1. 參數設定, 設定變數及初始值, 太陽及行星半徑、質量、遠日距、遠日點速率, 資料來源
"""
radius = {"Mercury": 2439700, "Venus": 6051800, "Earth": 6371000, "Mars": 3389500, "Sun": 696392000}
mass = {"Mercury": 0.33011E24, "Venus": 4.8675E24, "Earth": 5.9723E24, "Mars": 0.64171E24, "Sun": 1988500E24}
G = 6.67408E-11 # 重力常數
N = 5 # 將顯示的空間每邊切成 N 等份
"""
2. 產生行星類別, 回傳行星產生的重力場
"""
class Planet(sphere):
def g(self, pos):
return -G*self.m / mag2(pos-self.pos) * norm(pos-self.pos)
"""
3. 畫面設定
"""
# 產生動畫視窗
d = 2E7
L = 2*d
scene = canvas(title="Gravitatinoal Field", width=600, height=600, x=0, y=0, background=color.black, range=L)
# 產生地球及火星
earth = Planet(pos=vec(-d/2, 0, 0), radius=radius["Earth"], m=mass["Earth"], texture=textures.earth)
mars = Planet(pos=vec(d/2, 0, 0), radius=radius["Mars"], m=mass["Mars"], color=color.red)
# 計算畫箭頭的位置, 如果不在地球或火星內則加到串列 locations 當中
locations = []
for i in range(N+1):
for j in range(N+1):
for k in range(N+1):
location = vec(L/N*i - L/2, L/N*j - L/2, L/N*k - L/2)
if mag(location-earth.pos) > earth.radius and mag(location-mars.pos) > mars.radius:
locations.append(location)
# 依序讀取串列 locations 的元素, 在對應的位置產生箭頭
fields = []
for location in locations:
fields.append(arrow(pos=location, axis=vec(0, 0, 0), color=color.green))
# 更新箭頭的長度及方向, 長度乘以 1E6 才能看見, 記錄重力場強度最大值, 量值接近最大值偏紅色, 量值接近 0 偏綠色
fmax = 0
for field in fields:
field.axis = (earth.g(field.pos) + mars.g(field.pos))*1E6
if field.axis.mag >= fmax: fmax = field.axis.mag
for field in fields:
field.color = vec(field.axis.mag/fmax, 1 - field.axis.mag/fmax, 0)
```
<br />
程式 18-3 與 18-2 幾乎一模一樣，只是多畫了火星。由於要產生兩個行星，在此應該就能看出自訂類別較為方便之處。需要特別注意一點，兩個行星之間的距離並未按照真實數據繪製，否則畫面中什麼東西都看不到。
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/ztcJizC.png">
<div style="text-align:center">程式18-3畫面截圖</div>
<br />
## 程式 18-4. 電場, 2個球狀帶電體 （[取得程式碼](https://github.com/YiZheWangTw/VPythonTutorial/blob/master/18.%E9%87%8D%E5%8A%9B%E5%A0%B4%E5%8F%8A%E9%9B%BB%E5%A0%B4/18-4_electric_field.py)） （[GlowScript 網站動畫連結](http://www.glowscript.org/#/user/yizhe/folder/Public/program/18-4electricfield)）
```python=
"""
VPython教學: 18-4. 電場, 2個球狀帶電體
Ver. 1: 2018/3/2
Ver. 2: 2019/9/14 箭頭的顏色隨著量值改變, 改為不用繼承的類別
作者: 王一哲
"""
from vpython import *
"""
1. 參數設定, 設定變數及初始值
"""
size = 1 # 帶電球體半徑
d = 10 # 帶電球體連心線距離
L = 1.8*d # 畫面寬度
q1, c1 = 5, color.blue # 帶電球體1電量、顏色
q2, c2 = -5, color.red # 帶電球體2電量、顏色
ke = 8.988E9 # 靜電力常數
N = 6 # 將顯示的空間每邊切成 N 等份
"""
2. 產生帶電球體類別, 回傳帶電球體產生的電場
"""
class Ball:
def __init__(self, pos, radius, color, charge):
self.pos = pos
self.radius = radius
self.color = color
self.charge = charge
self.ball = sphere(pos=vec(self.pos), radius=self.radius, charge=self.charge, color=self.color)
def electric(self, pos2):
return ke*self.charge / mag2(pos2-self.pos) * norm(pos2-self.pos)
"""
3. 畫面設定
"""
# 產生動畫視窗
scene = canvas(title="Electric Field", width=600, height=600, x=0, y=0, background=color.black, range=L)
# 產生帶電球體1、2
b1 = Ball(pos=vec(-d/2, 0, 0), radius=size, color=c1, charge=q1)
b2 = Ball(pos=vec(d/2, 0, 0), radius=size, color=c2, charge=q2)
# 計算畫箭頭的位置, 如果不在帶電球體內則加到串列 locations 當中
locations = []
for i in range(N+1):
for j in range(N+1):
for k in range(N+1):
location = vec(L/N*i - L/2, L/N*j - L/2, L/N*k - L/2)
if mag(location-b1.pos) > 2*size and mag(location-b2.pos) > 2*size:
locations.append(location)
# 依序讀取串列 locations 的元素, 在對應的位置產生箭頭
fields = []
for location in locations:
fields.append(arrow(pos=location, axis=vec(0, 0, 0), color=color.green))
# 更新箭頭的長度及方向, 長度乘以 1E-9 才不會太長, 記錄電場強度最大值, 量值接近最大值偏紅色, 量值接近 0 偏綠色
fmax = 0
for field in fields:
field.axis = (b1.electric(field.pos) + b2.electric(field.pos))*1E-9
if field.axis.mag >= fmax: fmax = field.axis.mag
for field in fields:
field.color = vec(field.axis.mag/fmax, 1 - field.axis.mag/fmax, 0)
```
<br />
程式 18-4 與 18-3 幾乎一模一樣，以下只解釋不同之處：
1. 將球體半徑、畫面寬度大幅縮小，設定帶電球體的電量、靜電力常數。
2. 第 23 ~ 31 行，自訂類別 Ball 時不採用繼承的方式，產生此類別的物件時要輸入位置 pos、半徑 radius、顏色 color、電量 charge，產生一個球體，並且自訂計算電場用的函式 electric。
3. 若左側小球電量為 q1、右側小球電量為 q2，以下是 3 種不同組合的模擬結果。
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/GFQ8aD7.png">
<div style="text-align:center">q1 = 1, q2 = -5 模擬結果畫面截圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/0d5oPc1.png">
<div style="text-align:center">q1 = 5, q2 = -5 模擬結果畫面截圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/tvDP67Q.png">
<div style="text-align:center">q1 = 5, q2 = -1 模擬結果畫面截圖</div>
<br />
## 結語
這是我們第一次嘗試用 VPython 畫出向量場的示意圖，看起來效果還不錯，之後應該可以用類似的方法畫出載流導線或線圈產生的磁場示意圖，請參考〈[電流的磁效應](https://hackmd.io/@yizhewang/r1g8P0uf7)〉。
<br />
## 太陽系天體資料來源
1. **太陽** https://nssdc.gsfc.nasa.gov/planetary/factsheet/sunfact.html
2. **水星** https://nssdc.gsfc.nasa.gov/planetary/factsheet/mercuryfact.html
3. **金星** https://nssdc.gsfc.nasa.gov/planetary/factsheet/venusfact.html
4. **地球** https://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html
5. **火星** https://nssdc.gsfc.nasa.gov/planetary/factsheet/marsfact.html
<br />
## VPython官方說明書
1. **canvas**: http://www.glowscript.org/docs/VPythonDocs/canvas.html
2. **sphere**: http://www.glowscript.org/docs/VPythonDocs/sphere.html
3. **arrow**: http://www.glowscript.org/docs/VPythonDocs/arrow.html
<br />
---
## 補充 for 迴圈的使用方法
> 2020/6/10：3 層 for 迴圈
> 2021/12/23：for 迴圈基本語法
### 1 層 for 迴圈
若程式碼為
```python
for i in range(5): print(i)
```
執行時會輸出
```python
0
1
2
3
4
```
<br />
range是用來產生數列的函式，格式為
```python
range(起始值, 結束值, 增量)
```
其中起始值預設為0，增量預設為1，執行時不包含結束值。
<br />
若程式碼為
```python
for i in range(1, 10, 2): print(i)
```
執行時會輸出
```python
1
3
5
7
9
```
<br />
若程式碼為
```python
for i in range(10, 1, -2): print(i)
```
執行時會輸出
```python
10
8
6
4
2
```
<br />
### 2 層 for 迴圈
若程式碼為
```python
for i in range(4):
for j in range(4):
print("i = {:d} j = {:d}".format(i, j))
```
執行時會輸出
```python
i = 0 j = 0
i = 0 j = 1
i = 0 j = 2
i = 0 j = 3
i = 1 j = 0
i = 1 j = 1
i = 1 j = 2
i = 1 j = 3
i = 2 j = 0
i = 2 j = 1
i = 2 j = 2
i = 2 j = 3
i = 3 j = 0
i = 3 j = 1
i = 3 j = 2
i = 3 j = 3
```
<br />
如果想要印出九九乘法表可以這樣寫
```python
for i in range(2, 10):
for j in range(1, 10):
print("{:d} * {:d} = {:d}".format(i, j, i*j), end="\t")
print()
```
執行時會輸出
```python
2 * 1 = 2 2 * 2 = 4 2 * 3 = 6 2 * 4 = 8 2 * 5 = 10 2 * 6 = 12 2 * 7 = 14 2 * 8 = 16 2 * 9 = 18
3 * 1 = 3 3 * 2 = 6 3 * 3 = 9 3 * 4 = 12 3 * 5 = 15 3 * 6 = 18 3 * 7 = 21 3 * 8 = 24 3 * 9 = 27
4 * 1 = 4 4 * 2 = 8 4 * 3 = 12 4 * 4 = 16 4 * 5 = 20 4 * 6 = 24 4 * 7 = 28 4 * 8 = 32 4 * 9 = 36
5 * 1 = 5 5 * 2 = 10 5 * 3 = 15 5 * 4 = 20 5 * 5 = 25 5 * 6 = 30 5 * 7 = 35 5 * 8 = 40 5 * 9 = 45
6 * 1 = 6 6 * 2 = 12 6 * 3 = 18 6 * 4 = 24 6 * 5 = 30 6 * 6 = 36 6 * 7 = 42 6 * 8 = 48 6 * 9 = 54
7 * 1 = 7 7 * 2 = 14 7 * 3 = 21 7 * 4 = 28 7 * 5 = 35 7 * 6 = 42 7 * 7 = 49 7 * 8 = 56 7 * 9 = 63
8 * 1 = 8 8 * 2 = 16 8 * 3 = 24 8 * 4 = 32 8 * 5 = 40 8 * 6 = 48 8 * 7 = 56 8 * 8 = 64 8 * 9 = 72
9 * 1 = 9 9 * 2 = 18 9 * 3 = 27 9 * 4 = 36 9 * 5 = 45 9 * 6 = 54 9 * 7 = 63 9 * 8 = 72 9 * 9 = 81
```
<br />
### 3 層 for 迴圈
若程式碼為
```python=
N, counts = 4, 0
for i in range(N):
for j in range(N):
for k in range(N):
counts += 1
print("i = {:d}, j = {:d}, k = {:d}".format(i, j, k))
print("counts = {:d}".format(counts))
```
執行時會輸出
```python
i = 0, j = 0, k = 0
i = 0, j = 0, k = 1
i = 0, j = 0, k = 2
i = 0, j = 0, k = 3
i = 0, j = 1, k = 0
i = 0, j = 1, k = 1
i = 0, j = 1, k = 2
i = 0, j = 1, k = 3
i = 0, j = 2, k = 0
i = 0, j = 2, k = 1
i = 0, j = 2, k = 2
i = 0, j = 2, k = 3
i = 0, j = 3, k = 0
i = 0, j = 3, k = 1
i = 0, j = 3, k = 2
i = 0, j = 3, k = 3
i = 1, j = 0, k = 0
i = 1, j = 0, k = 1
i = 1, j = 0, k = 2
i = 1, j = 0, k = 3
i = 1, j = 1, k = 0
i = 1, j = 1, k = 1
i = 1, j = 1, k = 2
i = 1, j = 1, k = 3
i = 1, j = 2, k = 0
i = 1, j = 2, k = 1
i = 1, j = 2, k = 2
i = 1, j = 2, k = 3
i = 1, j = 3, k = 0
i = 1, j = 3, k = 1
i = 1, j = 3, k = 2
i = 1, j = 3, k = 3
i = 2, j = 0, k = 0
i = 2, j = 0, k = 1
i = 2, j = 0, k = 2
i = 2, j = 0, k = 3
i = 2, j = 1, k = 0
i = 2, j = 1, k = 1
i = 2, j = 1, k = 2
i = 2, j = 1, k = 3
i = 2, j = 2, k = 0
i = 2, j = 2, k = 1
i = 2, j = 2, k = 2
i = 2, j = 2, k = 3
i = 2, j = 3, k = 0
i = 2, j = 3, k = 1
i = 2, j = 3, k = 2
i = 2, j = 3, k = 3
i = 3, j = 0, k = 0
i = 3, j = 0, k = 1
i = 3, j = 0, k = 2
i = 3, j = 0, k = 3
i = 3, j = 1, k = 0
i = 3, j = 1, k = 1
i = 3, j = 1, k = 2
i = 3, j = 1, k = 3
i = 3, j = 2, k = 0
i = 3, j = 2, k = 1
i = 3, j = 2, k = 2
i = 3, j = 2, k = 3
i = 3, j = 3, k = 0
i = 3, j = 3, k = 1
i = 3, j = 3, k = 2
i = 3, j = 3, k = 3
counts = 64
```
<br />
---
###### tags:`VPython`
Sign in with Wallet
Connect another wallet