24小时热门版块排行榜

返回列表

当前主题已经存档。

zyj8119

木虫 (著名写手)

应助: 65 (初中生)
贵宾: 0.003
金币: 915.1
散金: 1440
红花: 35
帖子: 2936
在线: 1329.4小时
虫号: 664177
注册: 2008-11-29
性别: GG
专业: 理论和计算化学

[交流] 【求助】MS建立的结果MCM-41胺基铆接的FORTRAN程序问题已有1人参与

我自己建立了MCM-41,然后导出为PDB文件，自己写了个程序，编译都没有错误，就是运行不出来。
程序时这样的：    program main
integer natom,atomsn
parameter(natom=12348,nN=1000)
real xc,yc,zc,dl,coor(natom,3),zuo(nN,3)
character atomlabel(natom)*2,tlabel(natom)*2
   open(10,file='1.txt',status='old')
do 20 j=1,natom
read(10,*)j,atomlabel(j),coor(j,1),coor(j,2),coo
   & r(j,3),tlabel(j)
20 continue
   close(10)
C    generate the initiated coordinates for nitrogen atom
C    check the distance between nitrogen and atoms in framework
   do 25 i=1,nN
  do 26 j=1,3
zuo(i,j)=0
26 continue
25 continue
   do 30 i=1,nN
50 xc=RAN2(IDUM)*58.015
   yc=RAN2(IDUM)*133.80
zc=RAN2(IDUM)*39.168
do 40 j=1,natom
dl=sqrt((xc-coor(j,1))**2+(yc-coor(j,2))**2+(zc-coor(j,3))**2)
if(dl.lt.0.25)goto 50
40 continue
C    check the distance between nitrogen atoms
   if(i.ge.2)then
do 60 k=1,i-1
  dl=sqrt((xc-zuo(k,1))**2+(yc-zuo(k,2))**2+(zc-zuo(k,3))**2)
  if(dl.lt.0.25)goto 50
60 continue
   end if
   zuo(i,1)=xc
zuo(i,2)=yc
zuo(i,3)=zc
30 continue
C    ouput the file for next simulation
open(70,file='2.txt',access='append')
   NAME='N'
C    output the coordinates of nitrogen
do 120 j=1,natom
   write(70,666)j,atomlabel(j),coor(j,1),coor(j,2),
   & coor(j,3),tlabel(j)
120 continue
C    output the coordinates of atoms in framework
do 140 i=1,nN
atomsn=natom+i
write(70,666)atomsn,NAME,zuo(atomsn,1),zuo(a
   & tomsn,2),zuo(atomsn,3),NAME
140 continue
666 format(I5.1,A2,3X,3F5.3,3X,A2)
   close(70)
end

C    obtain the randam number between 0 and 1
FUNCTION RAN2(INUM)
INTEGER idum,IM1,IM2,IMM1,IA1,IA2,IQ1,IQ2,IR1,IR2,NTAB,NDIV
REAL ran2,AM,eps,RNMX
PARAMETER(IM1=214783563,IM2=214783399,AM=1./IM1,IMM1=IM1-1,
   *IA1=40014,IA2=40692,IQ1=53668,IQ2=52774,IR1=12211,IR2=3791,
   *NTAB=32,NDIV=1+IMM1/NTAB,EPS=1.2E-7,RNMX=1.-EPS)
INTEGER idum2,j,k,iv(NTAB),iy
  SAVE iv,iy,idum2
  DATA idum2/123456789/,iv/NTAB*0/,iy/0/
if(idum.le.0)then
idum=max(-idum,1)
idum2=idum
do 11 j=NTAB+8,1,-1
  k=idum/iq1
idum=IA1*(idum-k*IQ1)-k*IR1
if(idum.lt.0)idum=idum+IM1
if(j.le.NTAB)iv(j)=idum
11    continue
   iy=iv(1)
end if
k=idum/IQ1
idum=IA1*(idum-k*IQ1)-k*IR1
if(idum.lt.0)idum=idum+IM1
k=idum2/IQ2
idum2=IA2*(idum2-k*IQ2)-k*IR2
if(idum2.lt.0)idum2=idum2+IM2
j=1+iy/NDIV
iy=iv(j)-idum2
iv(j)=idum
if(iy.lt.1)iy=iy+IMM1
ran2=min(AM*iy,RNMX)
return
END
1.txt文件是这样的：
1 SI 15.65 37.172 15.303 Si
2 SI 15.762 56.005 14.655 Si
3 SI 14.211 74.922 14.369 Si
4 SI 16.512 94.068 14.554 Si
5 SI 14.894 131.378 15.386 Si
6 SI 34.373 0.61 16.75 Si
7 SI 53.776 17.43 16.264 Si
8 SI 53.474 114.575 15.138 Si
9 SI 11.633 4.363 14.668 Si
10 SI 14.412 24.617 13.956 Si
11 SI 13.28 43.362 13.062 Si
12 SI 13.565 82.058 12.901 Si
13 SI 12.393 121.022 12.86 Si
14 SI 33.626 102.581 12.622 Si
15 SI 52.292 82.29 13.707 Si
16 SI 18.889 32.386 18.925 Si
17 SI 18.202 70.681 17.455 Si
18 SI 18.403 126.612 18.51 Si
19 SI 38.236 32.885 18.612 Si
20 SI 36.352 71.09 19.168 Si
21 SI 36.798 129.172 0.007 Si
22 SI 56.176 12.407 20.212 Si
23 SI 57.832 52.169 19.696 Si
24 SI 18.002 31.644 13.86 Si
25 SI 16.688 71.132 12.441 Si
26 SI 17.71 127.019 13.594 Si
27 SI 36.757 32.941 13.735 Si
28 SI 56.708 52.054 14.185 Si
29 SI 8.206 20.724 16.217 Si
这里只截取前面几行，清高手指导下，说是数组越界，不知道为什么？

回复此楼