Fortran Coder

查看: 4435|回复: 2
打印 上一主题 下一主题

[其他行业算法] 牛顿法——求解程序

[复制链接]

23

帖子

9

主题

0

精华

熟手

F 币
124 元
贡献
79 点
跳转到指定楼层
楼主
发表于 2015-2-4 11:43:03 | 只看该作者 回帖奖励 |倒序浏览 |阅读模式
目的是想建立一个通用的牛顿求解算法的框架,就是说,对于不同的问题,只需要将初值A0,以及问题FX,导数DFX输入,调用牛顿算法,就可以求出结果。


我写的程序如下:有错误,请高手指点一下应该在怎样写?多谢多谢!!
[Fortran] 纯文本查看 复制代码
subroutine module_and_algorithm
implicit none
double precision,external :: cubic_roots,der_cubic_roots
integer,parameter :: length=4
double precision Co(length)
double precision :: x0 =1d0
double precision y
co=(/3d0,12d0,0d0,1d0/)
call newton(x0,cubic_roots,der_cubic_roots,y)
print *, "this is the roots near X0",y

end subroutine
subroutine cubic_roots(Co,length,y)
implicit none
integer length,i
double precision x
double precision y
double precision temp(length)
do i=1,length
temp(i)=co(i)*x**(i-1)
end do
y=sum(temp)
return
end subroutine
subroutine der_cubic_roots(Co,length,y)
implicit none
integer length,i
double precision x
double precision y
double precision temp(length-1)
do i=2,length
temp(i)=dble(i-1)*co(i)*x**(length-2)
end do
y=sum(temp)
return
end subroutine
subroutine newton(a,f,df,b)
implicit none
double precision :: epsilon=1d-6
double precision a,b,fb
double precision,external :: f,df

b=a-f(a)/df(a)
fb=f(b)

print *, "i am here (1)",fb
do while(abs(fb)>epsilon)
a=b
b=a-f(a)/df(a)
fb=f(b)
print *, "i am here (2)", fb
end do
return
end subroutine
program main
call module_and_algorithm
end

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

23

帖子

9

主题

0

精华

熟手

F 币
124 元
贡献
79 点
沙发
 楼主| 发表于 2015-2-4 14:08:44 | 只看该作者
解决了  麻烦各位了 。。

796

帖子

2

主题

0

精华

大宗师

F 币
3787 元
贡献
2266 点
板凳
发表于 2015-2-4 15:52:50 | 只看该作者
本帖最后由 li913 于 2015-2-4 15:54 编辑

初值A0,函数FX,导数DFX,这三者都作为输入参数。
subroutine Newton(a0,fx,dfx,...,...)这样的好处就是代码可不经修改(或少量修改),用于其他项目。
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

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

GMT+8, 2024-4-20 22:45

Powered by Tencent X3.4

© 2013-2024 Tencent

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