trial for gen_rand
[unres4.git] / source / unres / energy.F90
index 45dd795..9d59ebc 100644 (file)
           weights_(41)=wcatcat
           weights_(42)=wcatprot
           weights_(46)=wscbase
           weights_(41)=wcatcat
           weights_(42)=wcatprot
           weights_(46)=wscbase
-          weights_(47)=wscpho
-          weights_(48)=wpeppho
+          weights_(47)=wpepbase
+          weights_(48)=wscpho
+          weights_(49)=wpeppho
 !          wcatcat= weights(41)
 !          wcatprot=weights(42)
 
 !          wcatcat= weights(41)
 !          wcatprot=weights(42)
 
           wcatcat= weights(41)
           wcatprot=weights(42)
           wscbase=weights(46)
           wcatcat= weights(41)
           wcatprot=weights(42)
           wscbase=weights(46)
-          wscpho=weights(47)
-          wpeppho=weights(48)
+          wpepbase=weights(47)
+          wscpho=weights(48)
+          wpeppho=weights(49)
+!      welpsb=weights(28)*fact(1)
+!
+!      wcorr_nucl= weights(37)*fact(1)
+!     wcorr3_nucl=weights(38)*fact(2)
+!     wtor_nucl=  weights(35)*fact(1)
+!     wtor_d_nucl=weights(36)*fact(2)
+
         endif
         time_Bcast=time_Bcast+MPI_Wtime()-time00
         time_Bcastw=time_Bcastw+MPI_Wtime()-time00
         endif
         time_Bcast=time_Bcast+MPI_Wtime()-time00
         time_Bcastw=time_Bcastw+MPI_Wtime()-time00
       epeppho=0.0
       endif
 !      call ecatcat(ecationcation)
       epeppho=0.0
       endif
 !      call ecatcat(ecationcation)
-      print *,"after ebend", ebe_nucl
+!      print *,"after ebend", wtor_nucl 
 #ifdef TIMING
       time_enecalc=time_enecalc+MPI_Wtime()-time00
 #endif
 #ifdef TIMING
       time_enecalc=time_enecalc+MPI_Wtime()-time00
 #endif
       wtor=weights(13)*fact(1)
       wtor_d=weights(14)*fact(2)
       wsccor=weights(21)*fact(1)
       wtor=weights(13)*fact(1)
       wtor_d=weights(14)*fact(2)
       wsccor=weights(21)*fact(1)
-
+      welpsb=weights(28)*fact(1)
+      wcorr_nucl= weights(37)*fact(1)
+      wcorr3_nucl=weights(38)*fact(2)
+      wtor_nucl=  weights(35)*fact(1)
+      wtor_d_nucl=weights(36)*fact(2)
+      wpepbase=weights(47)*fact(1)
       return
       end subroutine rescale_weights
 !-----------------------------------------------------------------------------
       return
       end subroutine rescale_weights
 !-----------------------------------------------------------------------------
       integer :: i,iti1,iti,k,l
       real(kind=8) :: sin1,cos1,sin2,cos2,dwacos2,dwasin2,cost1,sint1,&
        sint1sq,sint1cub,sint1cost1,b1k,b2k,aux
       integer :: i,iti1,iti,k,l
       real(kind=8) :: sin1,cos1,sin2,cos2,dwacos2,dwasin2,cost1,sint1,&
        sint1sq,sint1cub,sint1cost1,b1k,b2k,aux
-       print *,"in set matrices"
+!       print *,"in set matrices"
 !
 ! Compute the virtual-bond-torsional-angle dependent quantities needed
 ! to calculate the el-loc multibody terms of various order.
 !
 ! Compute the virtual-bond-torsional-angle dependent quantities needed
 ! to calculate the el-loc multibody terms of various order.
 #endif
 #else
         if (i.gt. nnt+2 .and. i.lt.nct+2) then
 #endif
 #else
         if (i.gt. nnt+2 .and. i.lt.nct+2) then
