夏沫有约 发表于 2019-12-5 20:43:16

fortran代码运行出错请教

program main
         dimension in(2),x(2)
         data ntrap /20/
         integer i
       real*8 in,c,d
         in(1)=-sqrt(2.0)/2.0
         in(2)=sqrt(2.0)/2.0
         c=0.0
         d=-1.0
      do i=1,ntrap,1
         call external f()
                   c=x(1)
                   d=x(2)
         in(1)=in(1)-2.0*c/1.0*vdotn(in(1),x(1),1) !新方向的x坐标表达式
         in(2)=in(2)-2.0*d/1.0*vdotn(in(1),x(1),1) !新方向的y坐标表达式
                   write(*,*) 'x(1)=',x(1),'x(2)=',x(2)
                   write(*,*) 'in(1)=',in(1),'in(2)=',in(2)
            end do
          stop
          end
         
          subroutine external f()
          dimension x(2),y(2)
          data x/2*0.0/
          b=2.0
          n=2
          m=10
          eps=1.0e-5
          call dnmtc(x,n,b,m,eps,y)
          write(*,*)
          do 10 i=1,n
10          write(*,100) i,x(i)
          write(*,*)
100          format(5x,'x(',i2,1x,')=',e13.6)
          end subroutine
         
          function f(x,n)
          dimension x(n),in(n)
          f1=x(1)*x(1)+x(2)*x(2)-1.0
          f2=in(2)*x(1)-in(2)*c-in(1)*x(2)+in(1)*d
          f=sqrt(f1*f1+f2*f2)
          end

          subroutine dnmtc(x,n,b,m,eps,y)
          dimension x(n),y(n)
          double precision r
          real nrnd1
          a=b
          k=1
          r=1.0d0
          z=f(x,n)
10          if (a.ge.eps) then
            l=l+1
            do 20 i=1,n
20            y(i)=-a+2.0*a*nrnd1(r)+x(i)
            z1=f(y,n)
         k=k+1
            if (z1.ge.z) then
            if (k.gt.m) then
                k=1
                a=a/2.0
            end if
            goto 10
            else
            k=1
            do 30 i=1,n
30            x(i)=y(i)
            z=z1
            if (z.ge.eps) goto 10
            end if
          end if
          return
          end subroutine

          real function nrnd1(r)
          double precision s,u,v,r
          s=65536.0
          u=2053.0
          v=13849.0
          m=r/s
          r=r-m*s
          r=u*r+v
          m=r/s
          r=r-m*s
          nrnd1=r/s
          return
          end
我写了一个代码,主要是由一个初始点和初始方向确定的直线和一个圆相交,得到的下一个点继续作为初始点,又有一个新的方向作为新方向,继续和圆相交,一直循环。循环的次数为ntrap.我运行不知道为什么没有执行循环,输出in(1),in(2)也没有执行。
x(1)为x,x(2)为y,c为x初坐标,d为y初坐标,in(1)为直线斜率的x坐标,in(2)为直线斜率的y坐标。直线方程为:(x(1)-c)/in(1)=(x(2)-c)/in(2)
新人,还不太懂,代码可能写的有问题,没有报错,但结果不对,还请大神们指点!感谢感谢!

li913 发表于 2019-12-6 10:28:00

代码语法错误,根本就不能执行。program main
         dimension in(2),x(2)
         data ntrap /20/
         integer i
       real*8 in,c,d
         in(1)=-sqrt(2.0)/2.0
         in(2)=sqrt(2.0)/2.0
         c=0.0
         d=-1.0
      do i=1,ntrap,1
         !call external f()
         callf1()
                   c=x(1)
                   d=x(2)
         in(1)=in(1)-2.0*c/1.0*vdotn(in(1),x(1),1) !新方向的x坐标表达式
         in(2)=in(2)-2.0*d/1.0*vdotn(in(1),x(1),1) !新方向的y坐标表达式
                   write(*,*) 'x(1)=',x(1),'x(2)=',x(2)
                   write(*,*) 'in(1)=',in(1),'in(2)=',in(2)
            end do
          stop
          end
         
          subroutinef1()
          dimension x(2),y(2)
          data x/2*0.0/
          b=2.0
          n=2
          m=10
          eps=1.0e-5
          call dnmtc(x,n,b,m,eps,y)
          write(*,*)
          do 10 i=1,n
10          write(*,100) i,x(i)
          write(*,*)
100          format(5x,'x(',i2,1x,')=',e13.6)
          end subroutine
         
          function f(x,n)
          dimension x(n),in(n)
          f1=x(1)*x(1)+x(2)*x(2)-1.0
          f2=in(2)*x(1)-in(2)*c-in(1)*x(2)+in(1)*d
          f=sqrt(f1*f1+f2*f2)
          end

          subroutine dnmtc(x,n,b,m,eps,y)
          dimension x(n),y(n)
          double precision r
          real nrnd1
          a=b
          k=1
          r=1.0d0
          z=f(x,n)
10          if (a.ge.eps) then
            l=l+1
            do 20 i=1,n
20            y(i)=-a+2.0*a*nrnd1(r)+x(i)
            z1=f(y,n)
         k=k+1
            if (z1.ge.z) then
            if (k.gt.m) then
                k=1
                a=a/2.0
            end if
            goto 10
            else
            k=1
            do 30 i=1,n
30            x(i)=y(i)
            z=z1
            if (z.ge.eps) goto 10
            end if
          end if
          return
          end subroutine

          real function nrnd1(r)
          double precision s,u,v,r
          s=65536.0
          u=2053.0
          v=13849.0
          m=r/s
          r=r-m*s
          r=u*r+v
          m=r/s
          r=r-m*s
          nrnd1=r/s
          return
          end

夏沫有约 发表于 2019-12-6 19:36:04

我是在Windows下的Linux系统操作的,可以运行出来,得到的只有x(1),x(2)的值,感觉循环那块没有起作用,不知道要怎么改。
页: [1]
查看完整版本: fortran代码运行出错请教