24小时热门版块排行榜    

查看: 1379  |  回复: 7
当前主题已经存档。

zhangyatao

木虫 (正式写手)

[交流] [color=Blue]【求助】Bessel函数的Matlab解法[/color]

对于下面的Bessel函数,如何用Matlab求解,麻烦给出源程序或类似程序,多谢了!
D*x*J1(x)=J0(x),其中D是常数。

[ Last edited by zhangyatao on 2009-4-1 at 16:51 ]
回复此楼

» 猜你喜欢

资源共享,一起飞翔!
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

fspdlh

金虫 (正式写手)

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
zhangyatao(金币+8,VIP+0):程序有点错误,不过还是挺感谢你的! 4-1 15:41
sunxiao(金币+2,VIP+0):谢谢参与交流,欢迎常来仿真编程版 4-1 22:38
可以先画图像看看啊
fplot('x*Bessel(x,1)-Bessel(x,0)',[-10,100])
可以看出,方程无解,但函数随x的增大是趋近于零的
>> x=20

x =

    20

>> x*Bessel(x,1)-Bessel(x,0)

ans =

  7.7470e-024

>> x=10000

x =

       10000

>> x*Bessel(x,1)-Bessel(x,0)

ans =

     0

[ Last edited by fspdlh on 2009-4-1 at 15:19 ]
2楼2009-04-01 15:16:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

zhangyatao

木虫 (正式写手)

fplot('x*Bessel(x,1)-Bessel(x,0)',[-10,100])
应该是:
fplot('x*Bessel(1,x)-Bessel(0,x)',[-10,100]) 吧
资源共享,一起飞翔!
3楼2009-04-01 15:24:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

戴氓

新虫 (初入文坛)

★ ★ ★ ★ ★
zhangyatao(金币+5,VIP+0):十分感谢! 4-1 15:46
fplot('x*Bessel(1,x)-Bessel(0,x)',[-10,100])
有解的,但我没法找到零点
4楼2009-04-01 15:34:29
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

fspdlh

金虫 (正式写手)

★ ★
zhangyatao(金币+2,VIP+0):谢谢你! 4-1 15:42
呵呵,我不是学物理的,不太明白,这样是个波动的一个函数,好像周期是固定的啊,可以在每一段上用二分法来做。
5楼2009-04-01 15:34:46
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

zhangyatao

木虫 (正式写手)

呵呵,是的,就是不知道如何找到零点!
资源共享,一起飞翔!
6楼2009-04-01 15:37:28
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

fspdlh

金虫 (正式写手)

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
zhangyatao(金币+15,VIP+0):十分感谢!另外,如果方程里面有个常数D(见编辑后的帖子),还能否得到带有常数D的数值解,谢谢哈! 4-1 16:45
sunxiao(金币+4,VIP+0):谢谢参与交流,欢迎常来仿真编程版 4-1 22:37
function x=solvefun(a,b,tol)

%--------------------------------------------------------------------------
range=a:1:b;
if range(end)     range=[range,b];
end
tag=[abs(diff(sign(fun(range)))) 0];
range=range(find(tag>0))';
range=[range range+1];
%--------------------------------------------------------------------------
n=size(range,1);
tol=tol/10;
x=zeros(n,1);
for i=1:n
    xmin=range(i,1);
    xmax=range(i,2);
    xcur=(xmin+xmax)/2;
    while xmax-xmin>tol
        if sign(fun(xmin))*sign(fun(xcur))>0
            xmin=xcur;
        else
            xmax=xcur;
        end
        xcur=(xmin+xmax)/2;
    end
    x(i)=xcur;
end
%--------------------------------------------------------------------------
function y=fun(x)
y=x.*Bessel(1,x)-Bessel(0,x);


>> solvefun(0,100,1e-4)

ans =

    1.2558
    4.0795
    7.1558
   10.2710
   13.3984
   16.5312
   19.6667
   22.8040
   25.9422
   29.0812
   32.2207
   35.3606
   38.5007
   41.6411
   44.7817
   47.9223
   51.0631
   54.2040
   57.3450
   60.4860
   63.6271
   66.7682
   69.9094
   73.0506
   76.1918
   79.3331
   82.4744
   85.6157
   88.7570
   91.8984
   95.0398
   98.1811

>>
7楼2009-04-01 16:20:38
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

fspdlh

金虫 (正式写手)

★ ★ ★ ★
sunxiao(金币+4,VIP+0):谢谢参与交流,欢迎常来仿真编程版 4-1 22:37
function x=solvefun(a,b,D,tol)

%--------------------------------------------------------------------------
range=a:1:b;
if range(end)     range=[range,b];
end
tag=[abs(diff(sign(fun(range,D)))) 0];
range=range(find(tag>0))';
range=[range range+1];
%--------------------------------------------------------------------------
n=size(range,1);
tol=tol/10;
x=zeros(n,1);
for i=1:n
    xmin=range(i,1);
    xmax=range(i,2);
    xcur=(xmin+xmax)/2;
    while xmax-xmin>tol
        if sign(fun(xmin,D))*sign(fun(xcur,D))>0
            xmin=xcur;
        else
            xmax=xcur;
        end
        xcur=(xmin+xmax)/2;
    end
    x(i)=xcur;
end
%--------------------------------------------------------------------------
function y=fun(x,D)
y=D.*x.*Bessel(1,x)-Bessel(0,x);
8楼2009-04-01 16:55:01
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 zhangyatao 的主题更新
普通表情 高级回复 (可上传附件)
信息提示
请填处理意见