X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;ds=sidebyside;f=source%2Funres%2Fsrc_MD-M%2Fenergy_p_new_barrier.F;h=252cff03301a20d6c5e389837c9b2eab6d91d9f9;hb=08d089f8cf9bef0534c0fadc67230242fd2903bb;hp=e56e104ce46f8f20ca3675d59eb08e552b98131f;hpb=f1cd18c0a4c39624d6b667327bcc61a873f1bbee;p=unres.git

diff --git a/source/unres/src_MD-M/energy_p_new_barrier.F b/source/unres/src_MD-M/energy_p_new_barrier.F
index e56e104..252cff0 100644
--- a/source/unres/src_MD-M/energy_p_new_barrier.F
+++ b/source/unres/src_MD-M/energy_p_new_barrier.F
@@ -711,7 +711,7 @@ c      enddo
         do i=1,4*nres
           glocbuf(i)=gloc(i,icg)
         enddo
-#define DEBUG
+c#define DEBUG
 #ifdef DEBUG
       write (iout,*) "gloc_sc before reduce"
       do i=1,nres
@@ -720,7 +720,7 @@ c      enddo
        enddo
       enddo
 #endif
-#undef DEBUG
+c#undef DEBUG
         do i=1,nres
          do j=1,3
           gloc_scbuf(j,i)=gloc_sc(j,i,icg)
@@ -740,7 +740,7 @@ c      enddo
         call MPI_Reduce(gloc_scbuf(1,1),gloc_sc(1,1,icg),3*nres,
      &    MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
         time_reduce=time_reduce+MPI_Wtime()-time00
-#define DEBUG
+c#define DEBUG
 #ifdef DEBUG
       write (iout,*) "gloc_sc after reduce"
       do i=1,nres
@@ -749,7 +749,7 @@ c      enddo
        enddo
       enddo
 #endif
-#undef DEBUG
+c#undef DEBUG
 #ifdef DEBUG
       write (iout,*) "gloc after reduce"
       do i=1,4*nres
@@ -1436,30 +1436,39 @@ C we have the original box)
         yi=c(2,nres+i)
         zi=c(3,nres+i)
 C Return atom into box, boxxsize is size of box in x dimension
-  134   continue
-        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
-        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
+c  134   continue
+c        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
+c        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
 C Condition for being inside the proper box
-        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
-     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
-        go to 134
-        endif
-  135   continue
-        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
-        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
+c        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
+c     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
+c        go to 134
+c        endif
+c  135   continue
+c        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
+c        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
 C Condition for being inside the proper box
-        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
-     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
-        go to 135
-        endif
-  136   continue
-        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
-        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
+c        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
+c     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
+c        go to 135
+c        endif
+c  136   continue
+c        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
+c        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
 C Condition for being inside the proper box
-        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
-     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
-        go to 136
-        endif
+c        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
+c     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
+c        go to 136
+c        endif
+          xi=mod(xi,boxxsize)
+          if (xi.lt.0) xi=xi+boxxsize
+          yi=mod(yi,boxysize)
+          if (yi.lt.0) yi=yi+boxysize
+          zi=mod(zi,boxzsize)
+          if (zi.lt.0) zi=zi+boxzsize
+          xi=xi+xshift*boxxsize
+          yi=yi+yshift*boxysize
+          zi=zi+zshift*boxzsize
 
         dxi=dc_norm(1,nres+i)
         dyi=dc_norm(2,nres+i)
@@ -1505,31 +1514,36 @@ c           alf12=0.0D0
             yj=c(2,nres+j)
             zj=c(3,nres+j)
 C Return atom J into box the original box
-  137   continue
-        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
-        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
+c  137   continue
+c        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
+c        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
 C Condition for being inside the proper box
