Zmiany do WHAM-M dla pola Yi
authorAdam Kazimierz Sieradzan <adasko@sun1.chem.univ.gda.pl>
Mon, 3 Dec 2012 10:43:17 +0000 (05:43 -0500)
committerAdam Kazimierz Sieradzan <adasko@sun1.chem.univ.gda.pl>
Mon, 3 Dec 2012 10:43:17 +0000 (05:43 -0500)
bin/wham_multparm-D-sccor-oldparm
source/wham/src-M/COMMON.ALLPARM
source/wham/src-M/cinfo.f
source/wham/src-M/energy_p_new.F
source/wham/src-M/include_unres/COMMON.LOCAL
source/wham/src-M/parmread.F
source/wham/src-M/store_parm.F

index 2e94320..008abf2 100755 (executable)
Binary files a/bin/wham_multparm-D-sccor-oldparm and b/bin/wham_multparm-D-sccor-oldparm differ
index 8461a62..a65a246 100644 (file)
@@ -9,20 +9,33 @@
      & polthet_all(0:3,-ntyp:ntyp,max_parm),
      & gthet_all(3,-ntyp:ntyp,max_parm),theta0_all(-ntyp:ntyp,max_parm),
      & sig0_all(-ntyp:ntyp,max_parm),sigc0_all(-ntyp:ntyp,max_parm),
-     & aa0thet_all(maxthetyp1,maxthetyp1,maxthetyp1,max_parm),
-     & aathet_all(maxtheterm,maxthetyp1,maxthetyp1,maxthetyp1,max_parm),
-     & bbthet_all(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,
-     & maxthetyp1,max_parm),
-     & ccthet_all(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,
-     & maxthetyp1,max_parm),
-     & ddthet_all(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,
-     & maxthetyp1,max_parm),
-     & eethet_all(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,
-     & maxthetyp1,max_parm),
-     & ffthet_all(maxdouble,maxdouble,maxtheterm3,maxthetyp1,maxthetyp1,
-     &  maxthetyp1,max_parm),
-     & ggthet_all(maxdouble,maxdouble,maxtheterm3,maxthetyp1,maxthetyp1,
-     &  maxthetyp1,max_parm),
+     & aa0thet_all(-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,2,max_parm),
+     & aathet_all(maxtheterm,-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2,max_parm),
+     & bbthet_all(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2,max_parm),
+     & ccthet_all(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,2,max_parm),
+     & ddthet_all(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,
+     & maxthetyp1,2,max_parm),
+     & eethet_all(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,
+     & -maxthetyp1:maxthetyp1,2,max_parm),
+     & ffthet_all1(maxdouble,maxdouble,maxtheterm3,
+     & -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,
+     &  -maxthetyp1:maxthetyp1,max_parm),
+     & ggthet_all1(maxdouble,maxdouble,maxtheterm3,
+     &  -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,
+     &  -maxthetyp1:maxthetyp1,max_parm),
+     & ffthet_all2(maxdouble,maxdouble,maxtheterm3,
+     & -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,
+     &  -maxthetyp1:maxthetyp1,max_parm),
+     & ggthet_all2(maxdouble,maxdouble,maxtheterm3,
+     &  -maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,
+     &  -maxthetyp1:maxthetyp1,max_parm),
      & dsc_all(ntyp1,max_parm),bsc_all(maxlob,ntyp,max_parm),
      & censc_all(3,maxlob,-ntyp:ntyp,max_parm),
      & gaussc_all(3,3,maxlob,-ntyp:ntyp,max_parm),
@@ -77,7 +90,8 @@
      &nterm_sccor_all(-ntyp:ntyp,-ntyp:ntyp,max_parm)
       common /allparm/ ww_all,vbldp0_all,akp_all,vbldsc0_all,aksc_all,
      & abond0_all,aa0thet_all,aathet_all,bbthet_all,ccthet_all,
-     & ddthet_all,eethet_all,ffthet_all,ggthet_all,
+     & ddthet_all,eethet_all,ffthet_all1,ggthet_all1,
+     & ffthet_all2,ggthet_all2,
      & a0thet_all,athet_all,bthet_all,polthet_all,gthet_all,theta0_all,
      & sig0_all,sigc0_all,dsc_all,bsc_all,censc_all,gaussc_all,dsc0_all,
      & sc_parmin_all,
