1!
2! Copyright (C) 2002-2005 FPMD-CPV groups
3! This file is distributed under the terms of the
4! GNU General Public License. See the file `License'
5! in the root directory of the present distribution,
6! or http://www.gnu.org/copyleft/gpl.txt .
7!
8
9
10MODULE kohn_sham_states
11
12
13   IMPLICIT NONE
14   SAVE
15
16   PRIVATE
17
18   ! ...   print KS states to file KS.indx_ksout if ksout true
19   LOGICAL :: tksout
20   CHARACTER(LEN=2 ), PARAMETER :: ks_file       = 'KS'
21
22   INTEGER, ALLOCATABLE :: indx_ksout(:,:)  ! (state inds, spin indxs)
23   INTEGER, ALLOCATABLE :: n_ksout(:)       ! (spin indxs)
24
25   PUBLIC :: ks_states_init, ks_states_closeup
26   PUBLIC :: n_ksout, indx_ksout, tksout, print_all_states
27
28!  ----------------------------------------------
29CONTAINS
30!  ----------------------------------------------
31
32
33   SUBROUTINE ks_states_init( nspin, nprnks, iprnks )
34
35      INTEGER, INTENT(IN) :: nspin, nprnks(:)
36      INTEGER, INTENT(IN) :: iprnks(:,:)
37
38      INTEGER :: i, ip, k, nstates
39
40      ! ...   Tell the code which Kohn-Sham state should be printed to file
41      !
42      IF( ALLOCATED( n_ksout    ) ) DEALLOCATE( n_ksout )
43      IF( ALLOCATED( indx_ksout ) ) DEALLOCATE( indx_ksout )
44      !
45      tksout = ANY( nprnks > 0 )
46      !
47      IF( tksout ) THEN
48         nstates = MAXVAL( nprnks )
49         ALLOCATE( n_ksout( nspin ) )
50         ALLOCATE( indx_ksout( nstates, nspin) )
51         n_ksout( 1:nspin ) = nprnks( 1:nspin )
52         DO i = 1, nspin
53           DO k = 1, nprnks( i )
54              indx_ksout( k, i ) = iprnks( k, i )
55           END DO
56         END DO
57      END IF
58
59      RETURN
60   END SUBROUTINE ks_states_init
61
62!  ----------------------------------------------
63
64   SUBROUTINE ks_states_closeup()
65      IF( ALLOCATED( indx_ksout ) ) DEALLOCATE( indx_ksout )
66      IF( ALLOCATED( n_ksout ) ) DEALLOCATE( n_ksout )
67      tksout = .FALSE.
68      RETURN
69   END SUBROUTINE ks_states_closeup
70
71!  ----------------------------------------------
72!  ----------------------------------------------
73
74      SUBROUTINE print_all_states( ctot, iupdwn_tot, nupdwn_tot )
75
76        USE kinds,            ONLY : DP
77        USE mp_bands,         ONLY : intra_bgrp_comm
78        USE io_global,        ONLY : ionode, stdout
79        USE electrons_base,   ONLY : nupdwn, iupdwn, nspin
80
81        IMPLICIT NONE
82
83        ! ...   declare subroutine arguments
84        COMPLEX(DP), INTENT(IN) :: ctot(:,:)
85        INTEGER,     INTENT(IN) :: iupdwn_tot(2)
86        INTEGER,     INTENT(IN) :: nupdwn_tot(2)
87
88        ! ...   declare other variables
89        INTEGER ::  i, iss, iks, itot
90
91        CHARACTER(LEN=256) :: file_name
92        CHARACTER(LEN=10), DIMENSION(2) :: spin_name
93        CHARACTER (LEN=6), EXTERNAL :: int_to_char
94
95        CALL errore('print_all_states','unused since a long time, please check',1)
96        IF( tksout ) THEN
97
98          IF (ionode) THEN
99            WRITE( stdout,*)
100            WRITE( stdout,'( "   Kohn Sham state")')
101            WRITE( stdout,'( "   ---------------")')
102          END IF
103
104          IF( nspin == 2 ) THEN
105            spin_name(1) = '_UP_'
106            spin_name(2) = '_DW_'
107          ELSE
108            spin_name(1) = '_'
109            spin_name(2) = '_'
110          END IF
111
112          DO iss = 1, nspin
113            IF( tksout ) THEN
114              DO i = 1, n_ksout(iss)
115                iks = indx_ksout(i, iss)
116                IF( ( iks > 0 ) .AND. ( iks <= nupdwn( iss ) ) ) THEN
117                  itot = iks + iupdwn_tot(iss) - 1
118                  file_name = TRIM( ks_file ) // &
119                            & trim(spin_name(iss)) // trim( int_to_char( iks ) )
120                  CALL print_ks_states( ctot( :, itot ), file_name )
121                END IF
122              END DO
123            END IF
124          END DO
125
126        END IF
127
128        RETURN
129        ! ...
130      END SUBROUTINE print_all_states
131
132
133!  ----------------------------------------------
134!  ----------------------------------------------
135
136      SUBROUTINE print_ks_states( c, file_name )
137
138        USE kinds
139        USE mp, ONLY: mp_sum
140        USE io_global, ONLY: ionode, ionode_id, stdout
141        USE fft_base, ONLY: dfftp, dffts
142        USE fft_interfaces, ONLY: invfft
143        USE fft_helper_subroutines, ONLY: c2psi_gamma
144        USE fft_rho, ONLY: rho_r2g
145        USE mp_bands, ONLY: intra_bgrp_comm, inter_bgrp_comm, my_bgrp_id,&
146             root_bgrp_id, root_bgrp
147        USE gvect, ONLY: ig_l2g, mill, ecutrho
148        USE io_base, ONLY: write_rhog
149
150        IMPLICIT NONE
151
152        COMPLEX(DP),      INTENT(IN) :: c(:)
153        CHARACTER(LEN=*), INTENT(IN) :: file_name
154        COMPLEX(DP), ALLOCATABLE :: psi(:), rhog(:,:)
155        REAL(DP), ALLOCATABLE :: rhor(:,:)
156        INTEGER   ::  i
157        ! FIXME: reciprocal lattice vectors
158        REAL(DP) :: bogus1(3), bogus2(3), bogus3(3)
159        REAL(DP) :: charge
160
161        ALLOCATE( rhor(dfftp%nnr,1), psi(dfftp%nnr) )
162
163        CALL c2psi_gamma( dffts, psi, c )
164        CALL invfft( 'Wave', psi, dffts )
165
166        ! FIXME: not sure things will work in presence of a double grid
167        DO i = 1, dfftp%nnr
168           rhor( i,1 ) = DBLE( psi( i ) )**2
169        END DO
170        charge = SUM( rhor )
171        CALL mp_sum( charge, intra_bgrp_comm )
172
173        DEALLOCATE( psi )
174        ALLOCATE( rhog(dfftp%ngm,1) )
175        CALL rho_r2g( dfftp, rhor, rhog )
176        DEALLOCATE ( rhor )
177
178        ! Only the first band group collects and writes
179        IF ( my_bgrp_id == root_bgrp_id ) CALL write_rhog &
180                ( file_name, root_bgrp, intra_bgrp_comm, &
181                bogus1, bogus2, bogus3, .true., &
182                mill, ig_l2g, rhog, ecutrho )
183
184        IF ( ionode ) THEN
185          WRITE( stdout,'(3X,A15," integrated charge : ",F14.5)')  &
186     &      TRIM(file_name), charge / DBLE(dfftp%nr1*dfftp%nr2*dfftp%nr3)
187        END IF
188
189        DEALLOCATE( rhog )
190        ! ...
191        RETURN
192        ! ...
193      END SUBROUTINE print_ks_states
194
195!  ----------------------------------------------
196!
197END MODULE kohn_sham_states
198