created util_jacobinosort

jacobi diag that does NOT reorder the eigenvalues
this was need to get decent grid invariance

src/nwdft/grid/grid_invariant.F

@@ -1,6 +1,6 @@
       Subroutine grid_invariant
 c
-c$Id: grid_invariant.F,v 1.18 2002-12-18 01:13:34 edo Exp $
+c$Id: grid_invariant.F,v 1.19 2003-01-18 03:58:34 edo Exp $
 c
       implicit none
 c      
@@ -117,7 +117,7 @@ c
       end
       Subroutine grid_invar_mat(coord, charge)
 c
-c$Id: grid_invariant.F,v 1.18 2002-12-18 01:13:34 edo Exp $
+c$Id: grid_invariant.F,v 1.19 2003-01-18 03:58:34 edo Exp $
 c
       implicit none
 c
@@ -135,7 +135,7 @@ c     allocate space for M and T (rotational invariance;
 c     see B.G. Johnson et al, CPL, 220, 377 (1994), and auxiliary 
 c     species R_a - T and (R_a -T)(R_a - T)^T
 c
-      double precision tvec(3), ramt(3), ramtm(3,3), eigs(3), work(12)
+      double precision tvec(3), ramt(3), ramtm(3,3), eigs(3)
       double precision ztot, lramt,coss,dett
 c
       integer ictr, icart, jcart
@@ -210,10 +210,6 @@ c
             enddo
          enddo
       enddo
-      if(ga_nodeid().eq.0.and.oprint) then
-         write(LuOut,*)' nuclear charge qrot: '
-         call output(qrot, 1, 3, 1, 3, 3, 3, 1)
-      endif
       valratio=abs(qrot(1,2)+qrot(1,3)+qrot(2,3))/
      /     (qrot(1,1)+qrot(2,2)+qrot(3,3))
       if(ga_nodeid().eq.0.and.oprint) then
@@ -230,9 +226,14 @@ c
 c     diagonalize M
 c
       call dcopy(9,qrot,1,qin,1)
-      call dcopy(12,0d0,0,work, 1)
       call dcopy(3,0d0,0,eigs, 1)
-      call util_jacobi(3, qrot,3, eigs)
+      call util_jacobinosort(3, qrot,3, eigs)
+      if(ga_nodeid().eq.0.and.oprint) then
+         write(LuOut,*)' nuclear charge qrot: '
+         call output(qin, 1, 3, 1, 3, 3, 3, 1)
+         write(LuOut,*)'eig: ',eigs
+      endif
+      
       ddiag=abs(eigs(1)-qin(1,1))+abs(eigs(2)-qin(2,2))+
      +     abs(eigs(3)-qin(3,3))
       if(abs(ddiag).lt.eps) then
@@ -248,6 +249,7 @@ c
          if(ga_nodeid().eq.0.and.oprint) then
             write(LuOut,*)'wqrot: '
             call output(qrot, 1, 3, 1, 3, 3, 3, 1)
+            write(LuOut,*)'eig: ',eigs
          
             write(LuOut,*) ' WARNING: nonsense rotational matrix'
             write(LuOut,*) ' cosine      ',coss

src/stepper/stpr_partit.F

@@ -4,7 +4,7 @@
      &       ATMASS, CMASS, TENIN, ENERGY,
      &       CONVGE, CONVGG, CONVGGM,
      &       NINTER, NMODE, RPATH, TRACK, BCKSTP)
-c $Id: stpr_partit.F,v 1.10 2000-10-30 22:40:00 d3e129 Exp $
+c $Id: stpr_partit.F,v 1.11 2003-01-18 03:58:34 edo Exp $
 c
       IMPLICIT  REAL*8(A-H,O-Z), INTEGER(I-N)
       LOGICAL TRACK, BCKSTP, RPATH, TROUB
@@ -396,7 +396,8 @@ C
 *debug:      write(6,*)' hessian ninternal print'
 *debug:      call stpr_prntpd(hess,(ninter*(ninter+1)/2),ninter,6)
 *debug:      call stpr_datestp(6,'partit: before jacobi')
-      CALL stpr_sjacobi(NAT3,NINTER,HESS,HIEIGS,SCRSQ,IERR)
+      CALL stpr_sjacobi(NAT3,NINTER,HESS,HIEIGS,SCRSQ,IERR,
+     .     .true.)
 *debug:      call stpr_datestp(6,'partit: after  jacobi')
 *debug:      write(6,*)'hess stpr_partit after jacobi',nat3tr
 *debug:      do i=1,nat3tr

src/util/stpr_sjacobi.F

@@ -1,5 +1,5 @@
-      subroutine stpr_sjacobi(nm,n,a,w,eivr,ierr)
-c $Id: stpr_sjacobi.F,v 1.1 2000-08-28 23:27:41 edo Exp $
+      subroutine stpr_sjacobi(nm,n,a,w,eivr,ierr,sorting)
+c $Id: stpr_sjacobi.F,v 1.2 2003-01-18 03:58:31 edo Exp $
       IMPLICIT  REAL*8(A-H,O-Z), INTEGER(I-N)
       dimension a(n*(n+1)/2),w(n),eivr(nm,n)
 c
@@ -10,6 +10,7 @@ c     the routine is much faster than the old hdiag routine written
 c     at m.i.t. that uses the jacobi method but not the variable
 c     threshold technique that is applied here
 c
+      logical sorting
 c-----parameters---------------
       sq2inv = 1.d0/sqrt(2.d0)
       t1     = 1.d-12
@@ -176,6 +177,7 @@ c
         idiag = idiag + i
         w(i) = a(idiag)
   170 continue
+      if(sorting) then
 c
 c     Arrange eigenvalues & vectors in ascending order.
 c
@@ -198,5 +200,6 @@ c
   200   continue
         goto 180
       endif
+      endif
       return
       end

src/util/util_jacobi.F

@@ -1,13 +1,33 @@
       subroutine util_jacobi(n, a, lda, e)
 *
-* $Id: util_jacobi.F,v 1.2 1999-07-27 21:00:19 d3e129 Exp $
+* $Id: util_jacobi.F,v 1.3 2003-01-18 03:58:31 edo Exp $
 *
       implicit none
+      integer n, lda
+      double precision a(lda,n)
+      double precision e(n)
+      call util_jacobi0(n, a, lda, e,.true.)
+      return
+      end
+      subroutine util_jacobinosort(n, a, lda, e)
+*
+* $Id: util_jacobi.F,v 1.3 2003-01-18 03:58:31 edo Exp $
+*
+      implicit none
+      integer n, lda
+      double precision a(lda,n)
+      double precision e(n)
+      call util_jacobi0(n, a, lda, e,.false.)
+      return
+      end
+      subroutine util_jacobi0(n, a, lda, e,sorting)
+      implicit none
 #include "util.fh"
 #include "mafdecls.fh"
       integer n, lda
       double precision a(lda,n)
       double precision e(n)
+      logical sorting
 c
 c     Diagonalize A using Jacobi overwriting a with 
 c     the eigenvectors and returning in e the eigenvalues
@@ -31,7 +51,7 @@ c     unit matrix into a
 c
       call dfill(lda*n, 0d0, a, 1)
       call dfill(n, 1d0, a, lda+1)
-      call stpr_sjacobi(lda, n, dbl_mb(k_v), e, a, ierr)
+      call stpr_sjacobi(lda, n, dbl_mb(k_v), e, a, ierr,sorting)
       if (ierr .ne. 0) call errquit('util_jacobi: ierr ', ierr)
 c
       if (.not. ma_pop_stack(l_v)) call errquit('jacobi:ma',0)