24小时热门版块排行榜    

查看: 1768  |  回复: 8
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

jjgy_2002

木虫 (著名写手)

[求助] [求助] 帮看此第一类贝塞尔函数求解代码根据哪篇文献的计算方法编写?

[求助] 帮看此第一类贝塞尔函数求解代码根据哪篇文献的计算方法编写?
代码如下:
function Jx = bessel(n,x)
format long;
if n == 0
    a = [5.7568490574e10;-1.3362590354e10;6.516196407e8;
         -1.121442418e7;7.739233017e4;-1.849052456e2];
    b = [5.7568490411e10;1.029532985e9;9.494680718e6;
         5.927264853e4;2.678532712e2;1.0];
    r =[1.0;-0.1098628627e-2;0.2734510407e-4;
        -0.2073370639e-5;0.2093887211e-6];
    s =[-0.1562499995e-1;0.1430488765e-3;-0.6911147651e-5;
        0.7621095161e-6;-0.934945152e-7];
    if abs(x) < 8
        J0_u = a(1);
        J0_d = b(1);
        for i=2:6
            J0_u = J0_u + a(i)*power(x,2*(i-1));
            J0_d = J0_d + b(i)*power(x,2*(i-1));
        end
        J0 = J0_u/J0_d;
    else
        z = 8/abs(x);
        sita = abs(x) - pi/4;
        R0 = r(1);
        S0 = s(1);
        for i=2:5
            R0 = R0 + r(i)*power(z,2*(i-1));
            S0 = S0 + s(i)*power(z,2*(i-1));
        end
        J0 = (R0*cos(sita)-z*S0*sin(sita))*sqrt(2/pi/abs(x));
    end
    Jx = J0;
else
    if n == 1
        c = [7.2362614232e10;-7.895059235e9;2.423968531e8;
            -2.972611439e6;1.570448260e4;-3.016036606e1];
        d = [1.44725228443e11;2.300535178e9;1.858330474e7;
            9.944743994e4;3.769991397e2;1.0];
        p = [1.0;0.183105e-2;-0.3516396496e-4;
            0.2457520174e-5;-0.240337019e-6];
        q = [0.4687499995e-1;-0.2002690873e-3;
            0.8449199096e-5;-0.88228987e-6;
            0.105787412e-6];
        if abs(x) < 8
            J1_u = c(1);
            J1_d = d(1);
            for i=2:6
                J1_u = J1_u + c(i)*power(x,2*(i-1));
                J1_d = J1_d + d(i)*power(x,2*(i-1));
            end
            J1 = x*J1_u/J1_d;
        else
            z = 8/abs(x);
            fi = abs(x) - 3*pi/4;
            P1 = p(1);
            Q1 = q(1);
            for i=2:5
                P1 = P1 + p(i)*power(z,2*(i-1));
                Q1 = Q1 + q(i)*power(z,2*(i-1));
            end
            J1 = (P1*cos(fi)-z*Q1*sin(fi))*sqrt(2/pi/abs(x));
            if x<0
                J1 = -J1;
            end
        end
        Jx = J1;
    else
        if abs(x) > n
            y = 2/x;
            bess1 = bessel(0,x);
            bess2 = bessel(1,x);
            for j=1:n-1
                Jx = j*y*bess2 - bess1;
                bess1 = bess2;
                bess2 = Jx;
            end
        else
            M = 2*floor(((n + floor(sqrt(40*n))))/2);
            JS = zeros(M+2 ,1);
            JS(M+2) = 0;
            JS(M+1) = 1;
            y = 2/x;
            bSum = 0;
            for k=M:-1:1
                JS(k) = JS(k+1)*y*k - JS(k+2);
            end
            BK = JS(1);
            for k=1:floor(M/2)
                BK = BK + 2*JS(2*k+1);
            end
            Jx = JS(n+1)/BK;
        end
    end
end
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tong_xiao

新虫 (初入文坛)

就是一个visual basic 常用算法集里面有一个此种算法,不过是vb,跟你的matlab大同小异哦
8楼2013-11-19 21:06:47
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 9 个回答

jjgy_2002

木虫 (著名写手)

怎没人响应呢?
3楼2011-06-27 09:58:13
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

vivian5251

新虫 (初入文坛)

【答案】应助回帖

想帮助,但不会,希望有人能帮到你
大家好
4楼2011-07-03 01:10:06
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jjgy_2002

木虫 (著名写手)

引用回帖:
5楼: Originally posted by tong_xiao at 2013-11-12 21:14:14
金币有点少哈

呵呵,能帮忙吗?我多给金币就是了
6楼2013-11-13 13:02:23
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 328求调剂,英语六级551,有科研经历 +3 生物工程调剂 2026-03-17 7/350 2026-03-18 20:41 by Wangjingyue
[考研] 一志愿天津大学化学工艺专业(081702)315分求调剂 +10 yangfz 2026-03-17 10/500 2026-03-18 20:14 by walc
[考研] 材料专硕英一数二306 +4 z1z2z3879 2026-03-18 4/200 2026-03-18 20:00 by 楤哥
[考研] 化工学硕306求调剂 +10 42838695 2026-03-12 10/500 2026-03-18 14:42 by haxia
[考研] 材料与化工一志愿南昌大学327求调剂推荐 +8 Ncdx123456 2026-03-13 9/450 2026-03-18 14:40 by haxia
[考研] 297求调剂 +8 戏精丹丹丹 2026-03-17 8/400 2026-03-18 14:30 by laoshidan
[考研] 296求调剂 +5 大口吃饭 身体健 2026-03-13 5/250 2026-03-17 21:05 by 不惑可乐
[考研] 268求调剂 +8 一定有学上- 2026-03-14 9/450 2026-03-17 17:47 by laoshidan
[考研] 332求调剂 +6 Zz版 2026-03-13 6/300 2026-03-17 17:03 by ruiyingmiao
[考研] 285求调剂 +6 ytter 2026-03-12 6/300 2026-03-16 15:05 by njzyff
[考研] 机械专硕调剂 +3 笨笨兔子 2026-03-12 3/150 2026-03-15 20:02 by 栗子粥?
[考研] 本科南京大学一志愿川大药学327 +3 麦田耕者 2026-03-14 3/150 2026-03-14 20:04 by 外星文明
[考研] 一志愿哈工大材料324分求调剂 +5 闫旭东 2026-03-14 5/250 2026-03-14 14:53 by 木瓜膏
[考研] 311求调剂 +3 冬十三 2026-03-13 3/150 2026-03-13 20:41 by JourneyLucky
[考研] 0703化学求调剂 +7 绿豆芹菜汤 2026-03-12 7/350 2026-03-13 17:25 by njzyff
[考研] 295求调剂 +3 小匕仔汁 2026-03-12 3/150 2026-03-13 15:17 by vgtyfty
[考研] 土木第一志愿276求调剂,科研和技能十分丰富,求新兴方向的导师收留 +3 土木小天才 2026-03-12 3/150 2026-03-13 15:01 by JourneyLucky
[考研] 一志愿山大07化学 332分 四六级已过 本科山东双非 求调剂! +3 不想理你 2026-03-12 3/150 2026-03-13 14:18 by JourneyLucky
[考研] 290求调剂 +3 ADT 2026-03-13 3/150 2026-03-13 10:19 by peike
[考博] 26读博 +4 Rui135246 2026-03-12 10/500 2026-03-13 07:15 by gaobiao
信息提示
请填处理意见