Fortran Coder

查看: 19944|回复: 10
打印 上一主题 下一主题

[求助] 写了一个程序,结果输出的数组是一个值,求赐教啊

[复制链接]

55

帖子

17

主题

0

精华

熟手

F 币
261 元
贡献
169 点
跳转到指定楼层
楼主
发表于 2015-8-20 15:02:14 | 只看该作者 回帖奖励 |正序浏览 |阅读模式
[Fortran] 纯文本查看 复制代码
program solard
     implicit none
     REAL TMAX, TMIN, E, TTHIRTY
     INTEGER I, YEAR, DOY ,a 
     INTEGER ,parameter:: ndays=730
     REAL  B,delta,sita,R,womiga0,womiga1,womiga2,womiga3
     REAL  tuo0,alfa,C,b0,PAI,FAI,PZ,P0,Tfmax,q1
     real b1,b2,I0,D
     real, allocatable::Q0(:)
     real, allocatable::Qh(:)
     real, allocatable::TTMAX(:)
 !      NDAYS=730
       tuo0=0.87
       alfa=-0.000061
       C=1.5
       b0=0.031
       b1=0.201
       b2=0.185
       I0=1367.0
       D=86400
       PAI=3.14159
       FAI=38.33/(180*PAI)
       PZ=90810
       P0=101325
    DO I=1,NDAYS
     OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
      rewind(9)
     READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY
     enddo
    Allocate( Q0( ndays ) )
    Allocate( Qh( ndays ) )
    Allocate( TTMAX( ndays ) )
    DO I=1,NDAYS
    !TTHIRTY 30日平均温差
      delta=23.45*PAI/(180*sin(2*PAI*(284+i)/365))
      sita=2*PAI*i/365
      R=1.000109+0.033494*cos(sita)+0.001472*sin(sita)
      R=R+0.000768*cos(2*sita)+0.000079*sin(2*sita)
      womiga0=ABS(acos(-tan(delta)*tan(fai)))
      q1=womiga0*sin(fai)*sin(delta)+cos(fai)*cos(delta)*sin(womiga0)
      q0=q1*I0*D/(PAI*R*R)
      womiga3=0
        DO A=1,100
            womiga1=-womiga0+womiga0/100+womiga0*a/50
	        womiga2=sin(fai)*sin(delta)+cos(fai)*cos(delta)*cos(womiga1)
	        womiga3=womiga3+(womiga2*tuo0**(PZ/(P0*womiga2)))*womiga0/50
        ENDDO
        Ttmax=womiga3/(q0*2*PAI*R*R/I0/D)+alfa*e*1000;
       B=b0+b1*exp(-b2* TTHIRTY);
       Tfmax=1.0-0.9*exp(-B*(Tmax-Tmin)**(C));
       Qh=q0*Ttmax*Tfmax/1000000;
    ENDDO
     open(2,file='QH.dat')
       do I=1,ndays-1
	    write(2,'(i4,2x,f15.1,f6.1,f15.1)') I,Q0(I),TTMAX(I),QH(I)
       enddo 
    CLOSE(2)
    END

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

954

帖子

0

主题

0

精华

大师

F 币
184 元
贡献
75 点

规矩勋章元老勋章新人勋章水王勋章热心勋章

QQ
11#
发表于 2015-8-21 19:49:31 | 只看该作者
在我这里,并没有出现你截图中的 end-of-file 错误。而是浮点数错误。
究其原因,在于
womiga0 = abs(acos(-tan(delta)*tan(fai)))
此处的 tan(delta)*tan(fai) 无法保证在 [ -1 , 1 ] 区间内。当 delta 大于 -1.5049 时,tan(delta)*tan(fai) 的结果为 1.0288
至此,无法计算 acos
你可以在此句前面加上 write(*,*) delta,tan(delta),tan(delta)*tan(fai) 以便观察。

最后,对于这种运行时错误,建议你学会自行debug,可参考此文:http://debug.w.fcode.cn

55

帖子

17

主题

0

精华

熟手

F 币
261 元
贡献
169 点
10#
 楼主| 发表于 2015-8-21 10:47:08 | 只看该作者
数据文件

SOLARDFILE.dat

22.45 KB, 下载次数: 2

回复

使用道具 举报

2033

帖子

12

主题

5

精华

论坛跑堂

