rdmft.F90 5.82 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
!--------------------------------------------------------------------------------
! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

MODULE m_rdmft

CONTAINS

SUBROUTINE rdmft(eig_id,mpi,input,kpts,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
                 sphhar,sym,vTot,oneD,noco,results)

   USE m_types
   USE m_juDFT
16
   USE m_constants
17 18
   USE m_cdnval
   USE m_cdn_io
19 20 21 22 23
   USE m_cdncore
   USE m_qfix
#ifdef CPP_MPI
   USE m_mpi_bc_potden
#endif
24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45

   IMPLICIT NONE

   TYPE(t_mpi),           INTENT(IN)    :: mpi
   TYPE(t_input),         INTENT(IN)    :: input
   TYPE(t_kpts),          INTENT(IN)    :: kpts
   TYPE(t_banddos),       INTENT(IN)    :: banddos
   TYPE(t_cell),          INTENT(IN)    :: cell
   TYPE(t_atoms),         INTENT(IN)    :: atoms
   TYPE(t_enpara),        INTENT(IN)    :: enpara
   TYPE(t_stars),         INTENT(IN)    :: stars
   TYPE(t_vacuum),        INTENT(IN)    :: vacuum
   TYPE(t_dimension),     INTENT(IN)    :: dimension
   TYPE(t_sphhar),        INTENT(IN)    :: sphhar
   TYPE(t_sym),           INTENT(IN)    :: sym
   TYPE(t_potden),        INTENT(IN)    :: vTot
   TYPE(t_oneD),          INTENT(IN)    :: oneD
   TYPE(t_noco),          INTENT(IN)    :: noco
   TYPE(t_results),       INTENT(INOUT) :: results

   INTEGER,               INTENT(IN)    :: eig_id

46
   TYPE(t_cdnvalJob)                    :: cdnvalJob
47
   TYPE(t_potden)                       :: singleStateDen, overallDen
48 49 50
   TYPE(t_regionCharges)                :: regCharges
   TYPE(t_dos)                          :: dos
   TYPE(t_moments)                      :: moments
51 52
   INTEGER                              :: jspin, ikpt, iBand, jsp, jspmax
   REAL                                 :: fix
53
   LOGICAL                              :: converged
54
   CHARACTER(LEN=20)                    :: filename
55

56 57
   CALL juDFT_error('rdmft not yet implemented!', calledby = 'rdmft')

58 59 60
   converged = .FALSE.

   ! Calculate all single state densities
61 62 63 64 65 66 67 68 69
   CALL regCharges%init(input,atoms)
   CALL dos%init(input,atoms,dimension,kpts,vacuum)
   CALL moments%init(input,atoms)
   cdnvalJob%l_evp = .FALSE.
   cdnvalJob%nkptExtended = kpts%nkpt
   ALLOCATE(cdnvalJob%noccbd(kpts%nkpt))
   ALLOCATE(cdnvalJob%nStart(kpts%nkpt))
   ALLOCATE(cdnvalJob%nEnd(kpts%nkpt))
   ALLOCATE(cdnvalJob%weights(1,kpts%nkpt))
70
   DO jspin = 1, input%jspins
71
      jsp = MERGE(1,jspin,noco%l_noco)
72
      DO ikpt = 1, kpts%nkpt
73
         DO iBand = 1, results%neig(ikpt,jsp)
74
            ! Construct cdnvalJob object for this state
75 76 77 78 79 80 81 82 83 84 85 86 87 88 89
            ! (Reasonable parallelization is not yet done - should be placed over the loops enclosing this section)
            cdnvalJob%ikptStart = ikpt
            cdnvalJob%ikptIncrement = kpts%nkpt
            IF(mpi%irank.EQ.0) THEN
               cdnvalJob%noccbd = 1
               cdnvalJob%nStart = iBand
               cdnvalJob%nEnd = iBand
               cdnvalJob%weights = 0.0
               cdnvalJob%weights(1,ikpt) = 1.0
            ELSE
               cdnvalJob%noccbd = 0
               cdnvalJob%nStart = 1
               cdnvalJob%nEnd = 0
               cdnvalJob%weights = 0.0
            END IF
90 91

            ! Call cdnval to construct density
92 93 94 95 96
            WRITE(*,*) 'Note: some optional flags may have to be reset in rdmft before the cdnval call'
            WRITE(*,*) 'This is not yet implemented!'
            CALL singleStateDen%init(stars,atoms,sphhar,vacuum,input%jspins,noco%l_noco,POTDEN_TYPE_DEN)
            CALL cdnval(eig_id,mpi,kpts,jsp,noco,input,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
                        sphhar,sym,vTot,oneD,cdnvalJob,singleStateDen,regCharges,dos,results,moments)
97 98

            ! Store the density on disc (These are probably way too many densities to keep them in memory)
99 100
            filename = ''
            WRITE(filename,'(a,i1.1,a,i4.4,a,i5.5)') 'cdn-', jsp, '-', ikpt, '-', iBand
101 102 103 104
            IF (mpi%irank.EQ.0) THEN
               CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,CDN_ARCHIVE_TYPE_CDN_const,CDN_INPUT_DEN_const,&
                                 0,-1.0,0.0,.FALSE.,singleStateDen,TRIM(ADJUSTL(filename)))
            END IF
105 106 107 108 109 110 111
         END DO
      END DO
   END DO

   DO WHILE (.NOT.converged)

      ! Calculate overall density with current occupation numbers (don't forget core electron density)
112 113 114 115 116 117 118 119 120 121 122 123 124 125 126
      jspmax = input%jspins
      IF (noco%l_mperp) jspmax = 1
      DO jspin = 1,jspmax
         CALL cdnvalJob%init(mpi,input,kpts,banddos,noco,results,jspin)
         CALL cdnval(eig_id,mpi,kpts,jsp,noco,input,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
                     sphhar,sym,vTot,oneD,cdnvalJob,overallDen,regCharges,dos,results,moments)
      END DO
      CALL cdncore(results,mpi,dimension,oneD,input,vacuum,noco,sym,&
                   stars,cell,sphhar,atoms,vTot,overallDen,moments)
      IF (mpi%irank.EQ.0) THEN
         CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,overallDen,noco%l_noco,.TRUE.,.true.,fix)
      END IF
#ifdef CPP_MPI
      CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,overallDen)
#endif
127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145

      ! Calculate Coulomb potential for overall density (+including external potential)

      ! For all states calculate integral over Coulomb potential times single state density

      ! For all states calculate Integral over other potential contributions times single state density

      ! Construct exchange matrix in the basis of eigenstates

      ! Optimize occupation numbers

      ! Check convergence of occupation numbers and set "converged" flag

   END DO ! WHILE (.NOT.converged)

   ! Calculate final overall density

   ! Calculate total energy

146 147 148
END SUBROUTINE rdmft

END MODULE m_rdmft