correction in nucleic acid wham and in unres ions
[unres4.git] / source / unres / geometry.F90
index ae86ead..c89442f 100644 (file)
 !        print *,i,vbld(i),"vbld(i)"
         vbld_inv(i)=1.0d0/vbld(i)
         vbld(nres+i)=dist(nres+i,i)
-        if (itype(i,1).ne.10) then
+        if ((itype(i,1).ne.10).and.(molnum(i).ne.5)) then
           vbld_inv(nres+i)=1.0d0/vbld(nres+i)
         else
           vbld_inv(nres+i)=0.0d0
         dxi=dc_norm(1,nres+i)
         dyi=dc_norm(2,nres+i)
         dzi=dc_norm(3,nres+i)
+        print *,i,itypi,"sc_move"
         dsci_inv=dsc_inv(itypi)
 !
        do iint=1,nint_gr(i)
          if (itype(j,molnum(j)).eq.ntyp1_molec(molnum(j))) cycle
             ind=ind+1
             itypj=iabs(itype(j,molnum(j)))
+        print *,j,itypj,"sc_move"
+
             dscj_inv=dsc_inv(itypj)
             sig0ij=sigma(itypi,itypj)
             chi1=chi(itypi,itypj)
           +gloc(i-1,icg)*dphi(j,2,i+2)+ &
           gloc(i,icg)*dphi(j,1,i+3)+gloc(nres+i-4,icg)*dtheta(j,2,i+1)+ &
           gloc(nres+i-3,icg)*dtheta(j,1,i+2)
-          if(itype(i,1).ne.10) then
+          if((itype(i,1).ne.10).and.(molnum(nres-1).ne.5)) then
            gcart(j,i)=gcart(j,i)+gloc(ialph(i,1),icg)*dalpha(j,2,i)+ &
            gloc(ialph(i,1)+nside,icg)*domega(j,2,i)
           endif
           dphi(j,3,nres-1)+gloc(nres-3,icg)*dphi(j,2,nres) &
            +gloc(2*nres-6,icg)* &
            dtheta(j,2,nres-1)+gloc(2*nres-5,icg)*dtheta(j,1,nres)
-          if(itype(nres-2,1).ne.10) then
+          if((itype(nres-2,1).ne.10).and.(molnum(nres-1).ne.5)) then
               gcart(j,nres-2)=gcart(j,nres-2)+gloc(ialph(nres-2,1),icg)* &
              dalpha(j,2,nres-2)+gloc(ialph(nres-2,1)+nside,icg)* &
               domega(j,2,nres-2)
           endif
-          if(itype(nres-1,1).ne.10) then
+          if((itype(nres-1,1).ne.10).and.(molnum(nres-1).ne.5)) then
              gcart(j,nres-2)=gcart(j,nres-2)+gloc(ialph(nres-1,1),icg)* &
             dalpha(j,1,nres-1)+gloc(ialph(nres-1,1)+nside,icg)* &
              domega(j,1,nres-1)
          enddo
       endif 
 !  Settind dE/ddnres-1       
-!#define DEBUG
+#define DEBUG
 #ifdef DEBUG
           j=1
               write(iout,*)"in int to carta",nres-1,gcart(j,nres-1),gloc(nres-3,icg),dphi(j,3,nres), &
 
 #endif
 !#undef DEBUG
-        if(itype(nres-1,1).ne.10) then
+        if((itype(nres-1,1).ne.10).and.(molnum(nres-1).ne.5)) then
           gcart(j,nres-1)=gcart(j,nres-1)+gloc(ialph(nres-1,1),icg)* &
          dalpha(j,2,nres-1)+gloc(ialph(nres-1,1)+nside,icg)* &
           domega(j,2,nres-1)
               write(iout,*)"in int to cart",i, gxcart(j,i),gloc(ialph(i,1),icg),dalpha(j,3,i), &
               gloc(ialph(i,1)+nside,icg),domega(j,3,i)
 #endif
-!#undef DEBUG
+#undef DEBUG
             enddo
          endif     
        enddo                                                                                                                                                   
           endif
 !          if (itype(i+1,1).ne.10) then
           if ((itype(i+1,1).ne.10).and.&
-          (itype(i+1,molnum(i+1)).ne.ntyp1_molec(molnum(i+1)))) then
+          (itype(i+1,molnum(i+1)).ne.ntyp1_molec(molnum(i+1))).and.&
+           (molnum(i+1).ne.5)) then
           gxcart(j,i)=gxcart(j,i)-gloc_sc(3,i-1,icg) &
       *dtauangle(j,3,1,i+2)
            gcart(j,i)=gcart(j,i)+gloc_sc(3,i-1,icg) &
           endif
 !          if (itype(i-1,1).ne.10) then
           if ((itype(i-1,1).ne.10).and.&
-          (itype(i-1,molnum(i-1)).ne.ntyp1_molec(molnum(i-1)))) then
+          (itype(i-1,molnum(i-1)).ne.ntyp1_molec(molnum(i-1))).and.&
+         (molnum(i-1).ne.5)) then
            gcart(j,i)=gcart(j,i)+gloc_sc(1,i-2,icg)* &
            dtauangle(j,1,3,i+1)
           endif
 !          if (itype(i+1,1).ne.10) then
           if ((itype(i+1,1).ne.10).and.&
-          (itype(i+1,molnum(i+1)).ne.ntyp1_molec(molnum(i+1)))) then
+          (itype(i+1,molnum(i+1)).ne.ntyp1_molec(molnum(i+1)))&
+       .and. (molnum(i+1).ne.5)) then
            gcart(j,i)=gcart(j,i)+gloc_sc(2,i-1,icg)* &
            dtauangle(j,2,2,i+2)
 !          write(iout,*) "numer",i,gloc_sc(2,i-1,icg),