臭石头雪球

F 币
1641 元
贡献
709 点

美女勋章热心勋章星光勋章新人勋章贡献勋章管理勋章帅哥勋章爱心勋章规矩勋章元老勋章水王勋章

9#
发表于 2015-8-21 10:42:11 | 只看该作者
请给出数据文件(上传附件)

17

帖子

4

主题

0

精华

熟手

F 币
186 元
贡献
105 点
8#
发表于 2015-8-21 10:16:33 | 只看该作者
开头自己声明
integer::ierror
character(len=80)::msg
自己改
read(9,*) 改成 read(9,*,iostat=ierror,iomsg=msg)  
后面加一行
if (ierror/=0) write(*,*) msg
报什么错误自己调

评分

参与人数 1F 币 +9 贡献 +9 收起 理由
fcode + 9 + 9 很给力!

查看全部评分

55

帖子

17

主题

0

精华

熟手

F 币
261 元
贡献
169 点
7#
 楼主| 发表于 2015-8-21 10:15:40 | 只看该作者
vvt 发表于 2015-8-20 16:39
你改得不对

[mw_shl_code=fortran,true]

?可是你这样写读不出来啊

55

帖子

17

主题

0

精华

熟手

F 币
261 元
贡献
169 点
6#
 楼主| 发表于 2015-8-21 09:46:06 | 只看该作者
读不到那个位置是什么原因

微信截图_20150821095334.png (49.17 KB, 下载次数: 1)

微信截图_20150821095334.png

954

帖子

0

主题

0

精华

大师

F 币
184 元
贡献
75 点

规矩勋章元老勋章新人勋章水王勋章热心勋章

QQ
5#
发表于 2015-8-20 16:39:30 | 只看该作者
本帖最后由 vvt 于 2015-8-20 16:43 编辑

你改得不对

[Fortran] 纯文本查看 复制代码
program solard
     implicit none
     REAL TMAX, TMIN, E, TTHIRTY
     INTEGER I, YEAR, DOY ,a 
     INTEGER ,parameter:: ndays=730
     REAL  B,delta,sita,R,womiga0,womiga1,womiga2,womiga3
     REAL  tuo0,alfa,C,b0,PAI,FAI,PZ,P0,Tfmax,q1
     real b1,b2,I0,D
     real, allocatable::Q0(:)
     real, allocatable::Qh(:)
     real, allocatable::TTMAX(:)
 !      NDAYS=730
       tuo0=0.87
       alfa=-0.000061
       C=1.5
       b0=0.031
       b1=0.201
       b2=0.185
       I0=1367.0
       D=86400
       PAI=3.14159
       FAI=38.33/(180*PAI)
       PZ=90810
       P0=101325
    !!!DO I=1,NDAYS
     OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
    !!  rewind(9)
    !! READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY
    !! enddo
    Allocate( Q0( ndays ) )
    Allocate( Qh( ndays ) )
    Allocate( TTMAX( ndays ) )
    DO I=1,NDAYS
    !TTHIRTY 30日平均温差
      READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY !////// 加到这里
      delta=23.45*PAI/(180*sin(2*PAI*(284+i)/365))
      sita=2*PAI*i/365
      R=1.000109+0.033494*cos(sita)+0.001472*sin(sita)
      R=R+0.000768*cos(2*sita)+0.000079*sin(2*sita)
      womiga0=ABS(acos(-tan(delta)*tan(fai)))
      q1=womiga0*sin(fai)*sin(delta)+cos(fai)*cos(delta)*sin(womiga0)
      q0(i)=q1*I0*D/(PAI*R*R)
      womiga3=0
        DO A=1,100
            womiga1=-womiga0+womiga0/100+womiga0*a/50
         womiga2=sin(fai)*sin(delta)+cos(fai)*cos(delta)*cos(womiga1)
         womiga3=womiga3+(womiga2*tuo0**(PZ/(P0*womiga2)))*womiga0/50
        ENDDO
        Ttmax(i)=womiga3/(q0(i)*2*PAI*R*R/I0/D)+alfa*e*1000;
       B=b0+b1*exp(-b2* TTHIRTY);
       Tfmax=1.0-0.9*exp(-B*(Tmax-Tmin)**(C));
       Qh(i)=q0(i)*Ttmax(i)*Tfmax/1000000;
    ENDDO
     open(2,file='QH.dat')
       do I=1,ndays-1
     write(2,'(i4,2x,f15.1,f6.1,f15.1)') I,Q0(I),TTMAX(I),QH(I)
       enddo 
    CLOSE(2)
    END

