scipy計算三重積分
阿新 • • 發佈:2019-02-19
圓柱體直徑4cm,長20cm,兩端各有一個球,每個球直徑20cm,密度為5600kg/m3計算rotational inertia matrix .座標原點為圓柱的中心,沿長度的方向為z軸。
先在草稿紙上寫出直角座標的積分公式,將其寫入程式即可。計算三重積分用到的函式是tplquad,此外還有二重積分,定積分,n重積分。f為被積函式,cyl為求圓柱的積分,sphereup求上球的積分,spheredown求下球的積分,返回一個元組,第一個元素是積分值,第二個元素是精度。
from scipy import integrate f = lambda x, y, z: x**2+y**2 cyl=integrate.tplquad(f, -0.1, 0.1, lambda y: -0.02, lambda y: 0.02, lambda x, y: -np.sqrt(0.0004-y**2), lambda x, y: np.sqrt(0.0004-y**2)) sphereup=integrate.tplquad(f, 0.1, 0.3, lambda z: -np.sqrt(0.40*z - z**2 - 0.03), lambda z: np.sqrt(0.40*z - z**2 - 0.03), lambda z,y: -np.sqrt(0.40*z - z**2 - 0.03 - y**2), lambda z,y: np.sqrt(0.40*z - z**2 - 0.03 - y**2)) spheredown=integrate.tplquad(f, -0.3, -0.1, lambda z: -np.sqrt(-0.40*z - z**2 - 0.03), lambda z: np.sqrt(-0.40*z - z**2 - 0.03), lambda z,y: -np.sqrt(-0.40*z - z**2 -0.03 - y**2), lambda z,y: np.sqrt(-0.40*z - z**2 - 0.03 - y**2)) print((cyl[0]+sphereup[0]+spheredown[0])*5600)
當把被積函式改為1時,即為求圓柱體和球體的體積:
f=lambda x,y,z:1
此時計算的結果為
print(cyl[0],sphereup[0],spheredown[0])
0.0002513274461536556
0.004188790206670827
0.004188790206670827
用球體、圓柱體體積公式計算結果為:
In[5]:4/3*np.pi*0.001
Out[5]: 0.0041887902047863905
In[6]: np.pi*0.02**2*0.2
Out[6]: 0.0002513274122871835