24小时热门版块排行榜    

查看: 890  |  回复: 5

815292578

木虫 (著名写手)

[求助] MATLAB求解方程后没有报错,但无法绘制出曲线,求高手帮忙解决... 已有1人参与

最近利用MATLAB求解方程后没有报错,但无法绘制出曲线,求高手帮忙解决...
方程不是很复杂,具体方程如下:
clear;clc;
P=29.6; D=8245; R=0.021;  LL=0.021; x0=0.010;  pp=8960;  a=0.002;  b=0.002;  erfa=30*pi/180;  angle=atan(LL/R);  
n=10;  
h=a/n;  
N=x0/h;  %步长总数
k=1.2;  
t=1/(2*k);
Y=zeros(N,1);
L=zeros(N,1);
H=zeros(N,1);
L(1)=R*tan(angle);
H(1)=0;
x(1)=h;
y(1)=h*tan(erfa);
for i=1:1:N-1;
    x(i+1)=(i+1)*h;  %离散后xi横坐标
    y(i+1)=x(i+1)*tan(erfa);  %离散后yi纵坐标
    psaiI=atan((LL+(i+1)*h)/R);
    z=fsolve(@(z)tan(psaiI)*sqrt((1-(2*z-1)/z^2)*(1-1/((2*z-1)^t)^2))-1/((2*z-1)^t)-sqrt((2*z-1)/z^2),1);
    Dm=D*sqrt(z^2/(2*z-1));
    xx=fsolve(@(x)cos(x)/sin(psaiI-x)-sqrt(z^2/(2*z-1)),0);   
    X=[x,xx];
    L(i+1)=(L(i)*tan(X(i))*tan(psaiI)+tan(psaiI)*(R-H(i)))/(1+tan(X(i))*tan(psaiI));
    H(i+1)=(R*tan(X(i))*tan(psaiI)-L(i)*tan(X(i))+H(i))/(1+tan(X(i))*tan(psaiI));     
if H(i)>Y(i)
    Pm=z*P;
    else
    Pm=P;
    end
    if i==1
        vv(1)=(1/3)*pi*(h*tan(erfa))^2*h;  
        m(1)=pp*vv(1);
        v=Pm/m(1);
    elseif i>n
            vv(i)=(1/3)*pi*(i*h*tan(erfa))^2*i*h-vv(n)-(1/3)*pi*((i-n)*h*tan(erfa))^2*(i-n)*h;
            m(i)=pp*vv(i);
            v(i)=Pm/m(i);
    else
        vv(i)=(1/3)*pi*(i*h*tan(erfa))^2*i*h-vv(i-1);
        m(i)=pp*vv(i);
        v(i)=Pm/m(i);
    end
end
plot(x(i),v(i))
想绘制横坐标为x, 纵坐标为v的曲线。请熟悉MATLAB的高手给点建议,谢谢...
回复此楼

» 猜你喜欢

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

FMStation

至尊木虫 (知名作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
815292578: 金币+10, 有帮助 2016-08-15 06:45:41
plot(x(i),v(i))  => plot(x',v')

whos x v
  Name      Size            Bytes  Class     Attributes

  v         1x49              392  double              
  x         1x50              400  double              

維度不同
2楼2016-08-13 20:47:14
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

815292578

木虫 (著名写手)

引用回帖:
2楼: Originally posted by FMStation at 2016-08-13 20:47:14
plot(x(i),v(i))  => plot(x',v')

whos x v
  Name      Size            Bytes  Class     Attributes

  v         1x49              392  double              
  x         1x50              400   ...

谢谢您的回复。现在我将程序拆分出来分别计算求得x,z值; H值;Pm值;和微元质量m。最后绘制plot(x,v)。发现还是不行!
发现:这样可以得到离散的x(i)值,y(i)值,X(i)值。但是离散的z(i)值得不到!所以最后无法得到Pm=z*p值。请问给出意见如何修改?谢谢...

clear;clc;
P=29.6; D=8245; R=0.021;  LL=0.021; x0=0.010;  pp=8960;  a=0.002;  b=0.002;  erfa=30*pi/180;  angle=atan(LL/R);  
n=10;  
h=a/n;  
N=x0/h;  %步长总数
k=1.2;  
t=1/(2*k);
%求X值和z值
for i=1:N
    x(i)=i*h;
    y(i)=x(i)*tan(erfa);
    psaiI=atan((LL+i*h)/R);
    z=fsolve(@(z)tan(psaiI)*sqrt((1-(2*z-1)/z^2)*(1-1/((2*z-1)^t)^2))-1/((2*z-1)^t)-sqrt((2*z-1)/z^2),1);  
    xx=fsolve(@(x)cos(x)/sin(psaiI-x)-sqrt(z^2/(2*z-1)),0);   
    X(i)=xx;   
end
%求H值
L=zeros(N,1);
H=zeros(N,1);
L(1)=R*tan(angle);
H(1)=0;
for i=1:N-1
    L(i+1)=(L(i)*tan(X(i))*((LL+i*h)/R)+((LL+i*h)/R)*(R-H(i)))/(1+tan(X(i))*((LL+i*h)/R));
    H(i+1)=(R*tan(X(i))*((LL+i*h)/R)-L(i)*tan(X(i))+H(i))/(1+tan(X(i))*((LL+i*h)/R));   
end
%比较上面计算出的H和y值的大小后,求Pm值
if H(i)>y(i)
    Pm=z*P;
else Pm=P;
end
%求各个微元的质量
for i=1:N
    if i==1
        vv(1)=(1/3)*pi*(h*tan(erfa))^2*h;  %第一个微元体积
        m(1)=pp*vv(1);
        v(1)=Pm/m(1);
    elseif i>n
        vv(i)=(1/3)*pi*((i)*h*tan(erfa))^2*(i)*h-vv(n)-(1/3)*pi*((i-n)*h*tan(erfa))^2*(i-n)*h;
        m(i)=pp*vv(i);
        v(i)=Pm/m(i);
    else
        vv(i)=(1/3)*pi*((i)*h*tan(erfa))^2*(i)*h-vv(i);
        m(i)=pp*vv(i);
        v(i)=Pm/m(i);
    end
end
plot(x,v)
3楼2016-08-15 06:44:45
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

FMStation

至尊木虫 (知名作家)

【答案】应助回帖

>> vv(i)=(1/3)*pi*((i)*h*tan(erfa))^2*(i)*h-vv(i);
>> whos vv
  Name      Size            Bytes  Class     Attributes
  vv        1x1                 8  double              

??? Index exceeds matrix dimensions.
i =  2 => No vv(2)
4楼2016-08-15 09:56:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

815292578

木虫 (著名写手)

引用回帖:
4楼: Originally posted by FMStation at 2016-08-15 09:56:39
>> vv(i)=(1/3)*pi*((i)*h*tan(erfa))^2*(i)*h-vv(i);
>> whos vv
  Name      Size            Bytes  Class     Attributes
  vv        1x1                 8  double              

??? Ind ...

谢谢您的回复。vv(i)是指离散后各个微元的体积。由于求解公式不同,当i=1;1<i<=n;n<i<=N.三个不同公式求解的。
接下来该如何求解?能给点建议吗...谢谢
5楼2016-08-16 03:36:42
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

FMStation

至尊木虫 (知名作家)

【答案】应助回帖

當 i =  2 時
欲執行
vv(i)=(1/3)*pi*((i)*h*tan(erfa))^2*(i)*h-vv(i);
可是沒有已產生的 vv(2) ?

公式錯 or  vv(2)須先產生?
6楼2016-08-16 08:22:46
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 815292578 的主题更新
信息提示
请填处理意见