请教关于解方程画图

zdd13717421

一个比较复杂的方程，自己无法写出方程化简的结果，怎么编程画图
例如
a=0.2;
for x=0.1:0.01:0.4
y=1-exp(a+x.*y)
plot(x,y)
end
在.m文件里面matlab好像不能完成解方程的过程。如果用solve求解的话好像也只能对a，x取一个定值的时候才能求解，x是0.1到0.4的区间的时候求不出来。请教下如何处理才能画出最后的图形来

[ 本帖最后由 zdd13717421 于 2010-1-25 12:00 编辑 ]

messenger

可以改成隐函数形式用ezplot来画

1. 
   
2. a=0.2;
   
3. syms x y
   
4. f=1-exp(a+x*y)-y;
   
5. ezplot(f)
   
6. figure; ezplot(f,[0.1,0.4,-1,1])
   

复制代码

[ 本帖最后由 ChaChing 于 2010-1-26 00:07 编辑 ]

zdd13717421

谢谢楼上的答复。但还是有3个问题不明白
1.但是为什么我在syms x y以后输入隐函数方程以后matlab就提示出错。
   Error in ==> rat at 100
                C = [C*[d;1] C(:,1)];
   Error in ==> sym.sym>symr at 202
       [n,d] = rat(x,tol*abs(x));
   Error in ==> sym.sym at 121
            S(k).s = symr(double(x(k)));
   Error in ==> sym.plus at 8
   A = sym(A);
2.还有，定义的x范围是0.1到0.4，为什么第一幅图x取到负值了？还有第二幅图为什么y的值是一个定值了？
3.如果是两个类似的方程，如
   x=1-exp(-a*(x-0.1));
   y=1-exp(-a*(y-0.3));
   z=y-x;
   此时要画出z与x在x取01.到0.4的时候的图形，那这个时候该怎么写呢？ezplot不能写方程组的吧

messenger

第1个问题可能是因为你把x、y以前定义过为其他的类型；

第2、3个问题请help ezplot

zdd13717421

应该是将z=y-x换成y=x+z代入第二个方程，利用x，y与a的关系
ezplot('log(1-x)/(x+0.1)-log(1-x-z)/(x+z+0.3)',[0,0.4,0,1])这样写吧
但是如果两个式子的指数部分都有一个定义的参量值，如x=1-exp（-a*（b*x+0.1））这样的话，还是不行啊，每次计算我都得自己将b先换成数值才能画图。再次求救

应该能画出来，估计我哪个地方出错了，我再看看

[ 本帖最后由 ChaChing 于 2010-1-25 23:30 编辑 ]

messenger

最简单的办法是将需要画图的函数写到ezplot命令的外面，这样就可能传递参数的值了，比如：
syms x z
a=0.3;
b=0.1;
f=log(1-x)/(x+b)-log(1-x-z)/(x+z+a);
ezplot(f,[0,0.4,0,1])

zdd13717421

谢谢
我最后程序如下
function m(M0,GB)
syms x y
t=2*10^(-9);
q=1.6*10^(-19);
Idm=10^(-13);
deltaT=10^(-5);
td=2*10^(-7);
c=0.01;
Mg=10^8;
n=0.3;
N0=0.3;
Ntr=c*Mg/(1-c);
Ttr=M0/(2*pi*GB);
Ndm1=Idm*t/q;
Ndm2=Idm*M0*Ttr/q;
Pph=1-exp(-N0);
a=Ntr*(exp(t/td)-1)/(exp(deltaT/td)-1);
b=Ntr*(exp(Ttr/td)-1)/(exp(deltaT/td)-1);
f=log(1-y)/(Ndm2+Ndm1+y*(a+b))-log(1-x*0.2592-y)/(Ndm1+Ndm2+(x*0.2592+y)*(a+b)+0.09)
ezplot(f,[0,0.5,10^(-5),10^(-2)])
这样写运行的时候会提示??? Error using ==> y
Too many output arguments.
为什么我在另外一台机器上运行就不会出错了？
望指正

messenger

这段代码在我这里运行正常。

看你对函数的命名方式，function m(M0,GB)，谨慎怀疑你很可能已经命名了一个名为 y 的函数，当你在执行这段代码时，Matlab认为y是一个函数，而不是一个变量，你只要找到这个y.m文件将其删除或重命名即可。

另外，顺便说一句，要养成一个良好的函数命令习惯，不要把函数名称只写成一个字母。

zdd13717421

谢谢，确实如你所言，有命名为y函数的文件。
至于函数命名我以后会注意

再问个问题，如果我要将这个图画成对数图形（y轴为对数），该怎么改

[ 本帖最后由 ChaChing 于 2010-7-14 00:27 编辑 ]

ChaChing

help semilogy

zdd13717421

刚才学习了下，代码改成如下
function zuihou(M0,GB)
t=2*10^(-9);
q=1.6*10^(-19);
Idm=10^(-12);
deltaT=10^(-5);
td=2*10^(-7);
c=0.01;
Mg=10^8;
n=0.3;
N0=0.3;
Ntr=c*Mg/(1-c);
Ttr=M0/(2*pi*GB);
Ndm1=Idm*t/q;
Ndm2=Idm*M0*Ttr/q;
Pph=1-exp(-N0);
a=Ntr*(exp(t/td)-1)/(exp(deltaT/td)-1);
b=Ntr*(exp(Ttr/td)-1)/(exp(deltaT/td)-1);
[email=odefun=@(x,Pd)log(1-Pd)/(Ndm1+Ndm2+Pd*(a+b))-log(1-x*Pph-Pd)/(Ndm1+Ndm2+(x*Pph+Pd)*(a+b)+n*N0]odefun=@(x,Pd)log(1-Pd)/(Ndm1+Ndm2+Pd*(a+b))-log(1-x*Pph-Pd)/(Ndm1+Ndm2+(x*Pph+Pd)*(a+b)+n*N0[/email]);
x0=0;
x1=0.5;
Pd0=2.0097e-018;
[x,Pd]=ode45(odefun,[x0,x1],2.0097e-018);
semilogy(x,Pd)
但是Pd的范围过大了，我只想显示Pd在10^（-7）~10^(-2)内的图形，应该如何修改呢
（在运行的时候M0值为10，GB为3*10^10）

[ 本帖最后由 zdd13717421 于 2010-1-27 15:23 编辑 ]

ChaChing

help axis