!
! Copyright (C) 2001 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
!-----------------------------------------------------------------------
subroutine ns_adj
!-----------------------------------------------------------------------
! This routine tries to suggest to the code the right atomic orbital to 
! localize the charge on.
!
      use pwcom

      implicit none
      integer ::                                                        &
     &          na,nt,is,m1,m2,                                         &
     &          majs,mins,adjs,                                         &
     &          mol(5),nel,                                             &
     &          i,j,l,index(5)                                          

      real(kind=8) ::                                                   &
     &          totoc,                                                  &
     &          delta,lambda(5)                                         

      complex(kind=DP) ::                                               &
     &          vet(5,5), f(5,5),                                       &
     &          temp

      logical ::                                                        &
     &          nm,adjust

      do na = 1,nat
         nt = ityp(na)
         if (Hubbard_U(nt).ne.0.d0) then

            call tabd(nt,totoc)     
            nm = .false.

            if (starting_magnetization(nt).gt.0.d0) then
               majs = 1
               mins = 2
            else if (starting_magnetization(nt).lt.0.d0) then
               majs = 2
               mins = 1
            else
               nm = .true.
            end if

            if (.not.nm) then
               if (totoc.gt.5.d0) then
                  do m1 = 1,5
                     do m2 = 1,5
                        f(m1,m2) = nsnew(na,mins,m1,m2)
                     end do
                  end do
                  adjs = mins
               else if (totoc.lt.5.d0) then
                  do m1 = 1,5
                     do m2 = 1,5
                        f(m1,m2) = nsnew(na,majs,m1,m2)
                     end do
                  end do
                  adjs = majs
               end if
            else
               go to 100
            end if

            call cdiagh(5,f,5,lambda,vet)
            
            do i = 1,5
               mol(i) = 1
               index(i) = i
            end do

            do i = 1,4
               if (index(i).eq.i) then
                  do j = i+1,5
                     delta = dabs(lambda(i)-lambda(j))
                     if (delta.lt.5.d-4) then
                        mol(i) = mol(i) + 1
                        index(j) = i
                        mol(j) = 100
                     end if
                  end do
               end if
            end do
               
            if (totoc.ge.5.d0) then
               nel = nint(totoc-5.d0)
            else
               nel = nint(totoc)
            end if

            adjust = .false.
            do i = 5,1,-1
               if (mol(i).eq.nel) then
                  adjust = .true.
                  do j = 5,1,-1
                     if (index(j).eq.i) then
                        lambda(j) = 1.d0
                     else
                        lambda(j) = 0.d0
                     end if
                  end do
               end if
            end do

            if (adjust) then
               do i = 1,5
                  do j = i,5
                     temp = 0.d0
                     do l = 1,5
                        temp = temp + conjg(vet(i,l))*lambda(l)*vet(j,l)     
                     end do
                     nsnew(na,adjs,i,j) = dreal(temp)
                     nsnew(na,adjs,j,i) = nsnew(na,adjs,i,j)
                  end do
               end do
            end if

         end if
 100     continue
      end do          ! on na

      return
end subroutine ns_adj