introduction of PBC to prot_cat
authorAdam Sieradzan <adasko@piasek4.chem.univ.gda.pl>
Thu, 27 Sep 2018 11:06:12 +0000 (13:06 +0200)
committerAdam Sieradzan <adasko@piasek4.chem.univ.gda.pl>
Thu, 27 Sep 2018 11:06:12 +0000 (13:06 +0200)
source/unres/data/names.F90
source/unres/energy.F90
source/unres/io_base.F90

index c2d8b58..f901fc1 100644 (file)
@@ -15,7 +15,7 @@
 !-----------------------------------------------------------------------------
 !      block data nazwy
 !el      allocate(restyp(-ntyp1:ntyp1))        !(-ntyp1:ntyp1)
-        character(len=3),dimension(-ntyp1:ntyp1,maxmolec) :: restyp = reshape(&
+        character(len=3),dimension(-ntyp1:ntyp1,maxmolec) :: restyp = (&
         (/'DD ','DDX','DDY','DDZ','DAU','DAI','DDB','DSM','DPR','DLY', &
           'DAR','DHI','DAS',&
        'DGL','DSG','DGN','DSN','DTH',&
        '   ','   ','   ','   ','   ','   ',&
        'A  ','G  ','C  ','T  ','U  ','X  ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ',&
+       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ', &
        '   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ',&
+       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ', &
        '   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ',&
+       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ', &
        '   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ',&
        'NA+','MG2','K+ ','CA2','CL-','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
-       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
-       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
-       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
-       '   ','   ','   ','   ','   ','   ','   ','   ','   ','   ',&
        '   ','   ','   ','   ','   ','   ','   ','   ','   ','   '&
-       /),(/ntyp11,maxmolec/))
+       /))
 !el      allocate(onelet(-ntyp1:ntyp1))         !(-ntyp1:ntyp1)
         character(len=1),dimension(-ntyp1:ntyp1) :: onelet = reshape(&
         (/'z','z','z','z','z','z','z','z',&
index 6727cd8..4a763c2 100644 (file)
             vcatprm(k)=catprm(k,inum)
             enddo
             dASGL=catprm(7,inum)
-             do k=1,3
-                vcm(k)=(cm1(k)/cm1mag)*dASGL+c(k,i+nres)
-                valpha(k)=c(k,i)
-                vcat(k)=c(k,j)
-              enddo
-                      do k=1,3
+!             do k=1,3
+!                vcm(k)=(cm1(k)/cm1mag)*dASGL+c(k,i+nres)
+                vcm(1)=(cm1(1)/cm1mag)*dASGL+xi
+                vcm(2)=(cm1(2)/cm1mag)*dASGL+yi
+                vcm(3)=(cm1(3)/cm1mag)*dASGL+zi
+
+!                valpha(k)=c(k,i)
+!                vcat(k)=c(k,j)
+                if (subchap.eq.1) then
+                 vcat(1)=xj_temp
+                 vcat(2)=yj_temp
+                 vcat(3)=zj_temp
+                 else
+                vcat(1)=xj_safe
+                vcat(2)=yj_safe
+                vcat(3)=zj_safe
+                 endif
+                valpha(1)=xi-c(1,i+nres)+c(1,i)
+                valpha(2)=yi-c(2,i+nres)+c(2,i)
+                valpha(3)=zi-c(3,i+nres)+c(3,i)
+
+!              enddo
+        do k=1,3
           dx(k) = vcat(k)-vcm(k)
         enddo
         do k=1,3
         wquad1=wquad1/wh2o
         wquad2 = vcatprm(4)
         wquad2=wquad2/wh2o
-        wquad2p = 1-wquad2
+        wquad2p = 1.0d0-wquad2
         wvan1 = vcatprm(5)
         wvan2 =vcatprm(6)
         opt = dx(1)**2+dx(2)**2
         rsixp = rfourp*rsecp
         reight=rsixp*rsecp
         Ir = 1.0d0/rs
-        Irsecp = 1/rsecp
+        Irsecp = 1.0d0/rsecp
         Irthrp = Irsecp/rs
         Irfourp = Irthrp/rs
-        Irsixp = 1/rsixp
-        Ireight=1/reight
+        Irsixp = 1.0d0/rsixp
+        Ireight=1.0d0/reight
         Irtw=Irsixp*Irsixp
         Irthir=Irtw/rs
         Irfourt=Irthir/rs
             vcatprm(k)=catprm(k,inum)
             enddo
             dASGL=catprm(7,inum)
-             do k=1,3
-                vcm(k)=(cm1(k)/cm1mag)*dASGL+c(k,i+nres)
-                valpha(k)=c(k,i)
-                vcat(k)=c(k,j)
-              enddo
+!             do k=1,3
+!                vcm(k)=(cm1(k)/cm1mag)*dASGL+c(k,i+nres)
+!                valpha(k)=c(k,i)
+!                vcat(k)=c(k,j)
+!              enddo
+                vcm(1)=(cm1(1)/cm1mag)*dASGL+xi
+                vcm(2)=(cm1(2)/cm1mag)*dASGL+yi
+                vcm(3)=(cm1(3)/cm1mag)*dASGL+zi
+                if (subchap.eq.1) then
+                 vcat(1)=xj_temp
+                 vcat(2)=yj_temp
+                 vcat(3)=zj_temp
+                 else
+                vcat(1)=xj_safe
+                vcat(2)=yj_safe
+                vcat(3)=zj_safe
+                endif
+                valpha(1)=xi-c(1,i+nres)+c(1,i)
+                valpha(2)=yi-c(2,i+nres)+c(2,i)
+                valpha(3)=zi-c(3,i+nres)+c(3,i)
+
 
         do k=1,3
           dx(k) = vcat(k)-vcm(k)
             dscmag = 0.0d0
             do k=1,3
               dscvec(k) = c(k,i+nres)-c(k,i)
+! TU SPRAWDZ???
+!              dscvec(1) = xj
+!              dscvec(2) = yj
+!              dscvec(3) = zj
+
               dscmag = dscmag+dscvec(k)*dscvec(k)
             enddo
             dscmag3 = dscmag
            else
             rcal = 0.0d0
             do k=1,3
-              r(k) = c(k,j)-c(k,i+nres)
+!              r(k) = c(k,j)-c(k,i+nres)
+              r(1) = xj
+              r(2) = yj
+              r(3) = zj
               rcal = rcal+r(k)*r(k)
             enddo
             ract=sqrt(rcal)
index d07c1f0..e82c308 100644 (file)
        stop
       else if (molecule.eq.5) then
       do i=1,ntyp1_molec(molecule)
-        print *,i,restyp(i,molecule)(1:2)
+        print *,restyp(i,molecule)
+        print *,i,restyp(i,molecule)(1:2),nam(1:2)
         if (ucase(nam(1:2)).eq.restyp(i,molecule)(1:2)) then
           rescode=i
           return