-        if ((xj.gt.((0.5d0)*boxxsize)).or.
-     &       (xj.lt.((-0.5d0)*boxxsize))) then
-        go to 137
-        endif
-  138   continue
-        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
-        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
+c        if ((xj.gt.((0.5d0)*boxxsize)).or.
+c     &       (xj.lt.((-0.5d0)*boxxsize))) then
+c        go to 137
+c        endif
+c  138   continue
+c        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
+c        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
 C Condition for being inside the proper box
-        if ((yj.gt.((0.5d0)*boxysize)).or.
-     &       (yj.lt.((-0.5d0)*boxysize))) then
-        go to 138
-        endif
-  139   continue
-        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
-        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
+c        if ((yj.gt.((0.5d0)*boxysize)).or.
+c     &       (yj.lt.((-0.5d0)*boxysize))) then
+c        go to 138
+c        endif
+c  139   continue
+c        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
+c        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
 C Condition for being inside the proper box
-        if ((zj.gt.((0.5d0)*boxzsize)).or.
-     &       (zj.lt.((-0.5d0)*boxzsize))) then
-        go to 139
-        endif
-
+c        if ((zj.gt.((0.5d0)*boxzsize)).or.
+c     &       (zj.lt.((-0.5d0)*boxzsize))) then
+c        go to 139
+c        endif
+          xj=mod(xj,boxxsize)
+          if (xj.lt.0) xj=xj+boxxsize
+          yj=mod(yj,boxysize)
+          if (yj.lt.0) yj=yj+boxysize
+          zj=mod(zj,boxzsize)
+          if (zj.lt.0) zj=zj+boxzsize
             dxj=dc_norm(1,nres+j)
             dyj=dc_norm(2,nres+j)
             dzj=dc_norm(3,nres+j)
@@ -2871,7 +2885,7 @@ C 14/01/2014 TURN3,TUNR4 does no go under periodic boundry condition
         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1
      &  .or. itype(i+2).eq.ntyp1
      &  .or. itype(i+3).eq.ntyp1
-c     &  .or. itype(i-1).eq.ntyp1
+     &  .or. itype(i-1).eq.ntyp1
      &  .or. itype(i+4).eq.ntyp1
      &  ) cycle
         dxi=dc(1,i)
@@ -2884,30 +2898,36 @@ c     &  .or. itype(i-1).eq.ntyp1
         ymedi=c(2,i)+0.5d0*dyi
         zmedi=c(3,i)+0.5d0*dzi
 C Return atom into box, boxxsize is size of box in x dimension
