24小时热门版块排行榜    

查看: 1762  |  回复: 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的回帖

jjgy_2002

木虫 (著名写手)

引用回帖:
7楼: Originally posted by tong_xiao at 2013-11-19 21:03:13
我在一本书上看到过Vb编程的,可以吗?...

可以啊
9楼2013-11-20 15:51:13
已阅   回复此楼   关注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的回帖
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 085601专硕,总分342求调剂,地区不限 +3 share_joy 2026-03-16 3/150 2026-03-17 13:41 by houyaoxu
[考研] 0854可跨调剂,一作一项核心论文五项专利,省、国级证书40+数一英一287 +3 小李0854 2026-03-16 3/150 2026-03-17 13:40 by 热情沙漠
[考研] 材料专硕274一志愿陕西师范大学求调剂 +5 薛云鹏 2026-03-13 5/250 2026-03-17 10:15 by Sammy2
[考研] 278求调剂 +3 Yy7400 2026-03-13 3/150 2026-03-17 08:24 by laoshidan
[考研] 286求调剂 +3 lemonzzn 2026-03-16 5/250 2026-03-16 20:43 by lemonzzn
[考研] 0703化学调剂 +6 妮妮ninicgb 2026-03-15 9/450 2026-03-16 16:40 by houyaoxu
[考研] 中科院材料273求调剂 +4 yzydy 2026-03-15 4/200 2026-03-16 15:59 by Gaodh_82
[考研] 0856专硕279求调剂 +5 加油加油!? 2026-03-15 5/250 2026-03-15 11:58 by 2020015
[考研] 294求调剂 +3 Zys010410@ 2026-03-13 4/200 2026-03-15 10:59 by zhq0425
[考研] 297一志愿上交085600求调剂 +5 指尖八千里 2026-03-14 5/250 2026-03-14 17:26 by a不易
[考研] 331求调剂(0703有机化学 +5 ZY-05 2026-03-13 6/300 2026-03-14 10:51 by Jy?
[基金申请] 有必要更换申报口吗 20+3 fannyamoy 2026-03-11 3/150 2026-03-14 00:52 by zhanghaozhu
[考研] 材料工程,326分,求调剂 +6 KRSLSR 2026-03-10 6/300 2026-03-13 23:47 by JourneyLucky
[考研] 复试调剂 +9 Copy267 2026-03-10 9/450 2026-03-13 23:45 by userper
[考研] 0703,333分求调剂 一志愿郑州大学-物理化学 +3 李魔女斗篷 2026-03-11 3/150 2026-03-13 22:24 by JourneyLucky
[考研] (081700)化学工程与技术-298分求调剂 +12 11啦啦啦 2026-03-11 35/1750 2026-03-13 21:25 by JourneyLucky
[考研] 26调剂/材料科学与工程/总分295/求收留 +9 2026调剂侠 2026-03-12 9/450 2026-03-13 20:46 by 18595523086
[考研] 301求调剂 +6 Liyouyumairs 2026-03-11 6/300 2026-03-13 20:11 by JourneyLucky
[考研] 【0856】化学工程(085602)313 分,本科学科评估A类院校化学工程与工艺,诚求调剂 +7 小刘快快上岸 2026-03-11 7/350 2026-03-13 16:06 by ruiyingmiao
[考研] 材料专硕350 求调剂 +4 王金科 2026-03-12 4/200 2026-03-13 16:02 by ruiyingmiao
信息提示
请填处理意见