-         write(iout,*) "i,",molnum(i)
-         print *, "i,",molnum(i),i,itype(i-2,1)
+!         write(iout,*) "i,",molnum(i)
+!         print *, "i,",molnum(i),i,itype(i-2,1)
         if (molnum(i).eq.1) then
           iti = itype2loc(itype(i-2,1))
         else
         if (molnum(i).eq.1) then
           iti = itype2loc(itype(i-2,1))
         else
         else
           iti1=nloctyp
         endif
         else
           iti1=nloctyp
         endif
-        print *,i,iti
+!        print *,i,iti
         b1(1,i-2)=b(3,iti)
         b1(2,i-2)=b(5,iti)
         b2(1,i-2)=b(2,iti)
         b1(1,i-2)=b(3,iti)
         b1(2,i-2)=b(5,iti)
         b2(1,i-2)=b(2,iti)
 #endif
 !        print *, "before set matrices"
         call set_matrices
 #endif
 !        print *, "before set matrices"
         call set_matrices
-        print *,"after set martices"
+!        print *,"after set martices"
 #ifdef TIMING
         time_mat=time_mat+MPI_Wtime()-time01
 #endif
 #ifdef TIMING
         time_mat=time_mat+MPI_Wtime()-time01
 #endif
       integer :: i,j
       
       if(nres.lt.100) then
       integer :: i,j
       
       if(nres.lt.100) then
-        maxconts=nres
+        maxconts=10*nres
       elseif(nres.lt.200) then
       elseif(nres.lt.200) then
-        maxconts=0.8*nres      ! Max. number of contacts per residue
+        maxconts=10*nres      ! Max. number of contacts per residue
       else
       else
-        maxconts=0.6*nres ! (maxconts=maxres/4)
+        maxconts=10*nres ! (maxconts=maxres/4)
       endif
       maxcont=12*nres      ! Max. number of SC contacts
       maxvar=6*nres      ! Max. number of variables
       endif
       maxcont=12*nres      ! Max. number of SC contacts
       maxvar=6*nres      ! Max. number of variables
       enddo
 !      IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and.
        IF ( j.gt.i+1 .and.&
       enddo
 !      IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and.
        IF ( j.gt.i+1 .and.&
-          num_conti.le.maxconts) THEN
+          num_conti.le.maxcont) THEN
 !C
 !C Calculate the contact function. The ith column of the array JCONT will 
 !C contain the numbers of atoms that make contacts with the atom I (of numbers
 !C greater than I). The arrays FACONT and GACONT will contain the values of
 !C the contact function and its derivative.
 !C
 !C Calculate the contact function. The ith column of the array JCONT will 
 !C contain the numbers of atoms that make contacts with the atom I (of numbers
 !C greater than I). The arrays FACONT and GACONT will contain the values of
 !C the contact function and its derivative.
-        r0ij=2.20D0*sigma(itypi,itypj)
+        r0ij=2.20D0*sigma_nucl(itypi,itypj)
 !c        write (2,*) "ij",i,j," rij",1.0d0/rij," r0ij",r0ij
         call gcont(rij,r0ij,1.0D0,0.2d0/r0ij,fcont,fprimcont)
 !c        write (2,*) "fcont",fcont
 !c        write (2,*) "ij",i,j," rij",1.0d0/rij," r0ij",r0ij
         call gcont(rij,r0ij,1.0D0,0.2d0/r0ij,fcont,fprimcont)
 !c        write (2,*) "fcont",fcont
 
           if (num_conti.gt.maxconts) then
             write (iout,*) 'WARNING - max. # of contacts exceeded;',&
 
           if (num_conti.gt.maxconts) then
             write (iout,*) 'WARNING - max. # of contacts exceeded;',&
-                          ' will skip next contacts for this conf.'
+                          ' will skip next contacts for this conf.',maxconts
           else
             jcont_hb(num_conti,i)=j
 !c            write (iout,*) "num_conti",num_conti,
           else
             jcont_hb(num_conti,i)=j
 !c            write (iout,*) "num_conti",num_conti,