24小时热门版块排行榜    

查看: 712  |  回复: 3
当前主题已经存档。

阳光不锈7717

[交流] 【求助】一个计算库伦相互作用能的fortran程序的并行化(2)

c-------------------------------------------------
c       end if

c       if(myid==2)then

c*** calculate the energy of bond bond interaction
c----------------------------------------------      



      do in=natom+1,natombond         
         xic = x(in)
         yic = y(in)
         zic = z(in)
         xi = x(in) - xic
         yi = y(in) - yic
         zi = z(in) - zic
                   c2scale= 0.0d0  
          c4scale= 0.57d0
            c3scale=0.0d0            
         do j = natom+1, na
            cscale(j) = 0.57d0
         end do
         do j = 1, nbond12(in)
            cscale(bond12(j,in)) = c2scale
         end do       
          do j = 1, nbdatom12(in)
             cscale(bdatom12(j,in)) = c2scale
        if(atomlp(1,bdatom12(j,in)).ne.0)then
         cscale(atomlp(1,bdatom12(j,in)))
     &=c2scale
        endif
      if(atomlp(2,bdatom12(j,in)).ne.0) then
         cscale(atomlp(2,bdatom12(j,in)))
     &=c2scale
       endif
        if(atompie(1,bdatom12(j,in)).ne.0)
     &  cscale(atompie(1,bdatom12(j,in)))=c2scale     
      if(atompie(2,bdatom12(j,in)).ne.0)
     &  cscale(atompie(2,bdatom12(j,in)))=c2scale
          end do
          do j = 1, nbdatom13(in)       
             cscale(bdatom13(j,in)) = c3scale
        if(atomlp(1,bdatom13(j,in)).ne.0) then
         cscale(atomlp(1,bdatom13(j,in)))
     &=c3scale
         endif
      if(atomlp(2,bdatom13(j,in)).ne.0)then
         cscale(atomlp(2,bdatom13(j,in)))
     &=c3scale
       endif
        if(atompie(1,bdatom13(j,in)).ne.0)
     &  cscale(atompie(1,bdatom13(j,in)))=c3scale     
      if(atompie(2,bdatom13(j,in)).ne.0)
     &  cscale(atompie(2,bdatom13(j,in)))=c3scale
          end do                       
         do j = 1, nbond13(in)
            cscale(bond13(j,in)) = c3scale
         end do
            fi = f * pchg(in)  

           do kk = in+1, naTOMBOND+NLP
               xc = xic - x(kk)
               yc = yic - y(kk)
               zc = zic - z(kk)
               if (use_image)  call image (xc,yc,zc,0)
               rc2 = xc*xc + yc*yc + zc*zc               
                  xr = xc + xi - x(kk) + x(kk)
                  yr = yc + yi - y(kk) + y(kk)
                  zr = zc + zi - z(kk) + z(kk)
                  r2 = xr*xr + yr*yr + zr*zr
                  r = sqrt(r2)
                  fik = fi * pchg(kk) * cscale(kk)
c        WRITE(*,*) AM(IN),AM(KK),CSCALE(KK),R

        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                  fik =  fik*lamp2
        endif
                  e = fik / r
          if(kk.le.natombond)then
         if (molcule(bdatom12(1,in)) .ne. molcule(bdatom12(1,kk))) then
                 einter = einter + e
         end if
          else
         if (molcule(bdatom12(1,in)) .ne. molcule(lpatom(kk))) then
                 einter = einter + e
         end if
          endif       
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                 eintero = eintero + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).eq.obj))then
                 eintrao = eintrao + e
         end if
           if(abs(e).gt.800)then
           write(*,*) ii,kk,molcule(II),molcule(kk)
           endif                 
