24小时热门版块排行榜    

Znn3bq.jpeg
查看: 381  |  回复: 3
当前主题已经存档。

jianchaoyv2651

铁虫 (初入文坛)

[交流] 【求助】fortran程序出错了

一个计算积分的程序如下:
program integral_1
integer,parameter::nmax=12
real*8,parameter::pio2=6.2831852
external Func
real*8::A=0.0,B=pio2,s
write(*,'(1x,a)') 'integral of Func with 2^(n-1) points'
write(*,'(1x,t7,a,t16,a)')'n','Approx.integral'
do i=1,nmax
   call TRAPZD(Func,A,B,s,i)
   write(*,'(1x,i6,f12.6)')i,s
end do
end program integral_1


Function Func(x)
real*8::x
real*8::Func
    Func=sin(x)**4/x**3
end function Func



subroutine TRAPZD(Func,a,b,s,n)
integer,parameter::k1=selected_int_kind(9)
integer(kind=k1)::tnm
real*8::del,sum,x,a,b,s
integer::n
if(n==1)then
  s=0.5*(b-a)*(Func(a)+Func(b))
  else
  tnm=2**(n-1)
  del=(b-a)/tnm
  x=a
  sum=0.0
  do j=2,tnm
    x=x+del
    sum=sum+Func(x)
  end do
   s=0.5*(Func(a)+Func(b)+2.0*sum)*del
end if
end subroutine TRAPZD
运行结果如下:
integral of Func with 2^(n-1) points
      n        Approx.integral
      1 NaN
      2 NaN
      3 NaN
      4 NaN
      5 NaN
      6 NaN
      7 NaN
      8 NaN
      9 NaN
     10 NaN
     11 NaN
     12 NaN
Press any key to continue
我感觉是计算范围越界,请高手指点!!
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

老虎大王

木虫 (著名写手)


jianchaoyv2651(金币+1,VIP+0): 4-13 11:31
呵呵呵。小兄弟,明显的问题。在TRAPZD的最后,有 s=0.5*(Func(a)+Func(b)+2.0*sum)*del,但是a是0啊。    Func=sin(x)**4/x**3,这样,Func(a)就是无穷大喽。。。。

遇到输出NaN,而没有其它出错信息,首先要考虑是不是被零除了。
2楼2009-04-13 10:59:02
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jianchaoyv

金虫 (小有名气)

★ ★ ★ ★
jianchaoyv2651(金币+4,VIP+0): 4-13 11:31
怎样解决奇异点?
3楼2009-04-13 11:21:10
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jianchaoyv2651

铁虫 (初入文坛)

知道了,算法中有反常积分。谢谢了!
4楼2009-04-13 11:32:08
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 jianchaoyv2651 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 评审有感 +10 popular289 2026-05-18 15/750 2026-05-19 09:18 by jurkat.1640
[基金申请] 同样的基金本子,换个专家直接从C变A! +3 国自然国社科中 2026-05-19 3/150 2026-05-19 08:50 by Equinoxhua
[基金申请] 别被青基扩招骗了!26年科研内卷才刚刚开始 +3 国自然国社科中 2026-05-14 4/200 2026-05-19 08:48 by archvillain
[考博] 云南大学材料与能源学院解琳课题组钙钛矿博士招生 +3 光伏爱好者 2026-05-17 4/200 2026-05-19 01:19 by 风起于浮萍
[考博] 博士申请 +5 星…… 2026-05-18 6/300 2026-05-18 23:49 by 糊糊涂涂好
[基金申请] 面上本子正文33页,违规吗?会被低分嘛? +8 1234567wang 2026-05-17 10/500 2026-05-18 18:52 by zzahkj
[基金申请] 今年审到国自然15份,谈谈感受 +16 国自然国社科中 2026-05-17 16/800 2026-05-18 14:58 by gy116024
[基金申请] 重磅!青年科学基金项目(C类)资助增幅预计超过50% +7 水和泥不是水泥 2026-05-13 10/500 2026-05-18 07:50 by 水和泥不是水泥
[找工作] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +5 l7k6xnh0yc 2026-05-14 5/250 2026-05-17 19:39 by Equinoxhua
[公派出国] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +4 v9tggjlwd0 2026-05-15 5/250 2026-05-17 18:39 by Equinoxhua
[考博] 光量子物理方向 博士招生 1人(2026.09) +3 sandyworld 2026-05-15 4/200 2026-05-17 14:38 by sandyworld
[考博] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3 7hlccowb3h 2026-05-15 4/200 2026-05-17 07:46 by 11n4dfd8yn
[找工作] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3 l7k6xnh0yc 2026-05-14 6/300 2026-05-17 07:31 by 11n4dfd8yn
[找工作] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3 k37jurhrau 2026-05-16 3/150 2026-05-17 01:37 by ue3ir18jc3
[考研] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3 k37jurhrau 2026-05-16 3/150 2026-05-17 01:25 by ue3ir18jc3
[考研] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3 ky2p12rrjj 2026-05-15 4/200 2026-05-17 00:50 by ue3ir18jc3
[高分子] 本人最近太闲了,谁有问题可以提,每天会统一回复 +9 一切都是空工 2026-05-12 20/1000 2026-05-16 19:52 by Equinoxhua
[有机交流] 求有机合成大神指点三硫酸乙烯酯(CAS:2793408-99-6)的合成路线 30+3 Leekmid 2026-05-13 10/500 2026-05-16 16:37 by czyzsu
[考研] 售SCI一区T0P文章,我:8.O.5.5.1.O.5.4,科目齐全,可+急 +4 l7k6xnh0yc 2026-05-14 6/300 2026-05-16 11:29 by h3oerqvkv9
[教师之家] 教学课件你会给同学吗 +8 硕士研究生吗 2026-05-13 8/400 2026-05-14 22:23 by 常规沥青
信息提示
请填处理意见