24小时热门版块排行榜    

查看: 905  |  回复: 21
当前主题已经存档。
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

tpp001

金虫 (著名写手)

[交流] 【求助】新手请教一段程序的疑问....【已完结】

我的疑问是,在do i=1,natom, read (1,*) molname,atomname,index,x(i),y(i),z(i),这里do 做 一个循环,我后面的输出要求是输出,x(i-1),y(i-1),z(i-1), x(i),y(i),z(i),x(i+1),y(i+1),z(i+1)这3个坐标,这样输出可以吗??

open(1,file='benz-22.gro',status='unknown') !打开benz-co2.gro文件,并读取
        do istep=1, nstep
                read(1,*) sysname
                read(1,*) natom
                                do i=1, natom
                                        read(1,*) molname,atomname,index,x(i),y(i),z(i)
                                                if (atomname .eq.'Cb') then
                                                        r2=(x(i)-2.623)**2+(y(i)-2.554)**2
                                                if (r2 .LT. 0.01932) then
                                                        h2=(z(i)-2.608)**2
                                                if (h2 .LT. 0.5776 .AND. z(i) .NE. 2.608) then
                                                        open(20,file='aromaticco2.txt')
                                                        write(20,*)'Oa',x(i-1),y(i-1),z(i-1) !输出这个Cb原子的坐标和Cb编号的前面和后面的一个原子的坐标。
                                                        write(20,*)'Cb',x(i),y(i),z(i)
                                                        write(20,*)'Oc',x(i+1),y(i+1),z(i+1)

[ Last edited by 波不动 on 2009-11-25 at 23:52 ]
回复此楼
迷茫在知识的海洋里,需要你的指导。thankyou
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tpp001

金虫 (著名写手)

引用回帖:
Originally posted by nono2009 at 2009-11-25 22:25:
输入文件中的数据少了、或者read的过多了。

我贴一个小一点的输入文件的格式
Generated by trjconv : Pure Mix - Yummie! (and some water) t= 1500.00000 !!这是体系的名字sysname
   14 !!这是体系的原子个数
    1BEN     C1    1   3.265   3.036   3.125
    1BEN     C2    2   3.170   2.942   3.118
    1BEN     C3    3   3.038   2.983   3.131
    1BEN     C4    4   3.006   3.118   3.139
    1BEN     C5    5   3.106   3.216   3.143
    1BEN     C6    6   3.238   3.173   3.129
    1BEN     H7    7   3.368   3.000   3.126
    1BEN     H8    8   3.197   2.836   3.111
    1BEN     H9    9   2.957   2.912   3.141
    1BEN    H10   10   2.902   3.146   3.143
    1BEN    H11   11   3.317   3.247   3.126
    1BEN    C12   12   3.083   3.361   3.139
    1BEN    O13   13   2.979   3.418   3.111
    1BEN    H14   14   3.171   3.424   3.162
   6.25808   6.25808   6.25808 !!这是盒子的大小
Generated by trjconv : Pure Mix - Yummie! (and some water) t= 1500.19995
   14
    1BEN     C1    1   3.267   3.036   3.135
    1BEN     C2    2   3.173   2.940   3.131
    1BEN     C3    3   3.041   2.983   3.130
    1BEN     C4    4   3.004   3.117   3.124
    1BEN     C5    5   3.105   3.213   3.125
    1BEN     C6    6   3.237   3.172   3.135
    1BEN     H7    7   3.371   3.009   3.124
    1BEN     H8    8   3.196   2.836   3.109
    1BEN     H9    9   2.957   2.914   3.123
    1BEN    H10   10   2.904   3.152   3.101
    1BEN    H11   11   3.324   3.237   3.147
    1BEN    C12   12   3.074   3.358   3.109
    1BEN    O13   13   2.983   3.427   3.147
    1BEN    H14   14   3.161   3.406   3.059
   6.25808   6.25808   6.25808
迷茫在知识的海洋里,需要你的指导。thankyou
17楼2009-11-25 22:40:21
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 22 个回答

jjdg

版主 (知名作家)

优秀版主


小木虫(金币+0.5):给个红包,谢谢回帖交流
努力学习!以正当途径!获得需要的知识!
2楼2009-11-24 12:43:59
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

maomao1210

金虫 (正式写手)


小木虫(金币+0.5):给个红包,谢谢回帖交流
引用回帖:
Originally posted by tpp001 at 2009-11-24 09:25:
我的疑问是,在do i=1,natom, read (1,*) molname,atomname,index,x(i),y(i),z(i),这里do 做 一个循环,我后面的输出要求是输出,x(i-1),y(i-1),z(i-1), x(i),y(i),z(i),x(i+1),y(i+1),z(i+1)这3个坐标,这样输出可以 ...

很显然不行,i=1的时候,你仅仅读入了对应第一个原子的信息,何来i=2?。
3楼2009-11-24 16:49:36
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

maomao1210

金虫 (正式写手)


小木虫(金币+0.5):给个红包,谢谢回帖交流
open(1,file='benz-22.gro',status='unknown') !打开benz-co2.gro文件,并读取
        do istep=1, nstep
                read(1,*) sysname
                read(1,*) natom
          I_emuch=1
           do i=1, natom
          read(1,*) molname,atomname,index,x(i),y(i),z(i)
           if (atomname .eq.'Cb') then
             r2=(x(i)-2.623)**2+(y(i)-2.554)**2
            if (r2 .LT. 0.01932) then
             h2=(z(i)-2.608)**2
        if (h2 .LT. 0.5776 .AND. z(i) .NE. 2.608) then

             Index_emuch(I_emuch)=i
               endif
               endif
               endif
            I_emuch=I_emuch+1
          Enddo

         Do J_emuch=1,I_emuch
        open(20,file='aromaticco2.txt')
       if(ii.gt.1)write(20,*)'Oa',x(ii-1),y(ii-1),z(ii-1) !输出这个Cb原子的坐标和Cb编号的前面和后面的一个原子的坐标。
        write(20,*)'Cb',x(ii),y(ii),z(ii)
        write(20,*)'Oc',x(ii+1),y(ii+1),z(ii+1)        
           Enddo
          Enddo
4楼2009-11-24 17:04:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
普通表情 高级回复 (可上传附件)
信息提示
请填处理意见