牟合方盖的绘图与体积计算问题
沈其松 学号：[1**********]8
一：问题叙述：
魏晋时数学家刘徽在他的著作《九章算术注》中指出我国古代数学名著《九
9
章算术》中的球体积公式V =d 3（d 为球的直径）是错误的，错误的原因在
16
于误以为球和它的外切圆柱的体积的比是π∶4。为了纠正这一错误，刘徽在他的《九章算术注》中，提出一个独特的方法来计算球体的体积：他不直接求球体的体积，而是先计算另一个叫“牟合方盖”的立体的体积。
正方体内两轴互相垂直的内切圆柱面相交所围的空间立体。由于这个立体的外形如同两把上下对称的正方形雨伞，所以称它为牟合方盖。
刘徽通过计算，球体体积与“牟合方盖”的体积之比应为 π : 4；显然, 只要求出牟合方盖的体积，那么球体积便迎刃而解。可惜的是，刘徽功亏一篑，未能求出牟合方盖的体积。
所以本试验用MATLAB 画出牟合方盖，并用“祖暅方法”，“微积分方法”，“蒙特卡罗方法”，分别计算牟合方盖的体积，来实现刘徽的愿望。
二：问题分析：
1. 绘制牟合方盖
绘制柱面x2 + y2 = R2与柱面x2 + z2 =R2所围立体在x-y 平面上半部分曲面。由第二个方程解出z ，得
z =
， x =r cos t (0≤r ≤R ); y =r sin t (0≤t ≤2π);
则可以画出对应的曲面，当画四分之一，八分之一曲面时，只需设置r 与t 的范围就可以了。 2. 计算牟合方盖体积
2.1 祖暅方法：祖暅沿用了刘徽的思想，利用刘徽“牟合方盖”的理论去进行体积计算。由于没有微积分，祖暅用一种等效的方法来计算。
图 一
图 二
图 三
他的方法是将原来的“牟合方盖”平均分为八份，取它的八分之一（如图一），设 OP = h ，过 P 点作平面 PQRS 平行于 OABC 。又设内切球体的半径为 r ，则 OS = OQ = r ，由勾股定理有PS = PQ =r 2-h 2，故此正方形 PQRS 面积是 r 2 - h 2。如果将图一的立体放在一个边长为 r 的正立方体之内（如图二），不难证明图二中与图一等高处阴影部分的面积等于h 2。（如图三）设由方锥顶点至方锥截面的高度为h ，不难发现对于任何的h ，方锥截面面积也必为h 2。由此可知，在等高处，图二中阴影部分的面积与图三中倒立的正立方锥体的横切面的面积总相等。所以，有理由相信，虽然方锥跟小正立方体去掉小“牟合方盖”后的形状不同，但因它们的体积都可以用截面面积和高度来计算，而在等高处的截面面积总是相等的，所以它们的体积相等。所以V 牟=V正-V 锥。
2.2微积分方法：与祖暅的分析类似，牟合方盖的截面积为r 2 -h 2，高为h ，
r
体积为：⎰(r 2-h 2) δh ；由于MTALAB 内集成了以Maple 的内核开发了Matlab 的
符号计算工具箱。可以进行积分公式的符号计算，有了积分工具，可以用matlab 直接计算出牟合方盖的体积公式。
2.3蒙特卡罗方法：通过随机变量的统计试验求近似解，对被积函数变量区间进行随机均匀抽样，然后对被抽样点的函数值求平均，从而可以得到函数积分的近似值。此种方法的正确性是基于概率论的中心极限定理。当抽样点数为m 时，使用此种方法所得近似解的统计误差恒为1除于根号M, 不随积分维数的改变而改变。
实验中用rand （n,3) 产生n*3个0到 1之间均匀随机数（n 取值较大），随机数较均匀地分布在正方体内，随机变量X 落入某个小空间内的概率仅与小空间的体积有关, 而与小空间间位置无关，通过find （）函数，来统计落入牟合方盖点的个数为m 个，最终体积为V=8*m/n。
三：实验程序及注释
%%%%%%%%%%画牟合方盖的图形%%%%%%%%%%%%%%%%%%%%%%%%
%%牟合方盖的全图程序%%
t=(0:40)/40*pi;
x=cos(t); %%转换为极坐标 y=sin(t); z=y;
X=[1;1;1;1;1]*x; Y=[1;-1;-1;1;1]*y; Z=[1;1;-1;-1;1]*z;
figure(1),subplot(2,2,1),mesh(X,Y,Z);title('牟合方盖全图') %%牟合方盖的1/2图程序%%
h=2*pi/100; t=0:h:2*pi;r=0:0.05:1;
x=r'*cos(t); %%转换为极坐标 y=r'*sin(t); z=sqrt(1-x.^2);
figure(1),subplot(2,2,2),meshz(x,y,z); title('牟合方盖1/2图') colormap([0 0 1]) %%牟合方盖的1/4图程序%%
h=2*pi/100;t=0:h:pi;r=0:0.05:1;
x=r'*cos(t);y=r'*sin(t); %%转换为极坐标 zz=sqrt(1-x.^2);
figure(1),subplot(2,2,3),meshz(x,y,zz);title('牟合方盖的1/4图'); colormap([0 0 1]) %axis off view(120,34)
%%牟合方盖的1/8图程序%%
h=2*pi/100;t=0:h:pi/2;r=0:0.05:1;
x=r'*cos(t);y=r'*sin(t); %%转换为极坐标 zz=sqrt(1-x.^2);
figure(1),subplot(2,2,4),meshz(x,y,zz);title('牟合方盖的1/8图') colormap([0 0 1]) %axis off
view(120,34)
%%%%%%%%%%计算牟合方盖的体积%%%%%%%%%%%%%%%%%%%%%%%%
%%祖暅方法%%
format long
r=1; %%设定相交圆柱的半径为1 Vs=r^3; %%立方体的体积公式 Vl=(r^2)*r*1/3; %%棱锥的体积公式
V1=(Vs-Vl)*8; %%牟和方盖的面积等于正方体面积减棱锥面积 %%微积分方法%%
syms h ; %%设定积分变量 r=1;
f=[r^2 - h^2]; %%设定积分函数 Vn=int(f,h,0,r); %%求定积分
V2=8*Vn; %%求出牟和方盖的体积 %%蒙特卡罗方法%%
format long
n=1000000; %%随机分布点数 m=rand(n,3);
x=m(:,1); %%分布到正方体中 y=m(:,2); z=m(:,3);
II=find(x.^2+y.^2
V3=8*k/n； %%计算体积 Yn=(V3-16/3)/(16/3); %%计算体积误差
四：实验数据结果及分析
1．画出牟和方盖的图形
5．实验结论
本试验通过将直角坐标转换为极坐标，画出了牟和方盖的全图、1/2、1/4、1/8的三维网面图，并在1/8的牟和方盖上来计算体积，为了计算牟和方盖的体积，使用了祖暅方法、微积分方法、蒙特卡罗方法。
祖暅方法是通过一个计算一个等体积的正棱锥来计算牟和方盖的体积，是一种间接的巧妙方法。而且计算的结果是无误差的标准值。
微积分方法是通过MATLAB 内的积分函数，直接进行积分运算，得出体积公式，是一种直接的先进的方法。而且计算的结果是无误差的标准值。
蒙特卡罗方法是通过统计随机分布在牟和方盖内的点来计算方盖体积，是一种间接的，适应很广的方法，但是计算结果是有误差的，误差主要受随机点的个数的影响，随机点数越多，结果越精确。
对于半径为1的两个正交圆柱所截的牟和方盖的体积为16/3。所以实现了刘徽的愿望。
牟合方盖的绘图与体积计算问题
沈其松 学号：[1**********]8
一：问题叙述：
魏晋时数学家刘徽在他的著作《九章算术注》中指出我国古代数学名著《九
9
章算术》中的球体积公式V =d 3（d 为球的直径）是错误的，错误的原因在
16
于误以为球和它的外切圆柱的体积的比是π∶4。为了纠正这一错误，刘徽在他的《九章算术注》中，提出一个独特的方法来计算球体的体积：他不直接求球体的体积，而是先计算另一个叫“牟合方盖”的立体的体积。
正方体内两轴互相垂直的内切圆柱面相交所围的空间立体。由于这个立体的外形如同两把上下对称的正方形雨伞，所以称它为牟合方盖。
刘徽通过计算，球体体积与“牟合方盖”的体积之比应为 π : 4；显然, 只要求出牟合方盖的体积，那么球体积便迎刃而解。可惜的是，刘徽功亏一篑，未能求出牟合方盖的体积。
所以本试验用MATLAB 画出牟合方盖，并用“祖暅方法”，“微积分方法”，“蒙特卡罗方法”，分别计算牟合方盖的体积，来实现刘徽的愿望。
二：问题分析：
1. 绘制牟合方盖
绘制柱面x2 + y2 = R2与柱面x2 + z2 =R2所围立体在x-y 平面上半部分曲面。由第二个方程解出z ，得
z =
， x =r cos t (0≤r ≤R ); y =r sin t (0≤t ≤2π);
则可以画出对应的曲面，当画四分之一，八分之一曲面时，只需设置r 与t 的范围就可以了。 2. 计算牟合方盖体积
2.1 祖暅方法：祖暅沿用了刘徽的思想，利用刘徽“牟合方盖”的理论去进行体积计算。由于没有微积分，祖暅用一种等效的方法来计算。
图 一
图 二
图 三
他的方法是将原来的“牟合方盖”平均分为八份，取它的八分之一（如图一），设 OP = h ，过 P 点作平面 PQRS 平行于 OABC 。又设内切球体的半径为 r ，则 OS = OQ = r ，由勾股定理有PS = PQ =r 2-h 2，故此正方形 PQRS 面积是 r 2 - h 2。如果将图一的立体放在一个边长为 r 的正立方体之内（如图二），不难证明图二中与图一等高处阴影部分的面积等于h 2。（如图三）设由方锥顶点至方锥截面的高度为h ，不难发现对于任何的h ，方锥截面面积也必为h 2。由此可知，在等高处，图二中阴影部分的面积与图三中倒立的正立方锥体的横切面的面积总相等。所以，有理由相信，虽然方锥跟小正立方体去掉小“牟合方盖”后的形状不同，但因它们的体积都可以用截面面积和高度来计算，而在等高处的截面面积总是相等的，所以它们的体积相等。所以V 牟=V正-V 锥。
2.2微积分方法：与祖暅的分析类似，牟合方盖的截面积为r 2 -h 2，高为h ，
r
体积为：⎰(r 2-h 2) δh ；由于MTALAB 内集成了以Maple 的内核开发了Matlab 的
符号计算工具箱。可以进行积分公式的符号计算，有了积分工具，可以用matlab 直接计算出牟合方盖的体积公式。
2.3蒙特卡罗方法：通过随机变量的统计试验求近似解，对被积函数变量区间进行随机均匀抽样，然后对被抽样点的函数值求平均，从而可以得到函数积分的近似值。此种方法的正确性是基于概率论的中心极限定理。当抽样点数为m 时，使用此种方法所得近似解的统计误差恒为1除于根号M, 不随积分维数的改变而改变。
实验中用rand （n,3) 产生n*3个0到 1之间均匀随机数（n 取值较大），随机数较均匀地分布在正方体内，随机变量X 落入某个小空间内的概率仅与小空间的体积有关, 而与小空间间位置无关，通过find （）函数，来统计落入牟合方盖点的个数为m 个，最终体积为V=8*m/n。
三：实验程序及注释
%%%%%%%%%%画牟合方盖的图形%%%%%%%%%%%%%%%%%%%%%%%%
%%牟合方盖的全图程序%%
t=(0:40)/40*pi;
x=cos(t); %%转换为极坐标 y=sin(t); z=y;
X=[1;1;1;1;1]*x; Y=[1;-1;-1;1;1]*y; Z=[1;1;-1;-1;1]*z;
figure(1),subplot(2,2,1),mesh(X,Y,Z);title('牟合方盖全图') %%牟合方盖的1/2图程序%%
h=2*pi/100; t=0:h:2*pi;r=0:0.05:1;
x=r'*cos(t); %%转换为极坐标 y=r'*sin(t); z=sqrt(1-x.^2);
figure(1),subplot(2,2,2),meshz(x,y,z); title('牟合方盖1/2图') colormap([0 0 1]) %%牟合方盖的1/4图程序%%
h=2*pi/100;t=0:h:pi;r=0:0.05:1;
x=r'*cos(t);y=r'*sin(t); %%转换为极坐标 zz=sqrt(1-x.^2);
figure(1),subplot(2,2,3),meshz(x,y,zz);title('牟合方盖的1/4图'); colormap([0 0 1]) %axis off view(120,34)
%%牟合方盖的1/8图程序%%
h=2*pi/100;t=0:h:pi/2;r=0:0.05:1;
x=r'*cos(t);y=r'*sin(t); %%转换为极坐标 zz=sqrt(1-x.^2);
figure(1),subplot(2,2,4),meshz(x,y,zz);title('牟合方盖的1/8图') colormap([0 0 1]) %axis off
view(120,34)
%%%%%%%%%%计算牟合方盖的体积%%%%%%%%%%%%%%%%%%%%%%%%
%%祖暅方法%%
format long
r=1; %%设定相交圆柱的半径为1 Vs=r^3; %%立方体的体积公式 Vl=(r^2)*r*1/3; %%棱锥的体积公式
V1=(Vs-Vl)*8; %%牟和方盖的面积等于正方体面积减棱锥面积 %%微积分方法%%
syms h ; %%设定积分变量 r=1;
f=[r^2 - h^2]; %%设定积分函数 Vn=int(f,h,0,r); %%求定积分
V2=8*Vn; %%求出牟和方盖的体积 %%蒙特卡罗方法%%
format long
n=1000000; %%随机分布点数 m=rand(n,3);
x=m(:,1); %%分布到正方体中 y=m(:,2); z=m(:,3);
II=find(x.^2+y.^2
V3=8*k/n； %%计算体积 Yn=(V3-16/3)/(16/3); %%计算体积误差
四：实验数据结果及分析
1．画出牟和方盖的图形
5．实验结论
本试验通过将直角坐标转换为极坐标，画出了牟和方盖的全图、1/2、1/4、1/8的三维网面图，并在1/8的牟和方盖上来计算体积，为了计算牟和方盖的体积，使用了祖暅方法、微积分方法、蒙特卡罗方法。
祖暅方法是通过一个计算一个等体积的正棱锥来计算牟和方盖的体积，是一种间接的巧妙方法。而且计算的结果是无误差的标准值。
微积分方法是通过MATLAB 内的积分函数，直接进行积分运算，得出体积公式，是一种直接的先进的方法。而且计算的结果是无误差的标准值。
蒙特卡罗方法是通过统计随机分布在牟和方盖内的点来计算方盖体积，是一种间接的，适应很广的方法，但是计算结果是有误差的，误差主要受随机点的个数的影响，随机点数越多，结果越精确。
对于半径为1的两个正交圆柱所截的牟和方盖的体积为16/3。所以实现了刘徽的愿望。