zqllc 发表于 2007-6-13 10:57

符号解的问题

帮忙看一下这个程序哪儿有问题,一个超越方程要解一上午,我也不知道哪儿用错了,
而且最后出来是下面的形式,而我本意是要一个数,而不是表达式。是不是定义syms ws用错了?

另外,要在输出文件的每一行坐标前加上两个空格如何操作?谢谢:)

ws = 0.5313
wo =ws+2056510499726355/72057594037927936+7/36*pi
x =61379602534499038979226340382253/2535301200456458802993406410752*sin(ws+2056510499726355/72057594037927936+7/3

代码:

%this script is to calculate the face-gear face data
%清屏,清理内存
clear;
clc;
%start to input paramenter
Ns=20;
m=3;
as=20;
r=Ns*m*cos(as);
q2s=0.15;
w2=q2s;
N2=Ns/q2s;
os0=pi/(2*Ns)-inv(as);
syms ws;
%变量
y=-240;
os=0;
data=fopen('data.dat','w');
for i=0:7
os=0+5*i*pi/180;
wo=ws+os0+os;
f=r*(sin(q2s*ws)*(sin(ws+os0+os)-os*cos(wo))+N2*cos(q2s*ws)/(Ns*cos(wo)))-240;
ws=solve(f);
x=(r*(cos(w2)*(sin(wo)-os*cos(wo))-N2*sin(w2)/(Ns*cos(wo))));
z=(-r*(cos(wo)+os*sin(wo)));
%测试输出用
%x=10;
%z=20;
fprintf(data,'%6.4f %6.4f %6.4f\n',x,y,z);
end
fclose(data);

eight 发表于 2007-6-13 11:06

原帖由 zqllc 于 2007-6-13 10:57 发表 http://www.chinavib.com/forum/images/common/back.gif
帮忙看一下这个程序哪儿有问题,一个超越方程要解一上午,我也不知道哪儿用错了,
而且最后出来是下面的形式,而我本意是要一个数,而不是表达式。是不是定义syms ws用错了?

另外,要在输出文件的每一行坐 ...


fprintf(data,'%6.4f %6.4f %6.4f\n',x,y,z); 不行吗?试试吧,不知道

zqllc 发表于 2007-6-13 19:45

谢谢,这样加空格是可以的。
不过上面出表达式的问题是哪儿写错了能指点一下吗?:)

zqllc 发表于 2007-6-16 10:19

问题又太弱智了?:'(

eight 发表于 2007-6-16 10:59

原帖由 zqllc 于 2007-6-16 10:19 发表 http://www.chinavib.com/forum/images/common/back.gif
问题又太弱智了?:'(

help subs

zqllc 发表于 2007-6-16 20:05

回复 #5 eight 的帖子

谢谢
改成如下后:
%this script is to calculate the face-gear face data
%清屏,清理内存
clear;
clc;
%start to input paramenter
Ns=20;
m=3;
as=20;
r=Ns*m*cos(as);
q2s=0.15;
w2=q2s;
N2=Ns/q2s;
os0=pi/(2*Ns)-inv(as);
syms ws x z;
%变量
y=-240;
os=0;
data=fopen('data.dat','wt');
for y=240:10:300;
    fprintf(data,'%s %6.4f\n','when y =',y);
    for i=0:7
      os=0+5*i*pi/180;
      wo=ws+os0+os;
      fws=r*(sin(q2s*ws)*(sin(wo)-os*cos(wo))+cos(q2s*ws)/(q2s*cos(wo)))+y;
      subs(fws);
      %ws=solve(fws);
      ws=25;
      fx=x-r*(cos(w2)*(sin(wo)-os*cos(wo))-N2*sin(w2)/(Ns*cos(wo)));
      fz=z+r*(cos(wo)+os*sin(wo));
      subs(fx);
      subs(fz);
      solve(fx,x);
      solve(fz,z);
      fprintf(data,'%6.4f %6.4f %6.4f %12.8f\n',x,y,z,ws);
      %测试输出用
    end
end
fclose(data);
%M=load('data.dat');
%code end
fws =
6891893330697659/281474976710656*sin(3/20*ws)*sin(ws+2056510499726355/72057594037927936)+34459466653488295/211106232532992*cos(3/20*ws)/cos(ws+2056510499726355/72057594037927936)+240

报错如下,应该也是x和z仍然不是数值,fws,fx,fz已经代入,方程能解出ws。
??? Error using ==> fprintf
Function is not defined for 'sym' inputs.
Error in ==> facegear at 34
      fprintf(data,'%6.4f %6.4f %6.4f %12.8f\n',x,y,z,ws);

xjzuo 发表于 2007-6-17 09:11

程序感觉很糟糕:
1. fprintf(data,'%s %6.4f\n','when y =',y);-----分号要去掉;
2. subs(fx);subs(fz);----应该改为 fx=subs(fx); fz=subs(fz);
3. ...
自己再仔细检查一下吧,...

zqllc 发表于 2007-6-17 10:25

回复 #7 xjzuo 的帖子

谢谢。以前没用matlab编过程序,见笑:@(
1、如果去掉分号,commond窗口太多。
2、subs改过来了。

改成:
%this script is to calculate the face-gear face data
%清屏,清理内存
clear;
clc;
%start to input paramenter
syms ws;
Ns=25;            
m=3;               
q2s=1/6;         
as=20;            
r=Ns*m*cos(as);   
N2=Ns/q2s;         
w2=q2s*ws;
os0=pi/(2*Ns)-inv(as);
%变量
y=-240;
os=0;
data=fopen('point.dat','wt');
rotate=fopen('rotateangle.dat','wt');
fprintf(rotate,'%s\n','    Ns,    N2,    r,    wo,    w2,    ws');
for y=-290:10:-200;
    fprintf(data,'%s %6.4f\n','when y =',y);
    fprintf(rotate,'%s %6.4f\n','when y =',y);
    for i=0:1:7
      os=0+5*i*pi/180;
      wo=ws+os0+os;
      fws=r*(sin(w2)*(sin(wo)-os*cos(wo))+cos(w2)/(q2s*cos(wo)))+y;
      fws=subs(fws);
      wss=subs(solve(fws,'ws'));                %wss即是ws
      wss=wss(imag(wss)==0&real(wss)>0);      %只要大于零的实根
      %wss=0.253;                               %测试用
      wo=subs(wss+os0+os);
      w2=subs(q2s*wss);
      x=subs(r*(cos(w2)*(sin(wo)-os*cos(wo))-N2*sin(w2)/(Ns*cos(wo))));
      z=subs(r*(cos(wo)+os*sin(wo)));
      fprintf(data,' %6.4f %6.4f %6.4f\n',x,y,z);
      fprintf(rotate,'   %6.4f %6.4f %12.8f %12.8f %12.8f %12.8f\n',Ns,N2,r,wo,w2,wss);
      %测试输出用
    end
end
fclose(data);
%load point.dat content
points = textread('point.dat','%s','delimiter','\n','whitespace','')

这样就可以输出了。但速度还是不快。solve函数运算时间很长。

[ 本帖最后由 zqllc 于 2007-6-17 15:56 编辑 ]
页: [1]
查看完整版本: 符号解的问题