论坛 › K08：编程语言 › 请帮我看一下这段程序,谢谢!

zxjcool 发表于 2007-4-24 22:22:20

请帮我看一下这段程序,谢谢!

模型是这样的: 在一个边长为35的正方体晶格中心有一个大离子,其表面吸附有10个小离子,大离子半径为5,小离子半径为0.5,晶格其他空间内随机分布1000个溶剂粒子,半径为0.5. 它们间的约束条件为小离子-小离子之间的距离应大于1.122462,大离子-溶剂粒子之间距离大于5.622462,溶剂粒子-溶剂粒子之间距离大于1.122462.以下是我编的程序(其中一部分). 调试后出现的问题是它们之间的约束不正确. 谁帮我看一下,好吗? 谢谢!
g=1.122462
nparticle=1000!1000个溶剂粒子
r(1)=(4.5+g)*(4.5+g)!5.622462
r(2)=g*g
r(3)=g*g

xbox=35!晶格尺寸
ybox=35
zbox=35

x(1)=0
y(1)=0
z(1)=0
count=1

C   小离子在球面上的坐标
call seed(RND$TIMESEED)
do 200, i=1,10
count=count+1
Call random(ran)
Call random(rand)
x(count)=5.5*sin(3.1415926*ran)*cos(2*3.1415926*rand)
y(count)=5.5*sin(3.1415926*ran)*sin(2*3.1415926*rand)
z(count)=5.5*cos(3.1415926*ran)
200   continue

C溶剂粒子在晶格中随机分布
count=11
do 100, i=1, nparticle
                        count=count+1
                     Call random(ran)
   x(count)=35*ran
                      Call random(ran)
               y(count)=35*ran
                  Call random(ran)
                      z(count)=35*ran
100 continue            
      n=count

C开始约束距离
do s=3,11
140 do k=2,s-1
d4=(x(s)-x(k))**2+(y(s)-y(k))**2+(z(s)-z(k))**2
if (d4.le.r(2)) then
Call random(ran)
Call random(rand)
x(s)=5.5*sin(3.1415926*ran)*cos(2*3.1415926*rand)
y(s)=5.5*sin(3.1415926*ran)*sin(2*3.1415926*rand)
z(s)=5.5*cos(3.1415926*ran)
goto 140
else
endif
enddo
enddo

do i=12,n
170 do k=2,i-1
      d2=(x(i)-x(k))**2+(y(i)-y(k))**2+(z(i)-z(k))**2
d1=(x(i)-x(1))**2+(y(i)-y(1))**2+(z(i)-z(1))**2
if ((d2.le.r(2)).AND.(d1.le.r(1)))then
               Call random(ran)
   x(i)=35*ran
                      Call random(ran)
               y(i)=35*ran
               Call random(ran)
                      z(i)=35*ran
goto 170
else
endif
enddo
enddo

C   将系统坐标原点置于原子所在区域的中心
      
do 190, i=12,n
x(i)=x(i)-17.5
y(i)=y(i)-17.5
z(i)=z(i)-17.5
190 continue

C输出距离,发现距离还是没有约束好
   do i=1,n
   rxi=x(i)
   ryi=y(i)
   rzi=z(i)
   do j=i+1, n
               xr=rxi-x(j)
   yr=ryi-y(j)
   zr=rzi-z(j)
               rij_2=xr*xr+yr*yr+zr*zr
               rij = sqrt(rij_2)
if(rij_2.le.r(2)) then      !(这边我输出距离时发现没有约束)
write(*,*)rij_2,i,j
pause 'n'
endif
enddo
enddo

oasis_luo 发表于 2007-4-25 08:09:24

goto 140 和 goto 170的逻辑上不对。应该在goto 前加判断语句。因为你原来程序
的意思是如果不满足约束则重新产生新的位置。而且默认这样一次之后就满足
了约束（往往不可能）。所以，必须反复产生新位置和反复计算距离直到确认
满足距离约束位置再算下一个离子。

zxjcool 发表于 2007-4-26 11:10:40

我改了一下,怎么还不行? 是不是哪里又错了?帮我看下好吗?
do s=3,11
        do k=2,s-1
        d4=(x(s)-x(k))**2+(y(s)-y(k))**2+(z(s)-z(k))**2
        if (d4.le.r(2)) then
140   Call random(ran)
        Call random(rand)
      x(s)=5.5*sin(3.1415926*ran)*cos(2*3.1415926*rand)
        y(s)=5.5*sin(3.1415926*ran)*sin(2*3.1415926*rand)
        z(s)=5.5*cos(3.1415926*ran)
      do l=2,s-1
        d4=(x(s)-x(l))**2+(y(s)-y(l))**2+(z(s)-z(l))**2
        if (d4.le.r(2)) then
        goto 140
        endif
        enddo
        endif
        enddo
        enddo

        do i=12,n
        do k=2,i-1
      d2=(x(i)-x(k))**2+(y(i)-y(k))**2+(z(i)-z(k))**2
        d1=(x(i)-x(1))**2+(y(i)-y(1))**2+(z(i)-z(1))**2
        if ((d2.le.r(2)).AND.(d1.le.r(1)))then
170                      Call random(ran)
                                        x(i)=35*ran
                      Call random(ran)
                      y(i)=35*ran
                      Call random(ran)
                      z(i)=35*ran
      do l=2,i-1
      d2=(x(i)-x(l))**2+(y(i)-y(l))**2+(z(i)-z(l))**2
        d1=(x(i)-x(1))**2+(y(i)-y(1))**2+(z(i)-z(1))**2
      if ((d2.le.r(2)).AND.(d1.le.r(1)))then
        goto 170
        endif
        enddo
        endif
        enddo
        enddo

oasis_luo 发表于 2007-4-26 22:09:58

你最好在计算完小离子分布后输出一个数据文件，这样可以检查小离子分布部分是
否正确。就你的程序看好象问题不大，不过还需要你自己检查检查。

置于溶剂离子分布，两个if ((d2.le.r(2)).AND.(d1.le.r(1)))语句应该改成if ((d2.le.r(2)).OR.(d1.le.r(1))). 调试程序的时候一定要检查条件判断语句，这些语句是最容易出错的地方，而且
错误比较隐蔽不容易发现。至于数学计算部分的错误，很容易就能看出来。

耐心调试，也算是增长经验，人人都是这么过来的。

zxjcool 发表于 2007-4-28 15:02:25

程序已经调通了,谢谢.果然是那个错了.

查看完整版本: 请帮我看一下这段程序,谢谢!