c         write(11,*) 'bond-bond',e                                                                        
c-------------------------------------
          end do






           do kk = NATOMBOND+NLP+1, na
               xc = xic - x(kk)
               yc = yic - y(kk)
               zc = zic - z(kk)
               if (use_image)  call image (xc,yc,zc,0)
               rc2 = xc*xc + yc*yc + zc*zc               
                  xr = xc + xi - x(kk) + x(kk)
                  yr = yc + yi - y(kk) + y(kk)
                  zr = zc + zi - z(kk) + z(kk)
                  r2 = xr*xr + yr*yr + zr*zr
                  r = sqrt(r2)
                  fik = fi * pchg(kk) * cscale(kk)
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                  fik =  fik*lamp2
        endif
                  e = fik / r

         if (molcule(bdatom12(1,in)) .ne. molcule(PIEatom(kk))) then
                 einter = einter + e   
          endif       
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                 eintero = eintero + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).eq.obj))then
                 eintrao = eintrao + e
         end if
           if(abs(e).gt.800)then
           write(*,*) ii,kk,molcule(II),molcule(kk)
           endif                          
c            write(11,*) 'b-lp,pi',e                                                
                  ec1 = ec1 + e                    
         end do

      end do
c-----------------------------------------------------
      end if
      if(myid==3)then
c------------------------------------------------------

c         write(12,*) 'bond',ec1  
c*** calculate the energy of lp lpPIE interaction

      do in=natombond+1,natombond+nlp         
         xic = x(in)
         yic = y(in)
         zic = z(in)
         xi = x(in) - xic
         yi = y(in) - yic
         zi = z(in) - zic
                   c2scale= 0.0d0  
          c4scale= 0.57d0
            c3scale=0.0d0            
         do j = natomBOND+1, na
            cscale(j) = 0.57d0
         end do
         pi=lpatom(in)
          do j=1,2       
           if(atomlp(j,pi).ne.in.AND.NATOMLP(PI).NE.1)
     &cscale(atomlp(j,pi))=0.0d0
        enddo
           DO J=1,NATOMPIE(PI)
           CSCALE(ATOMPIE(J,PI))=0.0D0
         ENDDO

         do j=1,n12(pi)
        DO K=1,NATOMLP(I12(J,PI))
        cscale(atomlp(K,i12(j,pi)))=c2scale
        ENDDO
        DO K=1,NATOMPIE(I12(J,PI))
        cscale(atompie(K,i12(j,pi)))=c2scale
        ENDDO
           enddo

         do j=1,n13(pi)
        DO K=1,NATOMLP(I13(J,PI))
        cscale(atomlp(K,i13(j,pi)))=c3scale
      ENDDO
        DO K=1,NATOMPIE(I13(J,PI))
        cscale(atompie(K,i13(j,pi)))=c3scale
        ENDDO
           enddo

            fi = f * pchg(in)  
           do kk = in+1, naTOMBOND+NLP
               xc = xic - x(kk)
               yc = yic - y(kk)
               zc = zic - z(kk)
               if (use_image)  call image (xc,yc,zc,0)
               rc2 = xc*xc + yc*yc + zc*zc               
                  xr = xc + xi - x(kk) + x(kk)
                  yr = yc + yi - y(kk) + y(kk)
                  zr = zc + zi - z(kk) + z(kk)
                  r2 = xr*xr + yr*yr + zr*zr
                  r = sqrt(r2)
                  fik = fi * pchg(kk) * cscale(kk)
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                  fik =  fik*lamp2
        endif
                  e = fik / r
         if (molcule(lpatom(in)) .ne. molcule(lpatom(kk))) then
                 einter = einter + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                 eintero = eintero + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).eq.obj))then
                 eintrao = eintrao + e
         end if
           if(abs(e).gt.800)then
           write(*,*) ii,kk,molcule(II),molcule(kk)
           endif                                  
