楼主: 
YmemY (**米)   
2012-12-18 02:55:26请教大家,我在做一个一维扩散的模拟,
在我一步一步debug之后,知道问题在循环的部分,好像会无法结束??
可以成功编译成执行档,但一旦开始执行,就无法结束了,也无法输出结果..
以下是程式码,麻烦版有帮我看看哪边出错了??
      program diffusion   !1-D diffusion
      implicit none
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!宣告区
      real L,D,J,C,T,F,smallt
      !L=length,D=diffusivity,J=diffusion flux,C=initial concentration
      !F=Fourier's number, smallt=time interval
      real, parameter :: defaultL = 10
      real, parameter :: defaultD = 1.03
      real, parameter :: defaultJ = 1E-5
      real, parameter ::  defaultC = 1E-4
      real array1(-1:10001),array2(-1:10001)
      !the length is divided into 10000 pieces, 10001 nodes, 2 virtual nodes
      integer*8 p, pmax, x
      !p=pace, x=location
      character mode
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!输入设定参数
      write(*,*) "press ""1"" to select the default values, or""2"" to input
them"
      read(*,*) mode
      if(mode=="1") then
        L = defaultL
        D = defaultD
        J = defaultJ
        C = defaultC
      else
        write(*,*) "please input the length:"
        read(*,*)L
        write(*,*) "please input the diffusivity:"
        read(*,*)D
        write(*,*) "please input the diffusion flux:"
        read(*,*)J
        write(*,*) "please input the initial concentration:"
        read(*,*)C
      end if
      write(*,*) "please input the time:"
      read(*,*)T
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !计算用到的参数
      smallt = 0.25*L*L/D/10**8
      pmax = T/smallt
      F = D*smallt/(L*L)*10**8
      array1 = C
      array2 = C
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!这里就是出错的开始!!!
      do p=1,pmax
        do x=0,10000
                array1(-1) = J*L/10000/D + array1(0)
                array1(10001) = array1(10000) + array1(0) - array1(-1)
                array2(x)=F*(array1(x+1)+array1(x-1))+(1-2*F)*array1(x)
        end do
        array1 = array2
      end do
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!将阵列输出到文字档
      open(unit=10,file="diffusion1.txt")
      do x=-1,10001
        write(10,*) array1(x)
      end do
      STOP
      END
smallt = 0.25*L*L/D/10**8 检查你的smallt的值如果没有意外应该是0,导致pmax极大,循环跑不完改成smallt = 0.25*L*L/D/10.0**8.0 试试看不仅是这个地方,例如array1(-1) = J*L/10000/D这都会'让你的结果与预期的不同
作者: mk650   2012-02-18 13:56:00
t=10,pamx=41200000,x=412000000000,粉大的循环,如果参数没错,也只能丢著给它跑啦
作者: mk650   2012-02-18 14:00:00
呵三楼的精简化本人常干,不过它这个程式很简单,只有两行多余差别不大啦.重点在t要多大,t=0.0001够吗?t=1,十分钟跑完,I7 2600套装很便宜电脑,更好的电脑得会更快这样很快了,我的程式都要跑一周以上,有三周还没跑完的