index 0672985..d497941 100644 (file)
@@ -1,10 +1,10 @@
 C DO NOT EDIT THIS FILE - IT HAS BEEN GENERATED BY COMPINFO.C
-C 0 0 700
+C 0 0 702
       subroutine cinfo
       include 'COMMON.IOUNITS'
       write(iout,*)'++++ Compile info ++++'
-      write(iout,*)'Version 0.0 build 700'
-      write(iout,*)'compiled Tue Oct 30 05:53:04 2012'
+      write(iout,*)'Version 0.0 build 702'
+      write(iout,*)'compiled Mon Dec  3 05:37:30 2012'
       write(iout,*)'compiled by aks255@matrix.chem.cornell.edu'
       write(iout,*)'OS name:    Linux '
       write(iout,*)'OS release: 2.6.34.9-69.fc13.x86_64 '
index 93ae4b9..f9dbc12 100644 (file)
 c      write (iout,*) "ithetyp",(ithetyp(i),i=1,ntyp1)
       do i=ithet_start,ithet_end
         if (itype(i-1).eq.ntyp1) cycle
+        if (iabs(itype(i+1)).eq.20) iblock=2
+        if (iabs(itype(i+1)).ne.20) iblock=1
         dethetai=0.0d0
         dephii=0.0d0
         dephii1=0.0d0
         theti2=0.5d0*theta(i)
-        ityp2=ithetyp(iabs(itype(i-1)))
+        ityp2=ithetyp((itype(i-1)))
         do k=1,nntheterm
           coskt(k)=dcos(k*theti2)
           sinkt(k)=dsin(k*theti2)
@@ -3441,7 +3443,7 @@ c      write (iout,*) "ithetyp",(ithetyp(i),i=1,ntyp1)
 #else
           phii=phi(i)
 #endif
-          ityp1=ithetyp(iabs(itype(i-2)))
+          ityp1=ithetyp((itype(i-2)))
           do k=1,nsingle
             cosph1(k)=dcos(k*phii)
             sinph1(k)=dsin(k*phii)
@@ -3462,7 +3464,7 @@ c      write (iout,*) "ithetyp",(ithetyp(i),i=1,ntyp1)
 #else
           phii1=phi(i+1)
 #endif
-          ityp3=ithetyp(iabs(itype(i)))
+          ityp3=ithetyp((itype(i)))
           do k=1,nsingle
             cosph2(k)=dcos(k*phii1)
             sinph2(k)=dsin(k*phii1)
@@ -3478,7 +3480,7 @@ c      write (iout,*) "ithetyp",(ithetyp(i),i=1,ntyp1)
 c        write (iout,*) "i",i," ityp1",itype(i-2),ityp1,
 c     &   " ityp2",itype(i-1),ityp2," ityp3",itype(i),ityp3
 c        call flush(iout)
-        ethetai=aa0thet(ityp1,ityp2,ityp3)
+        ethetai=aa0thet(ityp1,ityp2,ityp3,iblock)
         do k=1,ndouble
           do l=1,k-1
             ccl=cosph1(l)*cosph2(k-l)
@@ -3500,11 +3502,12 @@ c        call flush(iout)
         enddo
         endif
         do k=1,ntheterm
-          ethetai=ethetai+aathet(k,ityp1,ityp2,ityp3)*sinkt(k)
-          dethetai=dethetai+0.5d0*k*aathet(k,ityp1,ityp2,ityp3)
+          ethetai=ethetai+aathet(k,ityp1,ityp2,ityp3,iblock)*sinkt(k)
+          dethetai=dethetai+0.5d0*k*aathet(k,ityp1,ityp2,ityp3,iblock)
      &      *coskt(k)
           if (lprn)
-     &    write (iout,*) "k",k," aathet",aathet(k,ityp1,ityp2,ityp3),
+     &    write (iout,*) "k",k,"
+     &      aathet",aathet(k,ityp1,ityp2,ityp3,iblock),
      &     " ethetai",ethetai
         enddo
         if (lprn) then
@@ -3523,24 +3526,24 @@ c        call flush(iout)
         endif
         do m=1,ntheterm2
           do k=1,nsingle