-  184   continue
-        if (xmedi.gt.((0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
-        if (xmedi.lt.((-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
-C Condition for being inside the proper box
-        if ((xmedi.gt.((0.5d0)*boxxsize)).or.
-     &       (xmedi.lt.((-0.5d0)*boxxsize))) then
-        go to 184
-        endif
-  185   continue
-        if (ymedi.gt.((0.5d0)*boxysize)) ymedi=ymedi-boxysize
-        if (ymedi.lt.((-0.5d0)*boxysize)) ymedi=ymedi+boxysize
+c  184   continue
+c        if (xmedi.gt.((0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
+c        if (xmedi.lt.((-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
 C Condition for being inside the proper box
-        if ((ymedi.gt.((0.5d0)*boxysize)).or.
-     &       (ymedi.lt.((-0.5d0)*boxysize))) then
-        go to 185
-        endif
-  186   continue
-        if (zmedi.gt.((0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
-        if (zmedi.lt.((-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
-C Condition for being inside the proper box
-        if ((zmedi.gt.((0.5d0)*boxzsize)).or.
-     &       (zmedi.lt.((-0.5d0)*boxzsize))) then
-        go to 186
-        endif
+c        if ((xmedi.gt.((0.5d0)*boxxsize)).or.
+c     &       (xmedi.lt.((-0.5d0)*boxxsize))) then
+c        go to 184
+c        endif
+c  185   continue
+c        if (ymedi.gt.((0.5d0)*boxysize)) ymedi=ymedi-boxysize
+c        if (ymedi.lt.((-0.5d0)*boxysize)) ymedi=ymedi+boxysize
+cC Condition for being inside the proper box
+c        if ((ymedi.gt.((0.5d0)*boxysize)).or.
+c     &       (ymedi.lt.((-0.5d0)*boxysize))) then
+c        go to 185
+c        endif
+c  186   continue
+c        if (zmedi.gt.((0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
+c        if (zmedi.lt.((-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
+cC Condition for being inside the proper box
+c        if ((zmedi.gt.((0.5d0)*boxzsize)).or.
+c     &       (zmedi.lt.((-0.5d0)*boxzsize))) then
+c        go to 186
+c        endif
+          xmedi=mod(xmedi,boxxsize)
+          if (xmedi.lt.0) xmedi=xmedi+boxxsize
+          ymedi=mod(ymedi,boxysize)
+          if (ymedi.lt.0) ymedi=ymedi+boxysize
+          zmedi=mod(zmedi,boxzsize)
+          if (zmedi.lt.0) zmedi=zmedi+boxzsize
         num_conti=0
         call eelecij(i,i+2,ees,evdw1,eel_loc)
         if (wturn3.gt.0.0d0) call eturn3(i,eello_turn3)
@@ -2918,6 +2938,8 @@ C Condition for being inside the proper box
      &    .or. itype(i+3).eq.ntyp1
      &    .or. itype(i+4).eq.ntyp1
      &    .or. itype(i+5).eq.ntyp1
+     &    .or. itype(i).eq.ntyp1
+     &    .or. itype(i-1).eq.ntyp1
      &                             ) cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
@@ -2929,30 +2951,36 @@ C Condition for being inside the proper box
         ymedi=c(2,i)+0.5d0*dyi
         zmedi=c(3,i)+0.5d0*dzi
 C Return atom into box, boxxsize is size of box in x dimension
-  194   continue
-        if (xmedi.gt.((0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
-        if (xmedi.lt.((-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
+c  194   continue
+c        if (xmedi.gt.((0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
+c        if (xmedi.lt.((-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
 C Condition for being inside the proper box
-        if ((xmedi.gt.((0.5d0)*boxxsize)).or.
-     &       (xmedi.lt.((-0.5d0)*boxxsize))) then
-        go to 194
-        endif
-  195   continue
-        if (ymedi.gt.((0.5d0)*boxysize)) ymedi=ymedi-boxysize
-        if (ymedi.lt.((-0.5d0)*boxysize)) ymedi=ymedi+boxysize
+c        if ((xmedi.gt.((0.5d0)*boxxsize)).or.
+c     &       (xmedi.lt.((-0.5d0)*boxxsize))) then
+c        go to 194
+c        endif
+c  195   continue
+c        if (ymedi.gt.((0.5d0)*boxysize)) ymedi=ymedi-boxysize
+c        if (ymedi.lt.((-0.5d0)*boxysize)) ymedi=ymedi+boxysize
 C Condition for being inside the proper box
-        if ((ymedi.gt.((0.5d0)*boxysize)).or.
-     &       (ymedi.lt.((-0.5d0)*boxysize))) then
-        go to 195
-        endif
-  196   continue
-        if (zmedi.gt.((0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
-        if (zmedi.lt.((-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
+c        if ((ymedi.gt.((0.5d0)*boxysize)).or.
+c     &       (ymedi.lt.((-0.5d0)*boxysize))) then
+c        go to 195
+c        endif
+c  196   continue
+c        if (zmedi.gt.((0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
+c        if (zmedi.lt.((-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
 C Condition for being inside the proper box
-        if ((zmedi.gt.((0.5d0)*boxzsize)).or.
-     &       (zmedi.lt.((-0.5d0)*boxzsize))) then
-        go to 196
-        endif
+c        if ((zmedi.gt.((0.5d0)*boxzsize)).or.
+c     &       (zmedi.lt.((-0.5d0)*boxzsize))) then
+c        go to 196
+c        endif
+          xmedi=mod(xmedi,boxxsize)
+          if (xmedi.lt.0) xmedi=xmedi+boxxsize
+          ymedi=mod(ymedi,boxysize)
+          if (ymedi.lt.0) ymedi=ymedi+boxysize
+          zmedi=mod(zmedi,boxzsize)
+          if (zmedi.lt.0) zmedi=zmedi+boxzsize
 
         num_conti=num_cont_hb(i)
         call eelecij(i,i+3,ees,evdw1,eel_loc)
@@ -2970,6 +2998,7 @@ c
       do i=iatel_s,iatel_e
         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1
      &  .or. itype(i+2).eq.ntyp1
+     &  .or. itype(i-1).eq.ntyp1
      &                ) cycle
         dxi=dc(1,i)
         dyi=dc(2,i)
@@ -2980,31 +3009,41 @@ c
         xmedi=c(1,i)+0.5d0*dxi
         ymedi=c(2,i)+0.5d0*dyi
         zmedi=c(3,i)+0.5d0*dzi
-C Return atom into box, boxxsize is size of box in x dimension
-  164   continue
-        if (xmedi.gt.((xshift+0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
-        if (xmedi.lt.((xshift-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
-C Condition for being inside the proper box
-        if ((xmedi.gt.((xshift+0.5d0)*boxxsize)).or.
-     &       (xmedi.lt.((xshift-0.5d0)*boxxsize))) then
-        go to 164
-        endif
-  165   continue
-        if (ymedi.gt.((yshift+0.5d0)*boxysize)) ymedi=ymedi-boxysize
-        if (ymedi.lt.((yshift-0.5d0)*boxysize)) ymedi=ymedi+boxysize
+          xmedi=mod(xmedi,boxxsize)
+          if (xmedi.lt.0) xmedi=xmedi+boxxsize
+          ymedi=mod(ymedi,boxysize)
+          if (ymedi.lt.0) ymedi=ymedi+boxysize
+          zmedi=mod(zmedi,boxzsize)
+          if (zmedi.lt.0) zmedi=zmedi+boxzsize
+          xmedi=xmedi+xshift*boxxsize
+          ymedi=ymedi+yshift*boxysize
+          zmedi=zmedi+zshift*boxzsize
+
+C Return tom into box, boxxsize is size of box in x dimension
+c  164   continue
+c        if (xmedi.gt.((xshift+0.5d0)*boxxsize)) xmedi=xmedi-boxxsize
+c        if (xmedi.lt.((xshift-0.5d0)*boxxsize)) xmedi=xmedi+boxxsize
 C Condition for being inside the proper box
-        if ((ymedi.gt.((yshift+0.5d0)*boxysize)).or.
-     &       (ymedi.lt.((yshift-0.5d0)*boxysize))) then
-        go to 165
-        endif
-  166   continue
-        if (zmedi.gt.((zshift+0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
-        if (zmedi.lt.((zshift-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
+c        if ((xmedi.gt.((xshift+0.5d0)*boxxsize)).or.
+c     &       (xmedi.lt.((xshift-0.5d0)*boxxsize))) then
+c        go to 164
+c        endif
+c  165   continue
+c        if (ymedi.gt.((yshift+0.5d0)*boxysize)) ymedi=ymedi-boxysize
+c        if (ymedi.lt.((yshift-0.5d0)*boxysize)) ymedi=ymedi+boxysize
 C Condition for being inside the proper box
-        if ((zmedi.gt.((zshift+0.5d0)*boxzsize)).or.
-     &       (zmedi.lt.((zshift-0.5d0)*boxzsize))) then
-        go to 166
-        endif
+c        if ((ymedi.gt.((yshift+0.5d0)*boxysize)).or.
+c     &       (ymedi.lt.((yshift-0.5d0)*boxysize))) then
+c        go to 165
+c        endif
+c  166   continue
+c        if (zmedi.gt.((zshift+0.5d0)*boxzsize)) zmedi=zmedi-boxzsize
+c        if (zmedi.lt.((zshift-0.5d0)*boxzsize)) zmedi=zmedi+boxzsize
+cC Condition for being inside the proper box
+c        if ((zmedi.gt.((zshift+0.5d0)*boxzsize)).or.
+c     &       (zmedi.lt.((zshift-0.5d0)*boxzsize))) then
+c        go to 166
+c        endif
 
 c        write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
         num_conti=num_cont_hb(i)
@@ -3012,6 +3051,7 @@ c        write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
 c          write (iout,*) i,j,itype(i),itype(j)
           if (itype(j).eq.ntyp1.or. itype(j+1).eq.ntyp1
      & .or.itype(j+2).eq.ntyp1
+     & .or.itype(j-1).eq.ntyp1
      &) cycle
           call eelecij(i,j,ees,evdw1,eel_loc)
         enddo ! j
@@ -3092,31 +3132,38 @@ C          zj=c(3,j)+0.5D0*dzj-zmedi
           xj=c(1,j)+0.5D0*dxj
           yj=c(2,j)+0.5D0*dyj
           zj=c(3,j)+0.5D0*dzj
+          xj=mod(xj,boxxsize)
+          if (xj.lt.0) xj=xj+boxxsize
+          yj=mod(yj,boxysize)
+          if (yj.lt.0) yj=yj+boxysize
+          zj=mod(zj,boxzsize)
+          if (zj.lt.0) zj=zj+boxzsize
+
 C        if ((i+3).lt.j) then !this condition keeps for turn3 and turn4 not subject to PBC
-  174   continue
-        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
-        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
+c  174   continue
+c        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
+c        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
 C Condition for being inside the proper box
-        if ((xj.gt.((0.5d0)*boxxsize)).or.
-     &       (xj.lt.((-0.5d0)*boxxsize))) then
-        go to 174
-        endif
-  175   continue
-        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
-        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
+c        if ((xj.gt.((0.5d0)*boxxsize)).or.
+c     &       (xj.lt.((-0.5d0)*boxxsize))) then
+c        go to 174
+c        endif
+c  175   continue
+c        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
+c        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
 C Condition for being inside the proper box
-        if ((yj.gt.((0.5d0)*boxysize)).or.
-     &       (yj.lt.((-0.5d0)*boxysize))) then
-        go to 175
-        endif
-  176   continue
-        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
-        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
+c        if ((yj.gt.((0.5d0)*boxysize)).or.
+c     &       (yj.lt.((-0.5d0)*boxysize))) then
+c        go to 175
+c        endif
+c  176   continue
+c        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
+c        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
 C Condition for being inside the proper box
-        if ((zj.gt.((0.5d0)*boxzsize)).or.
-     &       (zj.lt.((-0.5d0)*boxzsize))) then
-        go to 176
-        endif
+c        if ((zj.gt.((0.5d0)*boxzsize)).or.
+c     &       (zj.lt.((-0.5d0)*boxzsize))) then
+c        go to 176
+c        endif
 C        endif !endPBC condintion
         xj=xj-xmedi
         yj=yj-ymedi
@@ -4073,30 +4120,39 @@ cd    write (iout,*) 'iatscp_s=',iatscp_s,' iatscp_e=',iatscp_e
         yi=0.5D0*(c(2,i)+c(2,i+1))
         zi=0.5D0*(c(3,i)+c(3,i+1))
 C Return atom into box, boxxsize is size of box in x dimension
-  134   continue
-        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
-        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
-C Condition for being inside the proper box
-        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
-     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
-        go to 134
-        endif
-  135   continue
-        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
-        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
+c  134   continue
+c        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
+c        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
 C Condition for being inside the proper box
-        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
-     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
-        go to 135
-        endif
-  136   continue
-        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
-        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
+c        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
+c     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
+c        go to 134
+c        endif
+c  135   continue
+c        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
+c        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
 C Condition for being inside the proper box
-        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
-     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
-        go to 136
-        endif
+c        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
+c     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
+c        go to 135
+c c       endif
+c  136   continue
+c        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
+c        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
+cC Condition for being inside the proper box
+c        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
+c     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
+c        go to 136
+c        endif
+          xi=mod(xi,boxxsize)
+          if (xi.lt.0) xi=xi+boxxsize
+          yi=mod(yi,boxysize)
+          if (yi.lt.0) yi=yi+boxysize
+          zi=mod(zi,boxzsize)
+          if (zi.lt.0) zi=zi+boxzsize
+          xi=xi+xshift*boxxsize
+          yi=yi+yshift*boxysize
+          zi=zi+zshift*boxzsize
         do iint=1,nscp_gr(i)
 
         do j=iscpstart(i,iint),iscpend(i,iint)
@@ -4110,30 +4166,36 @@ C Uncomment following three lines for Ca-p interactions
           xj=c(1,j)
           yj=c(2,j)
           zj=c(3,j)
-  174   continue
-        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
-        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
-C Condition for being inside the proper box
-        if ((xj.gt.((0.5d0)*boxxsize)).or.
-     &       (xj.lt.((-0.5d0)*boxxsize))) then
-        go to 174
-        endif
-  175   continue
-        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
-        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
+c  174   continue
+c        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
+c        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
 C Condition for being inside the proper box
-        if ((yj.gt.((0.5d0)*boxysize)).or.
-     &       (yj.lt.((-0.5d0)*boxysize))) then
-        go to 175
-        endif
-  176   continue
-        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
-        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
+c        if ((xj.gt.((0.5d0)*boxxsize)).or.
+c     &       (xj.lt.((-0.5d0)*boxxsize))) then
+c        go to 174
+c        endif
+c  175   continue
+c        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
+c        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
+cC Condition for being inside the proper box
+c        if ((yj.gt.((0.5d0)*boxysize)).or.
+c     &       (yj.lt.((-0.5d0)*boxysize))) then
+c        go to 175
+c        endif
+c  176   continue
+c        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
+c        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
 C Condition for being inside the proper box
-        if ((zj.gt.((0.5d0)*boxzsize)).or.
-     &       (zj.lt.((-0.5d0)*boxzsize))) then
-        go to 176
-        endif
+c        if ((zj.gt.((0.5d0)*boxzsize)).or.
+c     &       (zj.lt.((-0.5d0)*boxzsize))) then
+c        go to 176
+          xj=mod(xj,boxxsize)
+          if (xj.lt.0) xj=xj+boxxsize
+          yj=mod(yj,boxysize)
+          if (yj.lt.0) yj=yj+boxysize
+          zj=mod(zj,boxzsize)
+          if (zj.lt.0) zj=zj+boxzsize
+c c       endif
           xj=xj-xi
           yj=yj-yi
           zj=zj-zi
@@ -4227,31 +4289,42 @@ cd    write (iout,*) 'iatscp_s=',iatscp_s,' iatscp_e=',iatscp_e
         xi=0.5D0*(c(1,i)+c(1,i+1))
         yi=0.5D0*(c(2,i)+c(2,i+1))
         zi=0.5D0*(c(3,i)+c(3,i+1))
+          xi=mod(xi,boxxsize)
+          if (xi.lt.0) xi=xi+boxxsize
+          yi=mod(yi,boxysize)
+          if (yi.lt.0) yi=yi+boxysize
+          zi=mod(zi,boxzsize)
+          if (zi.lt.0) zi=zi+boxzsize
+          xi=xi+xshift*boxxsize
+          yi=yi+yshift*boxysize
+          zi=zi+zshift*boxzsize
 C Return atom into box, boxxsize is size of box in x dimension
-  134   continue
-        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
-        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
+c  134   continue
+c        if (xi.gt.((xshift+0.5d0)*boxxsize)) xi=xi-boxxsize
+c        if (xi.lt.((xshift-0.5d0)*boxxsize)) xi=xi+boxxsize
 C Condition for being inside the proper box
-        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
-     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
-        go to 134
-        endif
-  135   continue
-        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
-        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
+c        if ((xi.gt.((xshift+0.5d0)*boxxsize)).or.
+c     &       (xi.lt.((xshift-0.5d0)*boxxsize))) then
+c        go to 134
+c        endif
+c  135   continue
+c          print *,xi,boxxsize,"pierwszy"
+
+c        if (yi.gt.((yshift+0.5d0)*boxysize)) yi=yi-boxysize
+c        if (yi.lt.((yshift-0.5d0)*boxysize)) yi=yi+boxysize
 C Condition for being inside the proper box
-        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
-     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
-        go to 135
-        endif
-  136   continue
-        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
-        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
+c        if ((yi.gt.((yshift+0.5d0)*boxysize)).or.
+c     &       (yi.lt.((yshift-0.5d0)*boxysize))) then
+c        go to 135
+c        endif
+c  136   continue
+c        if (zi.gt.((zshift+0.5d0)*boxzsize)) zi=zi-boxzsize
+c        if (zi.lt.((zshift-0.5d0)*boxzsize)) zi=zi+boxzsize
 C Condition for being inside the proper box
-        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
-     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
-        go to 136
-        endif
+c        if ((zi.gt.((zshift+0.5d0)*boxzsize)).or.
+c     &       (zi.lt.((zshift-0.5d0)*boxzsize))) then
+c        go to 136
+c        endif
         do iint=1,nscp_gr(i)
 
         do j=iscpstart(i,iint),iscpend(i,iint)
@@ -4265,30 +4338,36 @@ C Uncomment following three lines for Ca-p interactions
           xj=c(1,j)
           yj=c(2,j)
           zj=c(3,j)
-  174   continue
-        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
-        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
+          xj=mod(xj,boxxsize)
+          if (xj.lt.0) xj=xj+boxxsize
+          yj=mod(yj,boxysize)
+          if (yj.lt.0) yj=yj+boxysize
+          zj=mod(zj,boxzsize)
+          if (zj.lt.0) zj=zj+boxzsize
+c  174   continue
+c        if (xj.gt.((0.5d0)*boxxsize)) xj=xj-boxxsize
+c        if (xj.lt.((-0.5d0)*boxxsize)) xj=xj+boxxsize
 C Condition for being inside the proper box
-        if ((xj.gt.((0.5d0)*boxxsize)).or.
-     &       (xj.lt.((-0.5d0)*boxxsize))) then
-        go to 174
-        endif
-  175   continue
-        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
-        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
-C Condition for being inside the proper box
-        if ((yj.gt.((0.5d0)*boxysize)).or.
-     &       (yj.lt.((-0.5d0)*boxysize))) then
-        go to 175
-        endif
-  176   continue
-        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
-        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
+c        if ((xj.gt.((0.5d0)*boxxsize)).or.
+c     &       (xj.lt.((-0.5d0)*boxxsize))) then
+c        go to 174
+c        endif
+c  175   continue
+c        if (yj.gt.((0.5d0)*boxysize)) yj=yj-boxysize
+c        if (yj.lt.((-0.5d0)*boxysize)) yj=yj+boxysize
+cC Condition for being inside the proper box
+c        if ((yj.gt.((0.5d0)*boxysize)).or.
+c     &       (yj.lt.((-0.5d0)*boxysize))) then
+c        go to 175
+c        endif
+c  176   continue
+c        if (zj.gt.((0.5d0)*boxzsize)) zj=zj-boxzsize
+c        if (zj.lt.((-0.5d0)*boxzsize)) zj=zj+boxzsize
 C Condition for being inside the proper box
-        if ((zj.gt.((0.5d0)*boxzsize)).or.
-     &       (zj.lt.((-0.5d0)*boxzsize))) then
-        go to 176
-        endif
+c        if ((zj.gt.((0.5d0)*boxzsize)).or.
+c     &       (zj.lt.((-0.5d0)*boxzsize))) then
+c        go to 176
+c        endif
           xj=xj-xi
           yj=yj-yi
           zj=zj-zi