爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 6695|回复: 3

[求助] 雅克比行列式的子程序有什么问题呢,为什么数组溢出了

[复制链接]

新浪微博达人勋

发表于 2012-12-17 18:42:53 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
C==================== *SUB JACOBI*===========================
      SUBROUTINE JACOBI(N,A,D,V,EVL,B,Z)
      DIMENSION A(N,N),D(N),V(N,N),B(N),Z(N)
      LOGICAL EVL
      IF(.NOT.EVL) GO TO 10
      DO 70 IP=1,N
      DO 70 IQ=1,N
      IF (IP-IQ) 50,60,50
  60  V(IP,IP)=1.0
      GO TO 70
  50  V(IP,IQ)=0.0
  70  CONTINUE
  10  DO 80 IP=1,N
      D(IP)=A(IP,IP)
      B(IP)=D(IP)
  80  Z(IP)=0.0
      IROT=0
      DO 370 I=1,50
      SM=0.0
      NM1=N-1
      DO 100 IP=1,NM1
      IPP1=IP-1
      DO 100 IQ=IPP1,N
100  SM=SM+ABS(A(IP,IQ))
      IF (SM) 110,120,110
110  IF (I-4) 130,140,140
130  TRESH=0.2*SM/(FLOAT(N)*FLOAT(N))
      GO TO 150
140  TRESH=0.0
150  DO 160 IP=1,NM1
      IPP1=IP+1
      DO 160 IQ=IPP1,N
      G=100.0*ABS(A(IP,IQ))
      IF(I.GT.4.AND.ABS(D(IP))+G.EQ.ABS(D(IP)).AND.ABS(D(IQ))+G.EQ.
     *ABS(D(IQ))) GO TO 200
      IF(ABS(A(IP,IQ)).LE.TRESH) GO TO 160
      H=D(IQ)-D(IP)
      IF(ABS(H)+G.EQ.ABS(H)) GO TO 240
      THETA=0.5*H/A(IP,IQ)
      T=1.0/(ABS(THETA)+SQRT(1.0+THETA*THETA))
      IF(THETA.LT.0.0) T=-T
      GO TO 250
240  T=A(IP,IQ)/H
250  C=1.0/SQRT(1.0+T*T)
      S=T*C
      H=T*A(IP,IQ)
      Z(IP)=Z(IP)-H
      Z(IQ)=Z(IQ)+H
      D(IP)=D(IP)-H
      D(IQ)=D(IQ)+H
      A(IP,IQ)=0.0
      IPM1=IP-1
      IF (IPM1) 260,260,270
270  DO 280 J=1,IPM1
      G=A(J,IP)
      H=A(J,IQ)
      A(J,IP)=C*G-S*H
280  A(J,IQ)=S*G+C*H
260  IQM1=IQ-1
      IF (IQM1-IPP1) 300,290,290
290  DO 310 J=IPP1,IQM1
      G=A(IP,J)
      H=A(J,IQ)
      A(IP,J)=C*G-S*H
310  A(J,IQ)=S*G+C*H
300  IQP1=IQ+1
      IF (N-IQP1) 330,320,320
320  DO 340 J=IQP1,N
      G=A(IP,J)
      H=A(IQ,J)
      A(IP,J)=C*G-S*H
340  A(IQ,J)=S*G+C*H
330  IF  (.NOT.EVL) GO TO 350
      DO 360 J=1,N
      G=V(J,IP)
      H=V(J,IQ)
      V(J,IP)=C*G-S*H
360  V(J,IQ)=S*G+C*H
350  IROT=IROT+1
      GO TO 160
200  A(IP,IQ)=0.0
160  CONTINUE
      DO 370 IP=1,N
      B(IP)=B(IP)+Z(IP)
      D(IP)=B(IP)
370  Z(IP)=0.0
120   RETURN
       END
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 成长值: 0
发表于 2012-12-18 08:05:50 | 显示全部楼层
这玩意儿没法看,又是标号又是goto的,画地图一样···
建议debug吧,fortran里面的build选项下面有debug,点一下,程序会停在运行不下去的地方

ps:JCB主要中用来求解矩阵的,一般如果你没有修改子程序是不会有问题的,另外对于写的不好的程序,可能会先溢出之类的,可以再找其他的jcb子程序用着看看···
密码修改失败请联系微信:mofangbao

新浪微博达人勋

0
早起挑战累计收入
发表于 2012-12-18 08:46:33 | 显示全部楼层
可能是你调用时传进去的参数有问题
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2021-8-23 18:24:51 | 显示全部楼层
请问楼主最后问题解决了吗
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

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