第三单元用python学习微积分（二十）壳层法、圆盘法求体积 (下）

本文内容来自于学习麻省理工学院公开课：单变量微积分-壳层法、圆盘法求体积-网易公开课

开发环境准备：CSDN

考虑 y = x^2 的曲线在y轴右侧的部分，如果绕y轴旋转360度，会形成一个像水缸的体，而如上图中红色矩形部分，会形成一个空心的圆柱，一个圆柱形的壳。这个壳的厚度为dx(趋近于无穷小)，高度为a-y(注意最左侧图上y=a直线), 半径为x。

这个壳展开为一个厚片（最右侧图），这个厚片的厚度为dx，高度为a-y, 底边长为中间图中圆柱的底面圆形的周长 $2\pi x$ 。

$dV = (2 \pi x)(a-y)dx = 2\pi(ax - x^3)dx ( 0<x<\sqrt[]{a} , y\leq a , y = x^2 )$

$V = \int_{0}^{\sqrt a} 2\pi(ax - x^3)dx = \pi a x^2 - \frac{\pi x^4}{2}|_{0}^{\sqrt a} = \pi a^2 - \frac{\pi a^2}{2} = \frac{\pi}{2} a^2$

代码主体部分请参考Bullseye：第三单元用python学习微积分（二十）壳层法、圆盘法求体积 (上）

fig = plt.figure(figsize=(12, 12),
                 facecolor='lightyellow'
                )
# draw sphere
ax = fig.gca(fc='whitesmoke',
               projection='3d' 
              )
x = symbols('x')
a = 4
expr = x**2

DrawXY(0,2,100,expr,'b','',plt)

RotateByAxis(0,2,100,expr,-0.5*np.pi,'r',ax, EnumAxis.YAxis)
RotateByAxis(0,2,100,expr,0.5*np.pi,'r',ax, EnumAxis.YAxis)
RotateByAxis(0,2,100,expr,-np.pi,'b',ax, EnumAxis.YAxis)

DrawAxis(ax, 5,5,5)

   
plot_surface(-2, 4, -2, 4, 0, 4, ax)
ax.view_init(elev=30,    # 仰角
             azim=-20    # 方位角
            )
plt.show()

这里老师提示要注意单位！

设 a=100 cm

则体积为： $V= \frac{\pi}{2}(100)^2 (cm^3) = \frac{10\pi}{2}(1000) (cm^3) = \frac{10\pi}{2} (L)$

设 a=1 m

则体积为： $V= \frac{\pi}{2} (m^3) = \frac{\pi}{2}(100)^3 (cm^3) = \frac{1000\pi}{2} (L)$

所以都是1m的长度，体积的结果确实不同，为啥呢？老师的原话是公式 $\frac{\pi}{2} a^2$ 违反了比例原则。

其实是都是合理的，造成不同的原因是曲线的形状不同！如下：

import matplotlib.pyplot as plt
fig=plt.figure()

axes=[]
#define matplotlib figure and axis

def DrawXY(xFrom,xTo,steps,expr,color,label,ax, aspect):
    yarr = []
    xarr = np.linspace(xFrom ,xTo, steps) 
    for xval in xarr:
        yval = expr.subs(x,xval)
        yarr.append(yval)
    y_nparr = np.array(yarr) 
    ax.set_aspect(aspect)
    ax.plot(xarr, y_nparr, c=color, label=label)   

x = symbols('x')
expr = x**2

# add_subplot(row, col, row*col) 添加图在row行，col列 
axes.append( fig.add_subplot(1, 1,1) )
DrawXY(-10,10,100,expr,'green','y = x**2 a=100 cm ',axes[len(axes)-1], 0.8)
plt.legend(loc='upper right')
axes.append( fig.add_subplot(1, 3,3) )
DrawXY(-1,1,100,expr,'green','y = x**2 a=1 m ',axes[len(axes)-1], 1)
plt.legend(loc='upper right')

#set aspect ratio to 1
fig.tight_layout()
#display plot
plt.show()

由积分的性质，我们知道求取的是旋转曲线所得到的空间形状的体积，而由于两图形的x、y比例不同，所得到的体积也就不同。上图中左边的图形当x取10cm时，图形的高是100cm，而右边的图形，x取100cm时，图形的高是100cm。

所以当处理实际问题是，一定要考虑好比例关系。

下节课将要讲在像这样的坩埚中加热其中全部的水，需要多少热量。

2、壳层法（method of shells）求体积