[已解决] 求积分上限那出错,该怎么改呀

[复制链接]
一只小萌新^-^ 发表于 2021-5-18 21:31:45
10 财富积分
  1. R=10;
  2. r=8;
  3. rc=1;
  4. H=2;
  5. o=0;
  6. theata=[];
  7. T0=[];
  8. X=[];
  9. Y=[];
  10. while o~=91
  11.     syms t0;
  12.     eq=((R*cos(t0)-(r+rc)*sind(o))^2)*(1+(tan(t0)/sind(o))^2)-(rc*sind(o))^2;
  13.     t0=vpasolve(eq,t0,[0,pi/2]);%求出t0的值
  14.     x0=R*cos(t0)/sind(o);%利用t0求出x的值
  15.     y0=R*sin(t0);%利用t0求出y的值
  16.     h=rc.^2-(x0-r-rc).^2;
  17.     L=R+H-y0-h.^0.5+rc*asin((h.^0.5)/rc);%变形长度L
  18.     F=@(t1)integral(@(t)((cos(t)/sind(o)).^2+(sin(t)).^2).^0.5,t0,t1)-L;
  19.     t1=fsolve(F,0.5);
  20.     x=R*cos(t1)/sind(o);
  21.     theata=[theata,o];%theta角的合集
  22.     T0=[T0,t0];%t0角的合集
  23.     X=[X,x];%x的合集
  24.     Y=[Y,y0];%y的合集
  25.     o=o+1;
  26. end
  27. plot(theata,X);
复制代码


最佳答案


WarnerChang 发表于 2021-5-21 09:47:31

改好了,你看看吧:
R=10;
r=8;
rc=1;
H=2;
o=0;
theata=[];
T0=[];
X=[];
Y=[];
while o~=91
    syms t0;
    eq=((R*cos(t0)-(r+rc)*sind(o))^2)*(1+(tan(t0)/sind(o))^2)-(rc*sind(o))^2;
    t0=vpasolve(eq,t0,[0,pi/2]);%求出t0的值
    if ~isempty(t0)
        x0=R*cos(t0)/sind(o);%利用t0求出x的值
        y0=R*sin(t0);%利用t0求出y的值
        h=rc.^2-(x0-r-rc).^2;
        L=R+H-y0-h.^0.5+rc*asin((h.^0.5)/rc);%变形长度L
        F=@(t1)integral(@(t)((cos(t)/sind(o)).^2+(sin(t)).^2).^0.5,double(t0),t1)-double(L);
        t1=fsolve(F,0.5);
        if ~isempty(t1)
            x=R*cos(t1)/sind(o);
            theata=[theata,o];%theta角的合集
            T0=[T0,t0];%t0角的合集
            X=[X,x];%x的合集
            Y=[Y,y0];%y的合集
        end
    end
    o=o+1;
end
plot(theata,X);

5 条回复


WarnerChang 发表于 2021-5-18 21:36:05
你的t0都没有解呢

一只小萌新^-^ 发表于 2021-5-20 11:09:21

13行不是给出了t0的解吗?

WarnerChang 发表于 2021-5-20 11:39:34
一只小萌新^-^ 发表于 2021-5-20 11:09
13行不是给出了t0的解吗?

t0=vpasolve(eq,t0,[0,pi/2]);%求出t0的值
这一行解出来的是空值

一只小萌新^-^ 发表于 2021-5-20 12:04:56
WarnerChang 发表于 2021-5-20 11:39
t0=vpasolve(eq,t0,[0,pi/2]);%求出t0的值
这一行解出来的是空值

那我该怎么改呢?

WarnerChang 发表于 2021-5-21 09:41:13

方程没有解,你可以选择不要求解这个方程呢!那样的话倒是好改

WarnerChang 发表于 2021-5-21 09:47:31

改好了,你看看吧:
R=10;
r=8;
rc=1;
H=2;
o=0;
theata=[];
T0=[];
X=[];
Y=[];
while o~=91
    syms t0;
    eq=((R*cos(t0)-(r+rc)*sind(o))^2)*(1+(tan(t0)/sind(o))^2)-(rc*sind(o))^2;
    t0=vpasolve(eq,t0,[0,pi/2]);%求出t0的值
    if ~isempty(t0)
        x0=R*cos(t0)/sind(o);%利用t0求出x的值
        y0=R*sin(t0);%利用t0求出y的值
        h=rc.^2-(x0-r-rc).^2;
        L=R+H-y0-h.^0.5+rc*asin((h.^0.5)/rc);%变形长度L
        F=@(t1)integral(@(t)((cos(t)/sind(o)).^2+(sin(t)).^2).^0.5,double(t0),t1)-double(L);
        t1=fsolve(F,0.5);
        if ~isempty(t1)
            x=R*cos(t1)/sind(o);
            theata=[theata,o];%theta角的合集
            T0=[T0,t0];%t0角的合集
            X=[X,x];%x的合集
            Y=[Y,y0];%y的合集
        end
    end
    o=o+1;
end
plot(theata,X);
回复此楼
您需要登录后才可以回帖 登录 | 注册

本版积分规则

相关帖子
热门教程
站长推荐
快速回复 返回顶部 返回列表