Fortran Coder

查看: 4208|回复: 3
打印 上一主题 下一主题

[特殊函数] 还是同一个题,但是我换成勒让德函数就总出错,大神能.....

[复制链接]

15

帖子

6

主题

0

精华

入门

F 币
73 元
贡献
45 点
跳转到指定楼层
楼主
发表于 2016-4-20 19:23:05 | 只看该作者 回帖奖励 |正序浏览 |阅读模式
[Fortran] 纯文本查看 复制代码
module autoxinpusen
contains 
subroutine solve(func,s,a,b,tol,m)
implicit none
integer::m,i
external func
double precision s,a,b,tol,s1,del
m=10
do i=1,m/2.d0
call simp(func,s,a,b,m)
m=m*2
call simp(func,s1,a,b,m)
del=abs(s-s1)
if(del<=tol)  exit
end do
s=s1
end subroutine solve

subroutine simp(func,s,a,b,m)
implicit none
integer::m,k,n
external func
double precision s,a,b,t1,t2,f1,f2,h,f3,f4
s=0d0
h=(b-a)/(m*2d0)
call func(f1,a,n)
call func(f2,b,n)
s=f1+f2
!k=0 情况
call func(f1,a+h,n)
s=s+4d0*f1
do k=1,m-1
t1=a+(2d0*k+1)*h
t2=a+2d0*k*h
call func(f3,t1,n)
call func(f4,t2,n)
s=s+f3*4d0+f4*2d0  
end do
s=s*h/3d0
end subroutine simp

subroutine fun1(f,x,n)
implicit none
integer k,n
double precision f,x,p(0:n)

 p(0) = 1.d0 
 p(1) = x 
 do k = 2, n 
    p(k)=((2.d0*k-1.d0)*x*p(k-1)-(k-1.d0)*p(k-2))/k 
 end do 
 f=p(n) 
end subroutine fun1
end module autoxinpusen

program main
use autoxinpusen
implicit none
integer::n
double precision a,b,abs,s
write(6,*)'n'
read(5,*)n
write(6,*)'a'
read(5,*)a
write(6,*)'b'
read(5,*)b
write(6,*)'abs'
read(5,*)abs
open(unit=11,file='result.txt')
call solve(fun1,s,a,b,abs,n)
write(6,'(1x,i3,1x,d15.8)')n,s
write(11,*)n,s
pause
end program main
!  --------------------------------------------------
!  Silverfrost FTN95 for Microsoft Visual Studio
!  Free Format FTN95 Source File
!  --------------------------------------------------


分享到:  微信微信
收藏收藏 点赞点赞 点踩点踩

15

帖子

6

主题

0

精华

入门

F 币
73 元
贡献
45 点
地板
 楼主| 发表于 2016-4-20 22:00:09 | 只看该作者
本帖最后由 freedom 于 2016-4-22 09:47 编辑

                                                         ok了

15

帖子

6

主题

0

精华

入门

F 币
73 元
贡献
45 点
板凳
 楼主| 发表于 2016-4-20 21:55:28 | 只看该作者
ok,谢谢哈
回复

使用道具 举报

835

帖子

2

主题

0

精华

大宗师

F 币
3926 元
贡献
2334 点
沙发
发表于 2016-4-20 20:38:52 | 只看该作者
如果代码有误,给出错误提示或者截图。
subroutine simp 中,n未初始化。
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

捐赠本站|Archiver|关于我们 About Us|小黑屋|Fcode ( 京ICP备18005632-2号 )

GMT+8, 2024-12-26 17:21

Powered by Tencent X3.4

© 2013-2024 Tencent

快速回复 返回顶部 返回列表