-            aux=bbthet(k,m,ityp1,ityp2,ityp3)*cosph1(k)
-     &         +ccthet(k,m,ityp1,ityp2,ityp3)*sinph1(k)
-     &         +ddthet(k,m,ityp1,ityp2,ityp3)*cosph2(k)
-     &         +eethet(k,m,ityp1,ityp2,ityp3)*sinph2(k)
+            aux=bbthet(k,m,ityp1,ityp2,ityp3,iblock)*cosph1(k)
+     &         +ccthet(k,m,ityp1,ityp2,ityp3,iblock)*sinph1(k)
+     &         +ddthet(k,m,ityp1,ityp2,ityp3,iblock)*cosph2(k)
+     &         +eethet(k,m,ityp1,ityp2,ityp3,iblock)*sinph2(k)
             ethetai=ethetai+sinkt(m)*aux
             dethetai=dethetai+0.5d0*m*aux*coskt(m)
             dephii=dephii+k*sinkt(m)*(
-     &          ccthet(k,m,ityp1,ityp2,ityp3)*cosph1(k)-
-     &          bbthet(k,m,ityp1,ityp2,ityp3)*sinph1(k))
+     &          ccthet(k,m,ityp1,ityp2,ityp3,iblock)*cosph1(k)-
+     &          bbthet(k,m,ityp1,ityp2,ityp3,iblock)*sinph1(k))
             dephii1=dephii1+k*sinkt(m)*(
-     &          eethet(k,m,ityp1,ityp2,ityp3)*cosph2(k)-
-     &          ddthet(k,m,ityp1,ityp2,ityp3)*sinph2(k))
+     &          eethet(k,m,ityp1,ityp2,ityp3,iblock)*cosph2(k)-
+     &          ddthet(k,m,ityp1,ityp2,ityp3,iblock)*sinph2(k))
             if (lprn)
      &      write (iout,*) "m",m," k",k," bbthet",
-     &         bbthet(k,m,ityp1,ityp2,ityp3)," ccthet",
-     &         ccthet(k,m,ityp1,ityp2,ityp3)," ddthet",
-     &         ddthet(k,m,ityp1,ityp2,ityp3)," eethet",
-     &         eethet(k,m,ityp1,ityp2,ityp3)," ethetai",ethetai
+     &         bbthet(k,m,ityp1,ityp2,ityp3,iblock)," ccthet",
+     &         ccthet(k,m,ityp1,ityp2,ityp3,iblock)," ddthet",
+     &         ddthet(k,m,ityp1,ityp2,ityp3,iblock)," eethet",
+     &         eethet(k,m,ityp1,ityp2,ityp3,iblock)," ethetai",ethetai
           enddo
         enddo
         if (lprn)
@@ -3548,28 +3551,29 @@ c        call flush(iout)
         do m=1,ntheterm3
           do k=2,ndouble
             do l=1,k-1
