module math !----------------------------------------------------------------------------- use io_units, only:inp,iout implicit none !----------------------------------------------------------------------------- ! ! !----------------------------------------------------------------------------- contains !----------------------------------------------------------------------------- ! djacob.f !----------------------------------------------------------------------------- subroutine DJACOB(N,NMAX,MAXJAC,E,A,C,AII) ! IMPLICIT REAL*8 (A-H,O-Z) ! THE JACOBI DIAGONALIZATION PROCEDURE integer :: N,NMAX,MAXJAC ! COMMON INP,IOUT,IPN real(kind=8),DIMENSION(NMAX,N) :: A,C real(kind=8),DIMENSION(150) :: AJJ !el AII real(kind=8),DIMENSION(*) :: AII !el local variables integer :: l,i,j,k,IPIV,JPIV,IJAC,LIM,LT,IT,IN real(kind=8) :: e,SIN45,COS45,S45SQ,C45SQ real(kind=8) :: TEMPA,AIJMAX,TAIJ,TAII,TAJJ,TMT real(kind=8) :: ZAMMA,SINT,COST,SINSQ,COSSQ,AIIMIN real(kind=8) :: TAIK,TAJK,TCIK,TCJK,TEST,AMAX,GAMSQ,T SIN45 = .70710678 COS45 = .70710678 S45SQ = 0.50 C45SQ = 0.50 ! UNIT EIGENVECTOR MATRIX DO 70 I = 1,N DO 7 J = I,N A(J,I)=A(I,J) C(I,J) = 0.0 7 C(J,I) = 0.0 70 C(I,I) = 1.0 ! DETERMINATION OF SEARCH ARGUMENT, TEST AMAX = 0.0 DO 1 I = 1,N DO 1 J = 1,I TEMPA=DABS(A(I,J)) IF (AMAX-TEMPA) 2,1,1 2 AMAX = TEMPA 1 CONTINUE TEST = AMAX*E ! SEARCH FOR LARGEST OFF DIAGONAL ELEMENT DO 72 IJAC=1,MAXJAC AIJMAX = 0.0 DO 3 I = 2,N LIM = I-1 DO 3 J = 1,LIM TAIJ=DABS(A(I,J)) IF (AIJMAX-TAIJ) 4,3,3 4 AIJMAX = TAIJ IPIV = I JPIV = J 3 CONTINUE IF(AIJMAX-TEST)300,300,5 ! PARAMETERS FOR ROTATION 5 TAII = A(IPIV,IPIV) TAJJ = A(JPIV,JPIV) TAIJ = A(IPIV,JPIV) TMT = TAII-TAJJ IF(DABS(TMT/TAIJ)-1.0D-12) 60,60,6 60 IF(TAIJ) 10,10,11 6 ZAMMA=TAIJ/(2.0*TMT) 90 IF(DABS(ZAMMA)-0.38268)8,8,9 9 IF(ZAMMA)10,10,11 10 SINT = -SIN45 GO TO 12 11 SINT = SIN45 12 COST = COS45 SINSQ = S45SQ COSSQ = C45SQ GO TO 120 8 GAMSQ=ZAMMA*ZAMMA SINT=2.0*ZAMMA/(1.0+GAMSQ) COST = (1.0-GAMSQ)/(1.0+GAMSQ) SINSQ=SINT*SINT COSSQ=COST*COST ! ROTATION 120 DO 13 K = 1,N TAIK = A(IPIV,K) TAJK = A(JPIV,K) A(IPIV,K) = TAIK*COST+TAJK*SINT A(JPIV,K) = TAJK*COST-TAIK*SINT TCIK = C(IPIV,K) TCJK = C(JPIV,K) C(IPIV,K) = TCIK*COST+TCJK*SINT 13 C(JPIV,K) = TCJK*COST-TCIK*SINT A(IPIV,IPIV) = TAII*COSSQ+TAJJ*SINSQ+2.0*TAIJ*SINT*COST A(JPIV,JPIV) = TAII*SINSQ+TAJJ*COSSQ-2.0*TAIJ*SINT*COST A(IPIV,JPIV) = TAIJ*(COSSQ-SINSQ)-SINT*COST*TMT A(JPIV,IPIV) = A(IPIV,JPIV) DO 30 K = 1,N A(K,IPIV) = A(IPIV,K) 30 A(K,JPIV) = A(JPIV,K) 72 CONTINUE WRITE (IOUT,1000) AIJMAX 1000 FORMAT (/1X,'NONCONVERGENT JACOBI. LARGEST OFF-DIAGONAL ELE',& 'MENT = ',1PE14.7) ! ARRANGEMENT OF EIGENVALUES IN ASCENDING ORDER 300 DO 14 I=1,N 14 AJJ(I)=A(I,I) LT=N+1 DO 15 L=1,N LT=LT-1 AIIMIN=1.0E+30 DO 16 I=1,N IF(AJJ(I)-AIIMIN)17,16,16 17 AIIMIN=AJJ(I) IT=I 16 CONTINUE IN=L AII(IN)=AIIMIN AJJ(IT)=1.0E+30 DO 15 K=1,N 15 A(IN,K)=C(IT,K) DO 18 I=1,N IF(A(I,1))19,22,22 19 T=-1.0 GO TO 91 22 T=1.0 91 DO 18 J=1,N 18 C(J,I)=T*A(I,J) return end subroutine DJACOB !----------------------------------------------------------------------------- ! energy_p_new_barrier.F !----------------------------------------------------------------------------- subroutine vecpr(u,v,w) ! implicit real*8(a-h,o-z) real(kind=8),dimension(3) :: u,v,w w(1)=u(2)*v(3)-u(3)*v(2) w(2)=-u(1)*v(3)+u(3)*v(1) w(3)=u(1)*v(2)-u(2)*v(1) return end subroutine vecpr !----------------------------------------------------------------------------- real(kind=8) function scalar(u,v) !DIR$ INLINEALWAYS scalar !#ifndef OSF !DEC$ ATTRIBUTES FORCEINLINE::scalar !#endif ! implicit none real(kind=8),dimension(3) :: u,v !d double precision sc !d integer i !d sc=0.0d0 !d do i=1,3 !d sc=sc+u(i)*v(i) !d enddo !d scalar=sc scalar=u(1)*v(1)+u(2)*v(2)+u(3)*v(3) return end function scalar !----------------------------------------------------------------------------- ! sort.f !----------------------------------------------------------------------------- ! ! ! ################################################### ! ## COPYRIGHT (C) 1990 by Jay William Ponder ## ! ## All Rights Reserved ## ! ################################################### ! ! ######################################################### ! ## ## ! ## subroutine sort -- heapsort of an integer array ## ! ## ## ! ######################################################### ! ! ! "sort" takes an input list of integers and sorts it ! into ascending order using the Heapsort algorithm ! ! subroutine sort(n,list) ! implicit none integer :: i,j,k,n integer :: index,lists integer :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort !----------------------------------------------------------------------------- ! ! ! ############################################################## ! ## ## ! ## subroutine sort2 -- heapsort of real array with keys ## ! ## ## ! ############################################################## ! ! ! "sort2" takes an input list of reals and sorts it ! into ascending order using the Heapsort algorithm; ! it also returns a key into the original ordering ! ! subroutine sort2(n,list,key) ! implicit none integer :: i,j,k,n integer :: index,keys integer :: key(*) real(kind=8) :: lists real(kind=8) :: list(*) ! ! ! initialize index into the original ordering ! do i = 1, n key(i) = i end do ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) keys = key(k) else lists = list(index) keys = key(index) list(index) = list(1) key(index) = key(1) index = index - 1 if (index .le. 1) then list(1) = lists key(1) = keys return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) key(i) = key(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists key(i) = keys end do return end subroutine sort2 !----------------------------------------------------------------------------- ! ! ! ################################################################# ! ## ## ! ## subroutine sort3 -- heapsort of integer array with keys ## ! ## ## ! ################################################################# ! ! ! "sort3" takes an input list of integers and sorts it ! into ascending order using the Heapsort algorithm; ! it also returns a key into the original ordering ! ! subroutine sort3(n,list,key) ! implicit none integer :: i,j,k,n integer :: index integer :: lists integer :: keys integer :: list(*) integer :: key(*) ! ! ! initialize index into the original ordering ! do i = 1, n key(i) = i end do ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) keys = key(k) else lists = list(index) keys = key(index) list(index) = list(1) key(index) = key(1) index = index - 1 if (index .le. 1) then list(1) = lists key(1) = keys return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) key(i) = key(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists key(i) = keys end do return end subroutine sort3 !----------------------------------------------------------------------------- ! ! ! ################################################################# ! ## ## ! ## subroutine sort4 -- heapsort of integer absolute values ## ! ## ## ! ################################################################# ! ! ! "sort4" takes an input list of integers and sorts it into ! ascending absolute value using the Heapsort algorithm ! ! subroutine sort4(n,list) ! implicit none integer :: i,j,k,n integer :: index integer :: lists integer :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (abs(list(j)) .lt. abs(list(j+1))) j = j + 1 end if if (abs(lists) .lt. abs(list(j))) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort4 !----------------------------------------------------------------------------- ! ! ! ################################################################ ! ## ## ! ## subroutine sort5 -- heapsort of integer array modulo m ## ! ## ## ! ################################################################ ! ! ! "sort5" takes an input list of integers and sorts it ! into ascending order based on each value modulo "m" ! ! subroutine sort5(n,list,m) ! implicit none integer :: i,j,k,m,n integer :: index,smod integer :: jmod,j1mod integer :: lists integer :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then jmod = mod(list(j),m) j1mod = mod(list(j+1),m) if (jmod .lt. j1mod) then j = j + 1 else if (jmod.eq.j1mod .and. list(j).lt.list(j+1)) then j = j + 1 end if end if smod = mod(lists,m) jmod = mod(list(j),m) if (smod .lt. jmod) then list(i) = list(j) i = j j = j + j else if (smod.eq.jmod .and. lists.lt.list(j)) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort5 !----------------------------------------------------------------------------- ! ! ! ############################################################# ! ## ## ! ## subroutine sort6 -- heapsort of a text string array ## ! ## ## ! ############################################################# ! ! ! "sort6" takes an input list of character strings and sorts ! it into alphabetical order using the Heapsort algorithm ! ! subroutine sort6(n,list) ! implicit none integer :: i,j,k,n integer :: index character(len=256) :: lists character*(*) :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort6 !----------------------------------------------------------------------------- ! ! ! ################################################################ ! ## ## ! ## subroutine sort7 -- heapsort of text strings with keys ## ! ## ## ! ################################################################ ! ! ! "sort7" takes an input list of character strings and sorts it ! into alphabetical order using the Heapsort algorithm; it also ! returns a key into the original ordering ! ! subroutine sort7(n,list,key) ! implicit none integer :: i,j,k,n integer :: index integer :: keys integer :: key(*) character(len=256) :: lists character*(*) :: list(*) ! ! ! initialize index into the original ordering ! do i = 1, n key(i) = i end do ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) keys = key(k) else lists = list(index) keys = key(index) list(index) = list(1) key(index) = key(1) index = index - 1 if (index .le. 1) then list(1) = lists key(1) = keys return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) key(i) = key(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists key(i) = keys end do return end subroutine sort7 !----------------------------------------------------------------------------- ! ! ! ######################################################### ! ## ## ! ## subroutine sort8 -- heapsort to unique integers ## ! ## ## ! ######################################################### ! ! ! "sort8" takes an input list of integers and sorts it into ! ascending order using the Heapsort algorithm, duplicate ! values are removed from the final sorted list ! ! subroutine sort8(n,list) ! implicit none integer :: i,j,k,n integer :: index integer :: lists integer :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists ! ! remove duplicate values from final list ! j = 1 do i = 2, n if (list(i-1) .ne. list(i)) then j = j + 1 list(j) = list(i) end if end do if (j .lt. n) n = j return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort8 !----------------------------------------------------------------------------- ! ! ! ############################################################# ! ## ## ! ## subroutine sort9 -- heapsort to unique text strings ## ! ## ## ! ############################################################# ! ! ! "sort9" takes an input list of character strings and sorts ! it into alphabetical order using the Heapsort algorithm, ! duplicate values are removed from the final sorted list ! ! subroutine sort9(n,list) ! implicit none integer :: i,j,k,n integer :: index character(len=256) :: lists character*(*) :: list(*) ! ! ! perform the heapsort of the input list ! k = n/2 + 1 index = n do while (n .gt. 1) if (k .gt. 1) then k = k - 1 lists = list(k) else lists = list(index) list(index) = list(1) index = index - 1 if (index .le. 1) then list(1) = lists ! ! remove duplicate values from final list ! j = 1 do i = 2, n if (list(i-1) .ne. list(i)) then j = j + 1 list(j) = list(i) end if end do if (j .lt. n) n = j return end if end if i = k j = k + k do while (j .le. index) if (j .lt. index) then if (list(j) .lt. list(j+1)) j = j + 1 end if if (lists .lt. list(j)) then list(i) = list(j) i = j j = j + j else j = index + 1 end if end do list(i) = lists end do return end subroutine sort9 !----------------------------------------------------------------------------- ! pinorm.f !----------------------------------------------------------------------------- real(kind=8) function pinorm(x) ! implicit real*8 (a-h,o-z) ! use geometry_data, only: pi,dwapi ! this function takes an angle (in radians) and puts it in the range of ! -pi to +pi. ! integer :: n real(kind=8) :: x ! include 'COMMON.GEO' n = x / dwapi pinorm = x - n * dwapi if ( pinorm .gt. pi ) then pinorm = pinorm - dwapi else if ( pinorm .lt. - pi ) then pinorm = pinorm + dwapi end if return end function pinorm !----------------------------------------------------------------------------- ! minimize_p.F !----------------------------------------------------------------------------- subroutine xx2x(x,xx) ! implicit real*8 (a-h,o-z) use geometry_data use energy_data ! include 'DIMENSIONS' ! include 'COMMON.VAR' ! include 'COMMON.CHAIN' ! include 'COMMON.INTERACT' integer :: i,ij,ig,igall real(kind=8),dimension(6*nres) :: xx,x !(maxvar) (maxvar=6*maxres) do i=1,nvar x(i)=varall(i) enddo ig=0 igall=0 do i=4,nres igall=igall+1 if (mask_phi(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif enddo do i=3,nres igall=igall+1 if (mask_theta(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif enddo do ij=1,2 do i=2,nres-1 if (itype(i).ne.10) then igall=igall+1 if (mask_side(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif endif enddo enddo return end subroutine xx2x !----------------------------------------------------------------------------- !----------------------------------------------------------------------------- end module math