c         write(11,*) 'lp-lp',e                                                                                                                  
                  ec1 = ec1 + e                    
         end do

           do kk = naTOMBOND+NLP+1,NA
               xc = xic - x(kk)
               yc = yic - y(kk)
               zc = zic - z(kk)
               if (use_image)  call image (xc,yc,zc,0)
               rc2 = xc*xc + yc*yc + zc*zc               
                  xr = xc + xi - x(kk) + x(kk)
                  yr = yc + yi - y(kk) + y(kk)
                  zr = zc + zi - z(kk) + z(kk)
                  r2 = xr*xr + yr*yr + zr*zr
                  r = sqrt(r2)
                  fik = fi * pchg(kk) * cscale(kk)

        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                  fik =  fik*lamp2
        endif
                  e = fik / r
         if (molcule(lpatom(in)) .ne. molcule(PIEatom(kk))) then
                 einter = einter + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                 eintero = eintero + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).eq.obj))then
                 eintrao = eintrao + e
         end if
           if(abs(e).gt.800)then
           write(*,*) ii,kk,molcule(II),molcule(kk)
           endif       
c          write(11,*) 'lp-pie',e                                                                                                                                     
                  ec1 = ec1 + e                    
         end do

      end do
c--------------------------------------------------------
c         end if
c         if(myid==4)then
c--------------------------------------------------------

c         write(12,*) 'lp',ec1  
      do in=natombond+nlp+1,na-1         
         xic = x(in)
         yic = y(in)
         zic = z(in)
         xi = x(in) - xic
         yi = y(in) - yic
         zi = z(in) - zic
                   c2scale= 0.0d0  
          c4scale= 0.57d0
            c3scale=0.0d0            
         do j = natom+1, na
            cscale(j) = 0.57d0
         end do
         pi=pieatom(in)
        do j=1,2
           if(atompie(j,pi).ne.iN)cscale(atompie(j,pi))=0.0d0
        enddo
         do j=1,n12(pi)
        DO K=1,NATOMPIE(I12(J,PI))
        cscale(atompie(K,i12(j,pi)))=c2scale
      ENDDO
           enddo
         do j=1,n13(pi)
        DO K=1,NATOMPIE(I13(J,PI))                  
        cscale(atompie(K,i13(j,pi)))=c3scale
      ENDDO
           enddo
            fi = f * pchg(in)  
           do kk = in+1, na
               xc = xic - x(kk)
               yc = yic - y(kk)
               zc = zic - z(kk)
               if (use_image)  call image (xc,yc,zc,0)
               rc2 = xc*xc + yc*yc + zc*zc               
                  xr = xc + xi - x(kk) + x(kk)
                  yr = yc + yi - y(kk) + y(kk)
                  zr = zc + zi - z(kk) + z(kk)
                  r2 = xr*xr + yr*yr + zr*zr
                  r = sqrt(r2)
                  fik = fi * pchg(kk) * cscale(kk)

        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                  fik =  fik*lamp2
        endif
                  e = fik / r
         if (molcule(PIEatom(in)) .ne. molcule(PIEatom(kk))) then
                 einter = einter + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).ne.obj).or.
     &(molcule(in).ne.obj.and.molcule(kk).eq.obj))then
                 eintero = eintero + e
         end if
        if ((molcule(in).eq.obj.and.molcule(kk).eq.obj))then
                 eintrao = eintrao + e
         end if
           if(abs(e).gt.800)then
           write(*,*) ii,kk,molcule(II),molcule(kk)
           endif               
c                    write(11,*) 'pie-pie',e                                                                                                                                                
                  ec1 = ec1 + e   
         end do
      end do
c         write(12,*) 'pie',ec1   
c----------------------------------------
       endif
c 以下为计算力的程序不需要并行(由于字数要求已省略)
c---------------------------------------
       ecry=ecry+ec1
        if (.not. use_replica)  goto 111   
c*****calculate the interaction with other cells with abeem charges explictly dielec
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

czhui11

★
sunxiao(金币+1,VIP+0):谢谢参与,欢迎常来仿真编程版 5-19 02:14
你怎么发这么长啊,谁有那么长时间去看啊
2楼2009-05-18 20:44:19
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

