Vpython作業8

from vpython import *

"""
 1. 參數設定, 設定變數及初始值
"""
m = 1               # 小球質量
size = 0.2          # 小球半徑
L = 5               # 繩長
theta = radians(30) # 繩子與鉛垂線夾角
R = L*sin(theta)    # 錐動擺半徑
H = L*cos(theta)    # 錐動擺平面與懸掛點高度差
g = 9.8             # 重力加速度
v0 = 1*sqrt(g*R)               # 切線速度量值
T = 2.88415              # 週期理論值
Fc = m*g/cos(theta)               # 向心力量值
i = 0               # 小球經過週期次數
N = 5               # 小球經過週期次數上限
t = 0               # 時間
dt = 0.001          # 時間間隔

"""
 2. 畫面設定
"""
# 產生動畫視窗、天花板、小球、繩子
scene = canvas(title="Cone Pendulum", width=600, height=600, x=0, y=0, background=color.black, range=L)
roof = box(pos=vec(0, 0.5*(L+size), 0), size=vec(L, size, 0.5*L), color=color.blue)
ball = sphere(pos=roof.pos-vec(R, 0, 0), radius=size, color=color.red,
              make_trail=True, retain=200, v=vec(0, 0, -v0))
rope = cylinder(pos=roof.pos, axis=ball.pos, radius=0.1*size, color=color.yellow)
ver = curve(pos=[roof.pos, vec(0, roof.pos.y-H, 0)], rasius=0.1*size, color=color.magenta)

# 產生表示速度的箭頭
arrow_v = arrow(pos=ball.pos, axis=vec(0, 0, 0), shaftwidth=0.3*size, color=color.green)
arrow_a = arrow(pos=ball.pos, axis=vec(0, 0, 0), shaftwidth=0.3*size, color=color.cyan)

"""
 3. 物體運動部分, 重覆5個週期
"""
pos1 = ball.pos.z
print("週期理論值 T = %f s" % T)

while i < N:
    rate(1000)
# 計算小球加速度, 更新速度, 位置, 繩子
    r = rope.axis               # 圓心到小球位置的向量
    ball.a = -Fc/m*r.norm()               # 圓心到小球位置的單位向量
    ball.v += ball.a*dt
    ball.pos += ball.v*dt
    rope.axis = ball.pos-roof.pos
# 更新表示速度, 加速度的箭頭
    arrow_v.pos = ball.pos
    arrow_v.axis = 0.5*ball.v
    arrow_a.pos = ball.pos
    arrow_a.axis = 0.5*ball.a
# 檢驗小球是否經過一個週期
    pos2 = ball.pos.z
    if pos1>0 and pos2<0: # 由小球位置判斷是否通過出發點
        print("第 %d 次通過原點, t = %f s" %(i, t))
        i += 1
    pos1 = pos2
# 更新時間
    t += dt

第27行的向量不會算，初速度不知道這樣算對不對，45行求向量(應該是這樣吧)，46要求單位向量，所以要取norm()