24小时热门版块排行榜    

查看: 1223  |  回复: 3

herozhu

木虫 (正式写手)

[求助] 用symsum计算发散级数和无法转换为数值型

如题,需要计算一个发散级数的有限项和,级数喂,n取值范围为[1,1e18]。采用symsum计算。当n取值范围较小时,计算结果直接可用double转换为数值型带入后续计算,但当n取值范围很大时,此时的计算结果无法转换为数值型。
>> syms n x
>> ff=n^2*exp(-10*(1-1/n^2))

ff =

n^2*exp(10/n^2 - 10)

>> symsum(ff,n,1,1e18)

ans =

sum(n^2*exp(10/n^2 - 10), n == 1..1000000000000000000)

用double(ans),出现错误提示:
Error using mupadmex
Error in MuPAD command: DOUBLE cannot convert the input expression
into a double array.
If the input expression contains a symbolic variable, use the VPA
function instead.

Error in sym/double (line 827)
        Xstr = mupadmex('symobj::double', S.s, 0);


尝试用vpa转换,
>> vpa(ans)

ans =

numeric::sum(n^2*exp(10.0/n^2 - 10.0), n == 1..1000000000000000000)
得到的结果尝试进行后续数值计算,比如ans/1e5:

>> ans/1e5

ans =

numeric::sum(n^2*exp(10.0/n^2 - 10.0), n == 1..1000000000000000000)/100000


想求教是否有方法将发散级数大区间的求和结果用数值型表示?并且能够保留合适精度,带入后续计算过程?
不胜感谢!
回复此楼
最怕的即是让爱自己的人失望
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

hytao2012

铁杆木虫 (正式写手)

木头虫子

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
herozhu(nono2009代发): 金币+10 2015-09-09 06:04:44
貌似循环次数太多,建议分开算,分开也够多的。下面是我的想法
clear; clc;

uplim=1e6;
lowlim=1;

cycltime=uplim/1e5;
syms n
result=0;
for i=1:cycltime
    temp1=symsum(n^2*exp(10/n^2)/exp(10),n,(i-1)*1e5+1,i*1e5);
    temp2=vpa(temp1);
    result=result+temp2;
end
disp(result)
2楼2014-08-25 11:33:47
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

herozhu

木虫 (正式写手)

引用回帖:
2楼: Originally posted by hytao2012 at 2014-08-25 11:33:47
貌似循环次数太多,建议分开算,分开也够多的。下面是我的想法
clear; clc;

uplim=1e6;
lowlim=1;

cycltime=uplim/1e5;
syms n
result=0;
for i=1:cycltime
    temp1=symsum(n^2*exp(10/n^2)/exp(10), ...

多谢回复!

思路挺好的,大于1e5时vpa即不显示数值,单用兄台的方法可以循环计算出1e7对应的数值。正如你所述,循环太多也是个问题,1e7即耗时132s,我的1e18量级不知道得算到啥时去啊。

萨哈方程中的这个partition function肯定得有办法求的啊,再想想再想想
最怕的即是让爱自己的人失望
3楼2014-08-25 14:44:04
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

hytao2012

铁杆木虫 (正式写手)

木头虫子

引用回帖:
3楼: Originally posted by herozhu at 2014-08-25 14:44:04
多谢回复!

思路挺好的,大于1e5时vpa即不显示数值,单用兄台的方法可以循环计算出1e7对应的数值。正如你所述,循环太多也是个问题,1e7即耗时132s,我的1e18量级不知道得算到啥时去啊。

萨哈方程中的这个pa ...

嗯,不知道为什么循环过多或大于1e5就不行,是Matlab内存设置问题之类的吗?还是vpa函数自身限制,可以通过修改vpa这个函数解决吗?
也可以试试用其他办法解出来。
实在不行,把循环也分成若干次。。。就是略麻烦
4楼2014-08-25 15:56:31
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 herozhu 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考博] 澳大利亚 Murdoch University 全奖博士招生(3个名额)地质化工冶金领域 +38 AI8RaGaPaSCR 2026-08-07 40/2000 2026-08-13 23:45 by 9kfQ6ZJVaWOC
[基金申请] 关于Filecode分析方法 +9 majunge000 2026-08-10 12/600 2026-08-13 23:20 by iwuli
[基金申请] 咱们一起用铁证分析2026国家社科基金中标与否 +3 启萌科技 2026-08-12 11/550 2026-08-13 22:57 by 启萌科技
[基金申请] 静等基金结果 +8 gjjjzhong 2026-08-10 21/1050 2026-08-13 17:56 by 且听虎啸
[基金申请] 重要来源:本周末出结果 +10 瞬息宇宙 2026-08-12 10/500 2026-08-13 15:46 by likettle
[基金申请] Filecode 又变了,巨变 +3 WH3796 2026-08-12 4/200 2026-08-13 14:13 by 小木虫6752397
[基金申请] 我的国基提前知道中了,可是同事的操作让我实在接受不了,怎么会有这样的人 +9 家与远方 2026-08-10 14/700 2026-08-13 09:41 by 臭臭不臭01
[基金申请] 结合人工智能,周易传统文化,filecode打分制来了,3分以上希望很大。 +3 Tide man 2026-08-12 4/200 2026-08-13 08:35 by ZJTJZ
[基金申请] 好奇怪的filecode +5 布布和一二 2026-08-08 6/300 2026-08-12 16:11 by 云上清扬
[论文投稿] 职称评审,求友友推荐见刊最快的期刊 +4 工厂打螺丝 2026-08-08 4/200 2026-08-12 09:18 by hansi2025
[基金申请] 奇怪,两个人的filecode固定段从头到尾一模一样 +7 布布和一二 2026-08-10 10/500 2026-08-12 09:12 by 布布和一二
[基金申请] 是这周出结果还是下周出结果? +3 yuleib84 2026-08-11 3/150 2026-08-11 21:53 by jnhyjjm
[硕博家园] 读博的好处 +3 lnee 2026-08-11 3/150 2026-08-11 18:10 by 希望我好好的
[基金申请] 确定了,国自然21号放榜 +6 布布和一二 2026-08-10 7/350 2026-08-10 19:15 by 2000zf36392
[基金申请] 2026国自然放榜时间 +9 布布和一二 2026-08-08 9/450 2026-08-10 11:22 by xxxx2020
[基金申请] 面上项目filecode邪修 +5 西山十月 2026-08-09 7/350 2026-08-10 07:32 by 仁砚薪传
[基金申请] fileCode有新解读? +10 Tide man 2026-08-08 18/900 2026-08-09 12:55 by 仁砚薪传
[基金申请] 关于filecode,很负责任的告诉大家 +6 爱看书的可乐 2026-08-08 7/350 2026-08-08 22:13 by a_niu
[基金申请] 国基金的申报应该改成非等额制,评价高的钱多评价低的钱少,但是增加资助率 +7 a089 2026-08-07 7/350 2026-08-08 18:05 by gltch
[基金申请] 化学口download_prp&fileCode的固定段好像这几天一直没变,有变的大神么? +3 Tide man 2026-08-07 4/200 2026-08-07 22:39 by Tide man
信息提示
请填处理意见