| 查看: 3734 | 回复: 35 | ||||
wgdd木虫 (正式写手)
|
[求助]
求助一个数值积分问题，用matlab的quadgk函数来计算，谢谢！ 已有1人参与
|
|||
|
现求解一个半无限震荡积分问题： 被积函数为s*besselj(0,s*R)./(s.^2-k0^2), 积分区间为[0,inf], 其中R=6, k0=10, 均为已知常数。我打算用两种方法来计算这个问题。(解析解和数值解) 方法一，解析解 对于这样一个积分，可以查到数学手册中积分结果为pi*i/2*besselh(0,1,k0*R). 经过计算可以得到结果为-0.0744。 计算程序 % analytic method k0=10; R=6; p_inc=real(pi*1i/2*besselh(0,1,k0*R)); 方法二，采用matlab中的quadgk函数进行计算，由于k0为奇点，所以将积分区间分为两块。[0,k0]和[k0,inf]，计算结果为0.0390 计算程序 % numerical method k0=10; R=6; p_f=@(s)(s).*besselj(0,s*R)./(s.^2-k0^2); p1=quadgk(p_f,0,k0); p2=quadgk(p_f,k0,inf); p_incN=p1+p2; 求各位大牛支招，为啥两种方法差别这么大呀~正确的数值积分方法应该是怎么样的呢？十分感谢！ |
回复此楼» 猜你喜欢
小木虫没落了，除了祈祷帖子，几乎看不到有价值的帖子
已经有15人回复
-
面上的水很深
已经有18人回复
-
国社科通讯评审
已经有5人回复
-
5个%2F，非固定段22个字母
已经有9人回复
-
生命口C13面上会评时间
已经有11人回复
-
内心匮乏
已经有16人回复
-
生命口面上会评结束了吗大家有知道结果吗？现在还没有结果正常不？
已经有9人回复
-
投英文期刊如何查重
已经有6人回复
-
不要再数国自然申请书的 filecode 的分隔符个数了
已经有21人回复
-
伯胺上甲基求助
已经有7人回复
Mr__Right
专家顾问 (著名写手)
-
专家经验: +31 - 应助: 317 (大学生)
- 金币: 14456.3
- 散金: 500
- 红花: 54
- 帖子: 2716
- 在线: 950.6小时
- 虫号: 1972612
- 注册: 2012-09-04
- 性别: GG
- 专业: 应用数学方法
- 管辖: 程序语言
★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
jjdg: 金币+1, 感谢参与 2017-07-30 09:27:48
wgdd(月只蓝代发): 金币+20, 感谢帮助！ 2017-08-03 12:14:30
jjdg: 金币+1, 感谢参与 2017-07-30 09:27:48
wgdd(月只蓝代发): 金币+20, 感谢帮助！ 2017-08-03 12:14:30
|
如果你查到的解析正确，结果应该是复数-0.07439126817-0.1436835739*I, 所说的 数学手册 具体是哪本？ 此外，你是只要实部吗？ 积分如果得到的是复数，……一般的应用情况下没有继续纠结下去的必要了 这个例子，用matlab里面默认的方法似乎是不行的； 因为被积函数存在s=10的间断点； 对Bessel函数进行积分本来就波动不容易，加上存在被积函数为正负无穷大的间断点，又是在包含无穷大的区间上积分， 这个问题的难度系数：非常高。  |
» 本帖已获得的红花（最新10朵）
2楼2017-07-30 07:19:33
wgdd
木虫 (正式写手)
- 应助: 5 (幼儿园)
- 金币: 4176.2
- 散金: 100
- 红花: 3
- 帖子: 581
- 在线: 252.4小时
- 虫号: 2535934
- 注册: 2013-07-07
- 性别: GG
- 专业: 声学
送红花一朵 |
谢谢您的回复。 这个解析解的公式我是在一本声学手册的附录中找到的。 是这本书：https://www.win.tue.nl/~sjoerdr/papers/boek.pdf 附录中公式(D.69)。 关于这类公式更详尽的证明与推广，也有一篇论文进行了说明：http://aip.scitation.org/doi/full/10.1063/1.3596359 恩，我只要实部，因为算出来的是声压嘛，取实部才有意义~ 您说的积分得到复数，就没有继续纠结下去的必要了，，这个我不太理解哎？是算出来的是复数就不准吗？ 因为考虑到s=10的奇异性，我把积分区间分成了两段[0,k0]和[k0,inf]，我查了下，matlab是推荐采用quadgk时这么处理奇异积分的。 恩恩，Bessel函数在快速震荡，这个积分确实很难求解。。。不知道您有何高见呢？ |
3楼2017-07-30 08:58:37
Mr__Right
专家顾问 (著名写手)
-
专家经验: +31 - 应助: 317 (大学生)
- 金币: 14456.3
- 散金: 500
- 红花: 54
- 帖子: 2716
- 在线: 950.6小时
- 虫号: 1972612
- 注册: 2012-09-04
- 性别: GG
- 专业: 应用数学方法
- 管辖: 程序语言
4楼2017-07-30 13:11:22
wgdd
木虫 (正式写手)
- 应助: 5 (幼儿园)
- 金币: 4176.2
- 散金: 100
- 红花: 3
- 帖子: 581
- 在线: 252.4小时
- 虫号: 2535934
- 注册: 2013-07-07
- 性别: GG
- 专业: 声学
5楼2017-07-30 13:35:20
cooooldog
铁杆木虫 (著名写手)
ส็็็
- 应助: 237 (大学生)
- 金币: 6101.9
- 散金: 1114
- 红花: 39
- 帖子: 1380
- 在线: 553.8小时
- 虫号: 506699
- 注册: 2008-02-18
- 专业: 模式识别
6楼2017-07-30 16:31:37
wgdd
木虫 (正式写手)
- 应助: 5 (幼儿园)
- 金币: 4176.2
- 散金: 100
- 红花: 3
- 帖子: 581
- 在线: 252.4小时
- 虫号: 2535934
- 注册: 2013-07-07
- 性别: GG
- 专业: 声学
7楼2017-07-30 16:35:37
cooooldog
铁杆木虫 (著名写手)
ส็็็
- 应助: 237 (大学生)
- 金币: 6101.9
- 散金: 1114
- 红花: 39
- 帖子: 1380
- 在线: 553.8小时
- 虫号: 506699
- 注册: 2008-02-18
- 专业: 模式识别
8楼2017-07-30 17:14:59
wgdd
木虫 (正式写手)
- 应助: 5 (幼儿园)
- 金币: 4176.2
- 散金: 100
- 红花: 3
- 帖子: 581
- 在线: 252.4小时
- 虫号: 2535934
- 注册: 2013-07-07
- 性别: GG
- 专业: 声学
9楼2017-07-30 17:49:22
 |
10楼2017-08-17 17:24:56
10