超越方程求解中重根的去除问题（问题已经解决）
本帖最后由 zhoulon 于 2009-6-1 13:24 编辑现在有一个超越方程需要求解，为了验证fsolve命令的用法以及方程根的正确性，
专门找了一本书上附录里面的下面简单超越方程，进而进行求解对比，最后fsolve
求解出来的结果是对的，但是存在一个问题：求解出来好多重根，想问一下如何去除
这些重根或者是重根的显示？
matlab程序：
function y=fun(x)
y=x.*tan(x);
=fsolve(@fun,)
求解出来的解：
x =
Columns 1 through 9
0 0.0011 3.1416 3.1416 3.1416 6.2832 6.2832 6.2832 9.4248
Columns 10 through 16
9.4248 9.4248 12.5664 12.5664 12.5664 12.5664 15.7080
eval =
1.0e-005 *
Columns 1 through 9
0 0.1273 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000
Columns 10 through 16
-0.0000 -0.0000 0.0029 -0.0000 -0.0000 -0.0000 -0.0000
具体公式： 本帖最后由 zhoulon 于 2009-5-31 14:16 编辑
说明一下：我想求在0～5*pi之间的根，本来应该是6个根，但是从上面的计算结果来看有16个根，其中有几个应该算是重根吧（例如里面有好几个3.1416和6.2832），咋样去除呀，在此表示感谢。
另外，我试过options中调整eval的最大值的显示范围，可以将根调整到11个，也就是去除了5个，其它的就没有办法了。
其中options选项调整为：options=optimset( 'MaxIter',30000, 'MaxFunEvals',1.0e-6) unique(x) 首先表示感谢，但出来还是有一个重根
ans =
0 0.0011 3.1416 6.2832 9.4248 12.5664 12.5664 15.7080 你的数据显示有问题，没有把有效数字位数都显示出来，其实那两个显示相同的数是不一样的，采用long G，看到更多的有效位数就明白了。
ps. 其实，你这不是重根，只是局部收敛到同一个根。 另外一个上面计算出来的0.0011本身应该不是这个方程的一个根，但是却出现在计算结果里面，这个应该只能用我上面的options中的'MaxFunEvals',1.0e-6来调整吧，还有没有其它办法？ 5# messenger
好的，感谢，我也感觉应该是显示问题，这样就确认了，呵呵。 嗯，可以通过上面的options中的'MaxFunEvals',1.0e-6来调整，但不一定总是有效。
数值解就这样，数值解要根据物理意义来选择，
比如，0.0011这个解，数值计算中可能比其它解还靠近0，从数学上无法消除
另外，怀疑你的那个12.5664 12.5664 其实是一个根，你可以选择一下有效位数，如果计算精度要求不高，设置一个允许误差，将允许误差之内的根算一个根。 8# messenger
呵呵，非常感谢，这样问题就算解决的圆满了。 接着昨天的那个方程再提一问：
这两天在论坛搜索了一下以前的帖子，知道可以借助图形来大致了解一下方程根的分布情况（例如在一定范围内根的个数、根的大致取值大小等）。于是，我设置了一个函数：y=x.*tan(x)对自变量x绘制（0,5*pi）范围内的图象，结果不理想，这是我们都可以想到的，因为tan(x)当x=(2*n+1)*(pi/2)，（其中n是整数）的时候，是无穷大的，导致这个函数绘图出来的图形没有达到鉴定根的分布情况的效果。但我们知道tan(x)是以pi为周期的，难道只能按照pi这个周期一个周期一个周期进行绘制吗？例如，我要看在（0,5*pi）内的图象，就要绘制（0,pi）、（pi,2*pi）...的图象吗？这样是不是太麻烦了，呵呵。
如果非要一块绘制x在（0,5*pi）内的图象，是否可以这样来做，每当遇到x=(2*n+1)*(pi/2)，（其中n是整数）的时候，我们可以在取采样点的时候，不要取到x=(2*n+1)*(pi/2)这个点，可以用其它相近的点来代替这个点，但是这样跟实际有点不太相符，不过鉴于图形的幅面（不可能反映自变量和因变量在-inf到inf的所有值），我们所取的只是有限范围内的数值，这样是否可行？或者有其它更好的办法？ 借助此贴再提两个问题：
1）我现在要求一个无穷级数求和的问题，这个无穷级数里面有我在求解超越方程出来的无穷个特征值，但我们知道不可能让它无穷的计算下去，看了论坛以前的帖子，我知道在对函数在无穷区间（例如0到inf）上积分的时候，可以先设置一个很大的积分上限M，再设定一定的相对误差，计算积分函数在（0,M）和（M,2*M）上的积分，如果后者在整个积分中所占的比重小于我们设定的相对误差，那么就可以用这个积分函数在（0,M）的积分来代替（0,inf）的积分，如果不满足，就用2*M来代替M接计算，知道满足设定的相对误差为止。那么，对于我上面的无穷级数求和的问题是否也可以这样做：先取定M和2*M，计算在（1,M）和（M,2*M）范围内级数的和，同时象上面一样设定一个相对误差，如果后者在整个求和中所占的比重小于我们设定的相对误差，那么就可以用这个级数在（0,M）的和来代替（0,inf）的和，如果不满足，就用2*M来代替M接计算，知道满足设定的相对误差为止。个人感觉这样是可行的，因为无论是积分还是求级数和，都是一个逐渐逼近的过程，也都是一个求取极限的过程，不知道大家以为如何？
2）其实本身来说给定的超越函数应该是：x.*tan(x)=c，其中c是一个常数，我为了方便，把c取为0了，但是那本书的附录里面有一个c＝inf的解的情况。我想问一下对于这个情况是否也可以象上面1）中我说的求取（0,inf）上积分的方法，设定两个M和2*M以及一定的相对误差，看在这两个M条件下计算出来的x的数值是否满足这个相对误差，如果满足计算停止，如果不能满足，就把2*M赋值给M重新进行计算，重新进行判断，知道满足设定的相对误差为止。 11# zhoulon
看了下面这个以前的帖子，问题1）的解决方法正如我上面所说
http://forum.simwe.com/viewthread.php?tid=223556&highlight=%CE%DE%C7%EE%B4%F3
在此求另外两个问题的解决方法，在此提前表示感谢。 10# zhoulon
还是没明白为什么有无穷大的值就“导致这个函数绘图出来的图形没有达到鉴定根的分布情况的效果”，不过，如果你想去除无穷大值的干扰，可以用ylim命令限制一下函数的值域，使函数的值在0附近，这样更容易看一些。比如ylim([-50 50]) 11# zhoulon
2）答：应该不可以吧，这和第1）个问题完全是两回事 本帖最后由 zhoulon 于 2009-6-1 11:07 编辑
这个是绘制的0～5*pi之间图形的代码
>> x=0:pi/100:5*pi;
>> y=x.*tan(x);
>> plot(x,y)
下面是图形，可以看出来由于tan（x）在x=(2*n+1)*(pi/2)，（其中n是整数）
时的无穷大，导致图形上面显示不出来根的具体情况，而只能看到几个尖点和一条
y=的横线。 下面时ylim([-50 50])的结果，看起来是可行的，如果把y的值域限制的更小一点的话，
效果就更好了，感谢。 15# zhoulon
试试下面的代码：
x=0:pi/100:5*pi;
y=x.*tan(x);
plot(x,y,'red','LineWidth',1.5)
grid minor;
ylim([-50 50])
14# messenger
感谢一下，那么这个问题该如何来处理，那个附录上面的结果又是如何来计算出来的？ 18# zhoulon
想了一会才发现，这个问题真忽悠人~~~
我说以前怎么没遇到这种情况的方程求解呢
这种情况，一般方程两边取倒数，变为求解函数为零的情况。
不知道那本书为什么会把这个问题拿出来，真唬人，没留意的还真被他唬住了，让人觉得作者水平挺高的。 19# messenger
思维挺敏捷呀，我取了倒数以后的确问题解决了，计算结果总算是对上了，非常感谢。
我就是那个被唬的一伙啊，呵呵。
说明一下：我那本书附录上面只是给了一个c和x计算出来的数值的表格，计算过程根本就没有给。
页:
[1]