snoopyzhao

至尊木虫 (职业作家)

★ ★
sunxiao(金币+2,VIP+0):谢谢参与,欢迎常来仿真编程版 5-19 02:14
不仅长,而且没有注释,再有,也不能编译执行,呵呵……
3楼2009-05-18 22:48:40
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mickk

铁杆木虫 (职业作家)

★
gwdavid(金币+1,VIP+0):欢迎常来 5-31 19:49
不仅长,而且好像不完整,没法调试呀
4楼2009-05-30 23:09:02
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 阳光不锈7717 的主题更新
普通表情 龙 兔 虎 猫 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[职场人生] 40+中年人失业后,还能找到好工作吗? +10 xqy2009 2026-10-13 12/600 2026-10-15 19:24 by lpwang
[有机交流] 用NBS溴化时溶剂除了CCl4外还有什么 +3 18834002860 2026-10-12 6/300 2026-10-15 14:30 by 长沙摩尔信息科�
[基金申请] 申请科研项目纯靠实力能成功吗? 200+7 David1982 2026-10-13 11/550 2026-10-15 10:30 by zyqchem
[硕博家园] 极限拉扯,如何抉择 +11 Isak0320 2026-10-12 15/750 2026-10-14 21:55 by 猫鱼不是鱼
[论文投稿] 二审4个多月了还是Awaiting Reviewer Scores 10+3 纳斯达克ML 2026-10-11 4/200 2026-10-14 17:30 by babero
[有机交流] OpenAI又攻克数学难题了,离AI做有机合成还有多远? +8 wangyl_123 2026-10-10 13/650 2026-10-14 15:12 by zuesli
[教师之家] 课题组会帮你通过考核么? +4 好运来~9527 2026-10-10 7/350 2026-10-13 08:41 by xli1984
[公派出国] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +5 qR2AdCWTVtkw 2026-10-09 9/450 2026-10-11 22:34 by XQNyDCa1Tk5n
[公派出国] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +5 qR2AdCWTVtkw 2026-10-09 10/500 2026-10-11 22:26 by XQNyDCa1Tk5n
[找工作] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +5 qR2AdCWTVtkw 2026-10-09 8/400 2026-10-11 22:26 by XQNyDCa1Tk5n
[考博] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +5 qR2AdCWTVtkw 2026-10-09 9/450 2026-10-11 22:15 by XQNyDCa1Tk5n
[考博] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 qR2AdCWTVtkw 2026-10-09 7/350 2026-10-11 22:11 by XQNyDCa1Tk5n
[教师之家] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 (EPI+-1)(金币-50) +5 qR2AdCWTVtkw 2026-10-09 10/500 2026-10-11 14:25 by VOMexwnDzOFZ
[教师之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 (EPI+-1)(金币-50) +4 qR2AdCWTVtkw 2026-10-09 7/350 2026-10-11 14:05 by VOMexwnDzOFZ
[考研] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +5 qR2AdCWTVtkw 2026-10-09 8/400 2026-10-11 13:51 by VOMexwnDzOFZ
[硕博家园] 售SCI-T0P文章,我:8O.5.5.1.O.54,科目齐全,可+急 +5 qR2AdCWTVtkw 2026-10-09 8/400 2026-10-11 13:45 by VOMexwnDzOFZ
[找工作] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 qR2AdCWTVtkw 2026-10-09 7/350 2026-10-11 10:54 by VOMexwnDzOFZ
[找工作] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +4 qR2AdCWTVtkw 2026-10-09 7/350 2026-10-11 10:34 by VOMexwnDzOFZ
[考博] 招收博士生一名 +3 mandolin 2026-10-10 3/150 2026-10-11 10:21 by 不是山谷!
[公派出国] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +4 qR2AdCWTVtkw 2026-10-09 4/200 2026-10-11 03:14 by M0zrYVSAnxpE
信息提示
请填处理意见