-              aux=ffthet(l,k,m,ityp1,ityp2,ityp3)*cosph1ph2(l,k)+
-     &            ffthet(k,l,m,ityp1,ityp2,ityp3)*cosph1ph2(k,l)+
-     &            ggthet(l,k,m,ityp1,ityp2,ityp3)*sinph1ph2(l,k)+
-     &            ggthet(k,l,m,ityp1,ityp2,ityp3)*sinph1ph2(k,l)
+              aux=ffthet(l,k,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(l,k)+
+     &            ffthet(k,l,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(k,l)+
+     &            ggthet(l,k,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(l,k)+
+     &            ggthet(k,l,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(k,l)
               ethetai=ethetai+sinkt(m)*aux
               dethetai=dethetai+0.5d0*m*coskt(m)*aux
               dephii=dephii+l*sinkt(m)*(
-     &           -ffthet(l,k,m,ityp1,ityp2,ityp3)*sinph1ph2(l,k)-
-     &            ffthet(k,l,m,ityp1,ityp2,ityp3)*sinph1ph2(k,l)+
-     &            ggthet(l,k,m,ityp1,ityp2,ityp3)*cosph1ph2(l,k)+
-     &            ggthet(k,l,m,ityp1,ityp2,ityp3)*cosph1ph2(k,l))
+     &           -ffthet(l,k,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(l,k)-
+     &            ffthet(k,l,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(k,l)+
+     &            ggthet(l,k,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(l,k)+
+     &            ggthet(k,l,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(k,l))
               dephii1=dephii1+(k-l)*sinkt(m)*(
-     &           -ffthet(l,k,m,ityp1,ityp2,ityp3)*sinph1ph2(l,k)+
-     &            ffthet(k,l,m,ityp1,ityp2,ityp3)*sinph1ph2(k,l)+
-     &            ggthet(l,k,m,ityp1,ityp2,ityp3)*cosph1ph2(l,k)-
-     &            ggthet(k,l,m,ityp1,ityp2,ityp3)*cosph1ph2(k,l))
+     &           -ffthet(l,k,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(l,k)+
+     &            ffthet(k,l,m,ityp1,ityp2,ityp3,iblock)*sinph1ph2(k,l)+
+     &            ggthet(l,k,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(l,k)-
+     &            ggthet(k,l,m,ityp1,ityp2,ityp3,iblock)*cosph1ph2(k,l))
               if (lprn) then
               write (iout,*) "m",m," k",k," l",l," ffthet",
-     &            ffthet(l,k,m,ityp1,ityp2,ityp3),
-     &            ffthet(k,l,m,ityp1,ityp2,ityp3)," ggthet",
-     &            ggthet(l,k,m,ityp1,ityp2,ityp3),
-     &            ggthet(k,l,m,ityp1,ityp2,ityp3)," ethetai",ethetai
+     &            ffthet(l,k,m,ityp1,ityp2,ityp3,iblock),
+     &            ffthet(k,l,m,ityp1,ityp2,ityp3,iblock)," ggthet",
+     &            ggthet(l,k,m,ityp1,ityp2,ityp3,iblock),
+     &            ggthet(k,l,m,ityp1,ityp2,ityp3,iblock),
+     &            " ethetai",ethetai
               write (iout,*) cosph1ph2(l,k)*sinkt(m),
      &            cosph1ph2(k,l)*sinkt(m),
      &            sinph1ph2(l,k)*sinkt(m),sinph1ph2(k,l)*sinkt(m)
@@ -3935,7 +3939,7 @@ C     &   dc_norm(3,i+nres)
           y_prime(j) = (dc_norm(j,i) + dc_norm(j,i-1))*sinfac
         enddo
         do j = 1,3
-          z_prime(j) = -uz(j,i-1)
+          z_prime(j) = -uz(j,i-1)*dsign(1.0,dfloat(itype(i)))
         enddo     
 c       write (2,*) "i",i
 c       write (2,*) "x_prime",(x_prime(j),j=1,3)
@@ -3957,7 +3961,7 @@ c
         do j = 1,3
           xx = xx + x_prime(j)*dc_norm(j,i+nres)
           yy = yy + y_prime(j)*dc_norm(j,i+nres)
-          zz = zz + dsign(1.0,itype(i))*z_prime(j)*dc_norm(j,i+nres)
+          zz = zz + z_prime(j)*dc_norm(j,i+nres)
         enddo
 
         xxtab(i)=xx
@@ -4146,8 +4150,10 @@ c     &   (dC_norm(j,i-1),j=1,3)," vbld_inv",vbld_inv(i+1),vbld_inv(i)
          dZZ_Ci1(k)=0.0d0
          dZZ_Ci(k)=0.0d0
          do j=1,3
-           dZZ_Ci(k)=dZZ_Ci(k)-uzgrad(j,k,2,i-1)*dC_norm(j,i+nres)
-           dZZ_Ci1(k)=dZZ_Ci1(k)-uzgrad(j,k,1,i-1)*dC_norm(j,i+nres)
+           dZZ_Ci(k)=dZZ_Ci(k)-uzgrad(j,k,2,i-1)
+     & *dsign(1.0,dfloat(itype(i)))*dC_norm(j,i+nres)
+           dZZ_Ci1(k)=dZZ_Ci1(k)-uzgrad(j,k,1,i-1)
+     &  *dsign(1.0,dfloat(itype(i)))*dC_norm(j,i+nres)
          enddo
           
          dXX_XYZ(k)=vbld_inv(i+nres)*(x_prime(k)-xx*dC_norm(k,i+nres))
index 1594710..748e26a 100644 (file)
@@ -13,17 +13,23 @@ C Parameters of the side-chain probability distribution
      &    nlob(ntyp1)
 C Parameters of ab initio-derived potential of virtual-bond-angle bending
       integer nthetyp,ntheterm,ntheterm2,ntheterm3,nsingle,ndouble,
-     & ithetyp(ntyp1),nntheterm
-      double precision aa0thet(maxthetyp1,maxthetyp1,maxthetyp1),
-     & aathet(maxtheterm,maxthetyp1,maxthetyp1,maxthetyp1),
-     & bbthet(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,maxthetyp1),
-     & ccthet(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,maxthetyp1),
-     & ddthet(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,maxthetyp1),
-     & eethet(maxsingle,maxtheterm2,maxthetyp1,maxthetyp1,maxthetyp1),
-     & ffthet(maxdouble,maxdouble,maxtheterm3,maxthetyp1,maxthetyp1,
-     &  maxthetyp1),
-     & ggthet(maxdouble,maxdouble,maxtheterm3,maxthetyp1,maxthetyp1,
-     &  maxthetyp1)
+     & ithetyp(-ntyp1:ntyp1),nntheterm
+      double precision aa0thet(-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & aathet(maxtheterm,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & bbthet(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & ccthet(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & ddthet(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & eethet(maxsingle,maxtheterm2,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,-maxthetyp1:maxthetyp1,2),
+     & ffthet(maxdouble,maxdouble,maxtheterm3,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1, -maxthetyp1:maxthetyp1,2),
+     & ggthet(maxdouble,maxdouble,maxtheterm3,-maxthetyp1:maxthetyp1,
+     &-maxthetyp1:maxthetyp1,  -maxthetyp1:maxthetyp1,2)
       common /theta_abinitio/aa0thet,aathet,bbthet,ccthet,ddthet,eethet,
      &  ffthet,
      &  ggthet,ithetyp,nthetyp,ntheterm,ntheterm2,ntheterm3,nsingle,
index e033c78..b2b64d5 100644 (file)
@@ -299,46 +299,55 @@ C
      &  ntheterm3,nsingle,ndouble
       nntheterm=max0(ntheterm,ntheterm2,ntheterm3)
       read (ithep,*) (ithetyp(i),i=1,ntyp1)
-      do i=1,maxthetyp
-        do j=1,maxthetyp
-          do k=1,maxthetyp
-            aa0thet(i,j,k)=0.0d0
+      do i=-ntyp1,-1
+        ithetyp(i)=-ithetyp(-i)
+      enddo
+      do iblock=1,2
+      do i=-maxthetyp,maxthetyp
+        do j=-maxthetyp,maxthetyp
+          do k=-maxthetyp,maxthetyp
+            aa0thet(i,j,k,iblock)=0.0d0
             do l=1,ntheterm
-              aathet(l,i,j,k)=0.0d0
+              aathet(l,i,j,k,iblock)=0.0d0
             enddo
             do l=1,ntheterm2
               do m=1,nsingle
-                bbthet(m,l,i,j,k)=0.0d0
-                ccthet(m,l,i,j,k)=0.0d0
-                ddthet(m,l,i,j,k)=0.0d0
-                eethet(m,l,i,j,k)=0.0d0
+                bbthet(m,l,i,j,k,iblock)=0.0d0
+                ccthet(m,l,i,j,k,iblock)=0.0d0
+                ddthet(m,l,i,j,k,iblock)=0.0d0
+                eethet(m,l,i,j,k,iblock)=0.0d0
               enddo
             enddo
             do l=1,ntheterm3
               do m=1,ndouble
                 do mm=1,ndouble
-                 ffthet(mm,m,l,i,j,k)=0.0d0
-                 ggthet(mm,m,l,i,j,k)=0.0d0
+                 ffthet(mm,m,l,i,j,k,iblock)=0.0d0
+                 ggthet(mm,m,l,i,j,k,iblock)=0.0d0
                 enddo
               enddo
             enddo
           enddo
         enddo
       enddo
-      do i=1,nthetyp
-        do j=1,nthetyp
-          do k=1,nthetyp
-            read (ithep,'(3a)') res1,res2,res3
-            read (ithep,*) aa0thet(i,j,k)
-            read (ithep,*)(aathet(l,i,j,k),l=1,ntheterm)
+      enddo
+      do iblock=1,2
+      do i=0,nthetyp
+        do j=-nthetyp,nthetyp
+          do k=-nthetyp,nthetyp
+            read (ithep,'(6a)') res1
+            read (ithep,*) aa0thet(i,j,k,iblock)
+            read (ithep,*)(aathet(l,i,j,k,iblock),l=1,ntheterm)
             read (ithep,*)
-     &       ((bbthet(lll,ll,i,j,k),lll=1,nsingle),
-     &        (ccthet(lll,ll,i,j,k),lll=1,nsingle),
-     &        (ddthet(lll,ll,i,j,k),lll=1,nsingle),
-     &        (eethet(lll,ll,i,j,k),lll=1,nsingle),ll=1,ntheterm2)
+     &       ((bbthet(lll,ll,i,j,k,iblock),lll=1,nsingle),
+     &        (ccthet(lll,ll,i,j,k,iblock),lll=1,nsingle),
+     &        (ddthet(lll,ll,i,j,k,iblock),lll=1,nsingle),
+     &        (eethet(lll,ll,i,j,k,iblock),lll=1,nsingle)
+     &        ,ll=1,ntheterm2)
             read (ithep,*)
-     &      (((ffthet(llll,lll,ll,i,j,k),ffthet(lll,llll,ll,i,j,k),
-     &         ggthet(llll,lll,ll,i,j,k),ggthet(lll,llll,ll,i,j,k),
+     &      (((ffthet(llll,lll,ll,i,j,k,iblock),
+     &      ffthet(lll,llll,ll,i,j,k,iblock),
+     &         ggthet(llll,lll,ll,i,j,k,iblock)
+     &        ,ggthet(lll,llll,ll,i,j,k,iblock),
      &         llll=1,lll-1),lll=2,ndouble),ll=1,ntheterm3)
           enddo
         enddo
@@ -350,17 +359,54 @@ C
       do i=1,nthetyp
         do j=1,nthetyp
           do l=1,ntheterm
-            aathet(l,i,j,nthetyp+1)=aathet(l,i,j,1)
-            aathet(l,nthetyp+1,i,j)=aathet(l,1,i,j)
+            aathet(l,i,j,nthetyp+1,iblock)=0.0d0
+            aathet(l,nthetyp+1,i,j,iblock)=0.0d0
           enddo
-          aa0thet(i,j,nthetyp+1)=aa0thet(i,j,1)
-          aa0thet(nthetyp+1,i,j)=aa0thet(1,i,j)
+          aa0thet(i,j,nthetyp+1,iblock)=0.0d0
+          aa0thet(nthetyp+1,i,j,iblock)=0.0d0
         enddo
         do l=1,ntheterm
-          aathet(l,nthetyp+1,i,nthetyp+1)=aathet(l,1,i,1)
+          aathet(l,nthetyp+1,i,nthetyp+1,iblock)=0.0d0
         enddo
-        aa0thet(nthetyp+1,i,nthetyp+1)=aa0thet(1,i,1)
+        aa0thet(nthetyp+1,i,nthetyp+1,iblock)=0.0d0
       enddo
+      enddo
+C Substitution for D aminoacids from symmetry.
+      do iblock=1,2
+      do i=-nthetyp,0
+        do j=-nthetyp,nthetyp
+          do k=-nthetyp,nthetyp
+           aa0thet(i,j,k,iblock)=aa0thet(-i,-j,-k,iblock)
+           do l=1,ntheterm
+           aathet(l,i,j,k,iblock)=aathet(l,-i,-j,-k,iblock)
+           enddo
+           do ll=1,ntheterm2
+            do lll=1,nsingle
+            bbthet(lll,ll,i,j,k,iblock)=bbthet(lll,ll,-i,-j,-k,iblock)
+            ccthet(lll,ll,i,j,k,iblock)=-ccthet(lll,ll,-i,-j,-k,iblock)
+            ddthet(lll,ll,i,j,k,iblock)=ddthet(lll,ll,-i,-j,-k,iblock)
+            eethet(lll,ll,i,j,k,iblock)=-eethet(lll,ll,-i,-j,-k,iblock)
+            enddo
+          enddo
+          do ll=1,ntheterm3
+           do lll=2,ndouble
+            do llll=1,lll-1
+            ffthet(llll,lll,ll,i,j,k,iblock)=
+     &      ffthet(llll,lll,ll,-i,-j,-k,iblock)
+            ffthet(lll,llll,ll,i,j,k,iblock)=
+     &      ffthet(lll,llll,ll,-i,-j,-k,iblock)
+            ggthet(llll,lll,ll,i,j,k,iblock)=
+     &      -ggthet(llll,lll,ll,-i,-j,-k,iblock)
+            ggthet(lll,llll,ll,i,j,k,iblock)=
+     &      -ggthet(lll,llll,ll,-i,-j,-k,iblock)
+            enddo !ll
+           enddo  !lll  
+          enddo   !llll
+         enddo    !k
+        enddo     !j
+       enddo      !i
+      enddo       !iblock
+
 C
 C Control printout of the coefficients of virtual-bond-angle potentials
 C
@@ -372,7 +418,7 @@ C
               write (iout,'(//4a)')
      &         'Type ',onelett(i),onelett(j),onelett(k)
               write (iout,'(//a,10x,a)') " l","a[l]"
-              write (iout,'(i2,1pe15.5)') 0,aa0thet(i,j,k)
+              write (iout,'(i2,1pe15.5)') 0,aa0thet(i,j,k,iblock)
               write (iout,'(i2,1pe15.5)')
      &           (l,aathet(l,i,j,k),l=1,ntheterm)
             do l=1,ntheterm2
@@ -380,8 +426,8 @@ C
      &          "b",l,"c",l,"d",l,"e",l
               do m=1,nsingle
                 write (iout,'(i2,4(1pe15.5))') m,
-     &          bbthet(m,l,i,j,k),ccthet(m,l,i,j,k),
-     &          ddthet(m,l,i,j,k),eethet(m,l,i,j,k)
+     &          bbthet(m,l,i,j,k),ccthet(m,l,i,j,k,iblock),
+     &          ddthet(m,l,i,j,k),eethet(m,l,i,j,k,iblock)
               enddo
             enddo
             do l=1,ntheterm3
@@ -390,8 +436,10 @@ C
               do m=2,ndouble
                 do n=1,m-1
                   write (iout,'(i1,1x,i1,4(1pe15.5))') n,m,
-     &              ffthet(n,m,l,i,j,k),ffthet(m,n,l,i,j,k),
-     &              ggthet(n,m,l,i,j,k),ggthet(m,n,l,i,j,k)
+     &              ffthet(n,m,l,i,j,k,iblock),
+     &              ffthet(m,n,l,i,j,k,iblock),
+     &              ggthet(n,m,l,i,j,k,iblock),
+     &              ggthet(m,n,l,i,j,k,iblock)
                 enddo
               enddo
             enddo
index ad64f47..fa9d1d7 100644 (file)
@@ -81,35 +81,50 @@ c Store bond angle parameters
       nsingle_all(iparm)=nsingle
       ndouble_all(iparm)=ndouble
       nntheterm_all(iparm)=nntheterm
-      do i=1,ntyp1
+      do i=-ntyp,ntyp
         ithetyp_all(i,iparm)=ithetyp(i)
       enddo
-      do i=1,maxthetyp1
-        do j=1,maxthetyp1
-          do k=1,maxthetyp1
-            aa0thet_all(i,j,k,iparm)=aa0thet(i,j,k)
+      do iblock=1,2
+      do i=-maxthetyp1,maxthetyp1
+        do j=-maxthetyp1,maxthetyp1
+          do k=-maxthetyp1,maxthetyp1
+            aa0thet_all(i,j,k,iblock,iparm)=aa0thet(i,j,k,iblock)
             do l=1,ntheterm
-              aathet_all(l,i,j,k,iparm)=aathet(l,i,j,k)
+              aathet_all(l,i,j,k,iblock,iparm)=aathet(l,i,j,k,iblock)
             enddo
             do l=1,ntheterm2
               do m=1,nsingle
-                bbthet_all(m,l,i,j,k,iparm)=bbthet(m,l,i,j,k)
-                ccthet_all(m,l,i,j,k,iparm)=ccthet(m,l,i,j,k)
-                ddthet_all(m,l,i,j,k,iparm)=ddthet(m,l,i,j,k)
-                eethet_all(m,l,i,j,k,iparm)=eethet(m,l,i,j,k)
+                bbthet_all(m,l,i,j,k,iblock,iparm)=
+     & bbthet(m,l,i,j,k,iblock)
+                ccthet_all(m,l,i,j,k,iblock,iparm)=
+     &ccthet(m,l,i,j,k,iblock)
+                ddthet_all(m,l,i,j,k,iblock,iparm)=
+     &ddthet(m,l,i,j,k,iblock)
+                eethet_all(m,l,i,j,k,iblock,iparm)=
+     &eethet(m,l,i,j,k,iblock)
               enddo
             enddo
             do l=1,ntheterm3
               do m=1,ndouble
                 do mm=1,ndouble
-                 ffthet_all(mm,m,l,i,j,k,iparm)=ffthet(mm,m,l,i,j,k)
-                 ggthet_all(mm,m,l,i,j,k,iparm)=ggthet(mm,m,l,i,j,k)
+                if iblock.eq.1 then
+                 ffthet_all1(mm,m,l,i,j,k,iparm)=
+     &   ffthet(mm,m,l,i,j,k,iblock)
+                 ggthet_all1(mm,m,l,i,j,k,iparm)=
+     &ggthet(mm,m,l,i,j,k)
+                  else
+                 ffthet_all2(mm,m,l,i,j,k,iparm)=
+     &   ffthet(mm,m,l,i,j,k,iblock)
+                 ggthet_all2(mm,m,l,i,j,k,iparm)=
+     &ggthet(mm,m,l,i,j,k)
+                  endif
                 enddo
               enddo
             enddo
           enddo
         enddo
       enddo
+      enddo
 #endif
 #ifdef CRYST_SC
 c Store the sidechain rotamer parameters
@@ -337,35 +352,50 @@ c Restore bond angle parameters
       nsingle=nsingle_all(iparm)
       ndouble=ndouble_all(iparm)
       nntheterm=nntheterm_all(iparm)
-      do i=1,ntyp1
+      do i=-ntyp,ntyp
         ithetyp(i)=ithetyp_all(i,iparm)
       enddo
-      do i=1,maxthetyp1
-        do j=1,maxthetyp1
-          do k=1,maxthetyp1
-            aa0thet(i,j,k)=aa0thet_all(i,j,k,iparm)
+      do iblock=1,2
+      do i=-maxthetyp1,maxthetyp1
+        do j=-maxthetyp1,maxthetyp1
+          do k=-maxthetyp1,maxthetyp1
+            aa0thet(i,j,k,iblock)=aa0thet_all(i,j,k,iblock,iparm)
             do l=1,ntheterm
-              aathet(l,i,j,k)=aathet_all(l,i,j,k,iparm)
+              aathet(l,i,j,k,iblock)=aathet_all(l,i,j,k,iblock,iparm)
             enddo
             do l=1,ntheterm2
               do m=1,nsingle
-                bbthet(m,l,i,j,k)=bbthet_all(m,l,i,j,k,iparm)
-                ccthet(m,l,i,j,k)=ccthet_all(m,l,i,j,k,iparm)
-                ddthet(m,l,i,j,k)=ddthet_all(m,l,i,j,k,iparm)
-                eethet(m,l,i,j,k)=eethet_all(m,l,i,j,k,iparm)
+                bbthet(m,l,i,j,k,iblock)=
+     &bbthet_all(m,l,i,j,k,iblock,iparm)
+                ccthet(m,l,i,j,k,iblock)=
+     &ccthet_all(m,l,i,j,k,iblock,iparm)
+                ddthet(m,l,i,j,k,iblock)=
+     &ddthet_all(m,l,i,j,k,iblock,iparm)
+                eethet(m,l,i,j,k,iblock)=
+     &eethet_all(m,l,i,j,k,iblock,iparm)
               enddo
             enddo
             do l=1,ntheterm3
               do m=1,ndouble
                 do mm=1,ndouble
-                 ffthet(mm,m,l,i,j,k)=ffthet_all(mm,m,l,i,j,k,iparm)
-                 ggthet(mm,m,l,i,j,k)=ggthet_all(mm,m,l,i,j,k,iparm)
+                if iblock.eq.1 then
+                 ffthet(mm,m,l,i,j,k,iblock)=
+     &ffthet_all1(mm,m,l,i,j,k,iparm)
+                 ggthet(mm,m,l,i,j,k,iblock)=
+     &ggthet_all1(mm,m,l,i,j,k,iparm)
+                else
+                 ffthet(mm,m,l,i,j,k,iblock)=
+     &ffthet_all2(mm,m,l,i,j,k,iparm)
+                 ggthet(mm,m,l,i,j,k,iblock)=
+     &ggthet_all2(mm,m,l,i,j,k,iparm)
+                endif
                 enddo
               enddo
             enddo
           enddo
         enddo
       enddo
+      enddo
 #endif
 c Restore the sidechain rotamer parameters
 #ifdef CRYST_SC