Question

求助matlab程序题，球体容积量尺设计

一个直径为3米的球体。一根4米长的钢尺作为量尺进行测量（如下图）。

被油浸湿的标尺刻度就是储存罐中油的高度，一旦知道罐中油的高度“h”，则油罐中剩余油的体积“V”直接由如下公式算出：

其中 r是油罐的半径。
但这样仍然不直接，需设计一根量尺，使得油罐里油的体积可直接从量尺上读出来，我该怎样设计量尺上的刻度呢？
有完整的建模过程（如是程序设计，算法思想要写清楚），求助大神！





Answer 1

如果不使用积分, 分h=4米到许多小的间隔h(i )(i=1,n,n很大),r=1.5米,求多项式v(i),i=1,n ,
v(i)=3.14159*h(i)^2*(3*1.5-h(i))/3.

求尺上的刻度U(m)=∑(i=2,m) v(i)-v(i-1);m≤n

如果使用积分,
v=3.14159*h^2*(3*1.5-h)/3.
先找到函数v的导数
v'=2*3.14159*h*(3*1.5-h)/3.-[3.14159*h^2]/3.
%曲线长度=计算尺上的刻度U值=∫(起点,终点 )[(1+(v')^2)]^(1/2)dx

---------------------MATLAB:
myfun = @(h) 2*3.14159*h*(3*1.5-h)/3.-[3.14159*h^2]/3.
from = 0;%起点
to = 4;%终点
U = quad(myfun,from,to)
%使用循环为不同终点值to,计算尺上的刻度U值

追问

要求直接输出设计好的量尺图，可以吗？还有完整matlab程序。。。

追答

举例来说:
How about something like this--it basically zeros out the Y tick marks and tick lables and creates a thin plot window:
screenSize = get( 0, 'ScreenSize' );
screen_height = screenSize(4);
screen_width = screenSize(3);
width = screen_width * 0.5;
left = floor( (screen_width - width) / 2 );
height = 170;
bottom = floor( (screen_height - height) / 2 );
figure( 'Position', [ left bottom width height] );

x = cumsum( abs( randn( 10, 1 ) ) );
z = zeros( 10, 1 );

plot( x, z, 'r+', 'MarkerSize', 10, 'LineWidth', 1 );
set( gca, 'YTick', -1:1:1 );
set( gca, 'YTickLabel', {'','',''} );

hold on;
xlimits = xlim;
plot( xlimits, [ 0 0 ], 'Color', 'black' );
grid on;

追问

so sorry，我弄了半天，实在不理想，真的能否直接给出可行程序，以及算法思路，如果可以我再追加100，如果不好弄，可以选择把做好的m文件上传。

追答

%再追加100,谢谢
h=ones(1,10);
u=ones(1,10);

for i=1:10
h(i)=i*4./10.;
myfun = @(x) 2*3.14159*x*(3*1.5-x)/3.-[3.14159*x^2]/3.;
from = 0;%起点;
to = h(i);%终点 ;

%曲线长度=计算尺上的刻度u值=∫(起点,终点 )[(1+(v')^2)]^(1/2)dx;
%v=3.14159*h^2*(3*1.5-h)/3.
% v'=v的导数
u(i) = quad(myfun,from,to);
end
screenSize = get( 0, 'ScreenSize' );
screen_height = screenSize(4);
screen_width = screenSize(3);
width = screen_width * 0.5;
left = floor( (screen_width - width) / 2 );
height = 100;
bottom = floor( (screen_height - height) / 2 );
figure( 'Position', [ left bottom width height] );
z = zeros( 10, 1 );
plot( u, z, 'r+', 'MarkerSize', 10, 'LineWidth', 1 );
set( gca, 'YTick', -1:1:1 );
set( gca, 'YTickLabel', {'','',''} );
hold on;
xlimits = 4.0;
plot( xlimits, [ 0 0 ], 'Color', 'black' );
grid on;

追问

Error using  * 

Inner matrix dimensions must agree.

Error in @(x)2*3.14159*x*(3*1.5-x)/3.-[3.14159*x^2]/3.

Error in quad (line 72)

...

我这里好多错误提示，你把可运行的m文件直接上传行不，谢了！

追答

不允许[ ],更改为()
myfun = @(x) 2*3.14159*x*(4.5-x)/3.0-(3.14159*x^2)/3.0;
或者你也可以使用
myfun = inline('2*3.14159*x*(4.5-x)/3.0-(3.14159*x^2)/3.0')

追问

还是不行。。。哭了，最好直接上传可运行m文件，我用的是2012b

貌似不是括号问题，矩阵计算维数不相容。

追答

%sorry. 要有耐心

h=zeros(10,1);
u=zeros(10,1);

for i=1:10,
h(i)=i*4./10.;
myfun = @(x) (2.0*3.14159*x*(4.5-x)/3.0-(3.14159*x*x)/3.0);
from = 0;
to = h(i);
u(i) = quad(myfun,from,to);
%或者你也可以使用
%u(i)=quad(@(x) (2.0*3.14159*x*(4.5-x)/3.0-(3.14159*x*x)/3.0),0,h(i));
end
%或者你也可以使用
%f1=@(x) (2.0*3.14159*x*(4.5-x)/3.0);
%f2=@(x) ((3.14159*x*x)/3.0);
%myfun=f1+f2;

追问

。。。你电脑里没有matlab吗，没有验证吗。。。还是错的啊

追答

在家电脑里没有matlab,

h=zeros(10,1);
u=zeros(10,1);
for i=1:10,
h(i)=i*4./10.;
x=0:h(i)/100:h(i);
myfun = @(x) (2.0*3.14159*x*(4.5-x)/3.0-(3.14159*x*x)/3.0);
u(i) = trapz(x,myfun);
end

追问

Undefined function 'max' for input arguments of type 'function_handle'.

Error in trapz (line 43)
perm = [dim:max(ndims(y),dim) 1:dim-1];
还是不对，再追问。。。一条扣30财富了。。。要不你啥时能有matlab弄下...

追答

%the current folder would need to be added to the working directory,
%首先尝试简单的程序,看函数是不是太复杂
x=linespace(0,4./100.,4);
myfun(x) = @(x) (((2.0*3.14159*x*(4.5-x)/3.0-(3.14159*x*x)/3.0)^2+1.)^(1/2));
y=myfun(x);
u = trapz(x,y);

%如果上述程序成功，那么请尝试以下
u=zeros(10);
for i=1:10,
h=i*4./10.;
x=linespace(0,h/100,h);
myfun(x) = @(x) (((2.0*3.14159*x*(4.5-x)/3.0-(3.14159*x*x)/3.0)^2+1.)^(1/2));
y=myfun(x);
u = trapz(x,y);
end

追问

不行

追答

x=linspace(0,4,100) ,not linespace, sorry
x=linspace(0,h,100)