55

帖子

17

主题

0

精华

熟手

F 币
261 元
贡献
169 点
地板
 楼主| 发表于 2015-8-20 16:29:46 | 只看该作者
vvt 发表于 2015-8-20 16:04
    DO I=1,NDAYS
     OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
...

[Fortran] 纯文本查看 复制代码
program solard
     implicit none
     REAL TMAX, TMIN, E, TTHIRTY
     INTEGER I, YEAR, DOY ,a 
     INTEGER ,parameter:: ndays=730
     REAL  B,delta,sita,R,womiga0,womiga1,womiga2,womiga3
     REAL  tuo0,alfa,C,b0,PAI,FAI,PZ,P0,Tfmax,q1
     real b1,b2,I0,D
     real, allocatable::Q0(:)
     real, allocatable::Qh(:)
     real, allocatable::TTMAX(:)
 !      NDAYS=730
       tuo0=0.87
       alfa=-0.000061
       C=1.5
       b0=0.031
       b1=0.201
       b2=0.185
       I0=1367.0
       D=86400
       PAI=3.14159
       FAI=38.33/(180*PAI)
       PZ=90810
       P0=101325
    !38.33是样地纬度
     OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
     DO I=1,NDAYS
      rewind(9)
     READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY
     enddo
     close(9)
    Allocate( Q0( ndays ) )
    Allocate( Qh( ndays ) )
    Allocate( TTMAX( ndays ) )
    DO I=1,NDAYS
    !TTHIRTY 30日平均温差
      delta=23.45*PAI/(180*sin(2*PAI*(284+i)/365))
      sita=2*PAI*i/365
      R=1.000109+0.033494*cos(sita)+0.001472*sin(sita)
      R=R+0.000768*cos(2*sita)+0.000079*sin(2*sita)
      womiga0=ABS(acos(-tan(delta)*tan(fai)))
      q1=womiga0*sin(fai)*sin(delta)+cos(fai)*cos(delta)*sin(womiga0)
      q0=q1*I0*D/(PAI*R*R)
      womiga3=0
        DO A=1,100
            womiga1=-womiga0+womiga0/100+womiga0*a/50
                womiga2=sin(fai)*sin(delta)+cos(fai)*cos(delta)*cos(womiga1)
                womiga3=womiga3+(womiga2*tuo0**(PZ/(P0*womiga2)))*womiga0/50
        ENDDO
        Ttmax=womiga3/(q0*2*PAI*R*R/I0/D)+alfa*e*1000;
       B=b0+b1*exp(-b2* TTHIRTY);
       Tfmax=1.0-0.9*exp(-B*(Tmax-Tmin)**(C));
       Qh=q0*Ttmax*Tfmax/1000000;
    ENDDO
     open(2,file='QH.dat')
       do I=1,ndays-1
            write(2,'(i4,2x,f15.1,f6.1,f15.1)') I,Q0(I),TTMAX(I),QH(I)
       enddo 
    CLOSE(2)
    END

按您说的改了,可是结果还是一样的

954

帖子

0

主题

0

精华

大师

F 币
184 元
贡献
75 点

规矩勋章元老勋章新人勋章水王勋章热心勋章

QQ
板凳
发表于 2015-8-20 16:04:20 | 只看该作者
本帖最后由 vvt 于 2015-8-20 16:07 编辑

[Fortran] 纯文本查看 复制代码
    DO I=1,NDAYS
     OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
      rewind(9)
     READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY
     enddo

这样没有意义。你读了 NDAYS 次第一行。所以每次结果都一样。

正确的做法是:
[Fortran] 纯文本查看 复制代码
OPEN (UNIT = 9, FILE = '\B90V3\DAT\SOLARDFILE.DAT')
DO I=1,NDAYS  
    !TTHIRTY 30日平均温差 
     READ (9, *) YEAR, DOY, TMAX, TMIN, E, TTHIRTY
    ! 计算
END DO
Close( 9 )
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

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

GMT+8, 2024-12-23 09:41

Powered by Tencent X3.4

© 2013-2024 Tencent

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