Fortran Coder

标题: 算出的y2结果错了,请求帮忙! [打印本页]

作者: 马媛    时间: 2015-7-6 16:20
标题: 算出的y2结果错了,请求帮忙!
算出来的y1是正确的,但是算y2时算错了。找不到原因了。请大家帮我看看!

T1NWB4G_Z]$HUAWPES21F@E.png (4.22 KB, 下载次数: 223)

T1NWB4G_Z]$HUAWPES21F@E.png

作者: fcode    时间: 2015-7-6 17:42
请给出 ess.dat 文件
作者: vvt    时间: 2015-7-6 17:55
OCR 工具,嘿嘿

[Fortran] 纯文本查看 复制代码
Program qmzy11
  Implicit None
  Integer , parameter :: N = 500
  Integer :: i
  Real(Kind=8) :: x0(N), y0(N) , x1 , y1 , x2 , y2
  Open (11, File='ess.dat')
  Open (12, File='dt.dat')
  Do i = 1, 50
    Read (11, *) x0(i), y0(i)
  End Do
  x1 = 0.2655D0
  Call lagr(x0, y0, N, x1, y1)
  x2 = 0.3886D0
  Call lagr(x0, y0, N, x2, y2)
  Write (12, 1) x1, x2
  Write (*, 1) y1, y2
  1 Format (7X, 2E15.7)
End Program qmzy11

Subroutine lagr(x0, y0, n, x, y)
  Implicit None
  Integer :: i , j , n
  Real(Kind=8) :: x0(n), y0(n) , p , x , y
  y = 0.D0
  Do i = 1, n
    p = 1.0D0
    Do j = 1, n
      If (i==j) cycle
      p = p*(x-x0(j))/(x0(i)-x0(j))
    End Do
    y = y + p*y0(i)
  End Do
End Subroutine lagr





欢迎光临 Fortran Coder (http://bbs.fcode.cn/) Powered by Discuz! X3.2