cdngen.F90 9.29 KB
Newer Older
1 2 3 4 5
!--------------------------------------------------------------------------------
! Copyright (c) 2016 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.
!--------------------------------------------------------------------------------
6 7 8 9
MODULE m_cdngen
CONTAINS

SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
10
                  dimension,kpts,atoms,sphhar,stars,sym,&
11
                  enpara,cell,noco,vTot,results,oneD,coreSpecInput,&
12
                  archiveType, xcpot,outDen,EnergyDen)
13 14 15 16 17 18 19

   !*****************************************************
   !    Charge density generator
   !    calls cdnval to generate the valence charge and the
   !    core routines for the core contribution
   !*****************************************************

20
   USE m_types
21
   USE m_constants
22
   USE m_juDFT
23 24
   USE m_prpqfftmap
   USE m_cdnval
25
   USE m_plot
26 27 28 29
   USE m_cdn_io
   USE m_wrtdop
   USE m_cdntot
   USE m_qfix
30
   USE m_genNewNocoInp
31
   USE m_xmlOutput
32
   USE m_magMoms
33
   USE m_orbMagMoms
34
   USE m_resMoms
35
   USE m_cdncore
36 37
   USE m_doswrite
   USE m_Ekwritesl
38
   USE m_banddos_io
39
   USE m_metagga
40
   USE m_unfold_band_kpts
41
   USE m_denMultipoleExp
42
#ifdef CPP_MPI
43
   USE m_mpi_bc_potden
44
#endif
45

46
   IMPLICIT NONE
47

48 49 50 51
#ifdef CPP_MPI
   INCLUDE 'mpif.h'
#endif

52 53 54 55 56 57 58 59 60 61
   ! Type instance arguments
   TYPE(t_results),INTENT(INOUT)    :: results
   TYPE(t_mpi),INTENT(IN)           :: mpi
   TYPE(t_dimension),INTENT(IN)     :: dimension
   TYPE(t_oneD),INTENT(IN)          :: oneD
   TYPE(t_enpara),INTENT(INOUT)     :: enpara
   TYPE(t_banddos),INTENT(IN)       :: banddos
   TYPE(t_sliceplot),INTENT(IN)     :: sliceplot
   TYPE(t_input),INTENT(IN)         :: input
   TYPE(t_vacuum),INTENT(IN)        :: vacuum
62
   TYPE(t_noco),INTENT(INOUT)       :: noco
63 64 65 66 67 68 69 70
   TYPE(t_sym),INTENT(IN)           :: sym
   TYPE(t_stars),INTENT(IN)         :: stars
   TYPE(t_cell),INTENT(IN)          :: cell
   TYPE(t_kpts),INTENT(IN)          :: kpts
   TYPE(t_sphhar),INTENT(IN)        :: sphhar
   TYPE(t_atoms),INTENT(IN)         :: atoms
   TYPE(t_coreSpecInput),INTENT(IN) :: coreSpecInput
   TYPE(t_potden),INTENT(IN)        :: vTot
71 72
   CLASS(t_xcpot),INTENT(INOUT)     :: xcpot
   TYPE(t_potden),INTENT(INOUT)     :: outDen, EnergyDen
73 74

   !Scalar Arguments
75
   INTEGER, INTENT (IN)             :: eig_id, archiveType
76 77

   ! Local type instances
78
   TYPE(t_noco)          :: noco_new
Matthias Redies's avatar
Matthias Redies committed
79 80 81 82
   TYPE(t_regionCharges) :: regCharges, fake_regCharges
   TYPE(t_dos)           :: dos, fake_dos
   TYPE(t_moments)       :: moments, fake_moments
   TYPE(t_results)       :: fake_results
83 84
   TYPE(t_mcd)           :: mcd
   TYPE(t_slab)          :: slab
85
   TYPE(t_orbcomp)       :: orbcomp
86
   TYPE(t_cdnvalJob)     :: cdnvalJob
87
   TYPE(t_potden)        :: val_den, core_den
88

89

90
   !Local Scalars
91
   REAL                  :: fix, qtot, dummy,eFermiPrev
92
   INTEGER               :: jspin, jspmax, ierr
Matthias Redies's avatar
Matthias Redies committed
93
   INTEGER               :: dim_idx
94

95 96 97
#ifdef CPP_HDF
   INTEGER(HID_T)        :: banddosFile_id
#endif
98
   LOGICAL               :: l_error, perform_MetaGGA
99

100 101
   CALL regCharges%init(input,atoms)
   CALL dos%init(input,atoms,dimension,kpts,vacuum)
102
   CALL moments%init(mpi,input,sphhar,atoms)
103 104 105
   CALL mcd%init1(banddos,dimension,input,atoms,kpts)
   CALL slab%init(banddos,dimension,atoms,cell,input,kpts)
   CALL orbcomp%init(input,banddos,dimension,atoms,kpts)
106

107 108
   CALL outDen%init(stars,    atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_DEN)
   CALL EnergyDen%init(stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_EnergyDen)
109

Matthias Redies's avatar
Matthias Redies committed
110
   IF (mpi%irank == 0) CALL openXMLElementNoAttributes('valenceDensity')
111

112
   !In a non-collinear calcuation where the off-diagonal part of the
113 114
   !density matrix in the muffin-tins is calculated, the a- and
   !b-coef. for both spins are needed at once. Thus, cdnval is only
115
   !called once and both spin directions are calculated in a single run.
116
   results%force=0.0
117 118 119
   jspmax = input%jspins
   IF (noco%l_mperp) jspmax = 1
   DO jspin = 1,jspmax
120 121
      CALL cdnvalJob%init(mpi,input,kpts,noco,results,jspin)
      IF (sliceplot%slice) CALL cdnvalJob%select_slice(sliceplot,results,input,kpts,noco,jspin)
122
      CALL cdnval(eig_id,mpi,kpts,jspin,noco,input,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
123
                  sphhar,sym,vTot,oneD,cdnvalJob,outDen,regCharges,dos,results,moments,coreSpecInput,mcd,slab,orbcomp)
124
   END DO
125
   call val_den%copyPotDen(outDen)
126

127 128
   ! calculate kinetic energy density for MetaGGAs
   if(xcpot%exc_is_metagga()) then
Matthias Redies's avatar
Matthias Redies committed
129
      CALL calc_EnergyDen(eig_id, mpi, kpts, noco, input, banddos, cell, atoms, enpara, stars,&
130
                             vacuum, DIMENSION, sphhar, sym, vTot, oneD, results, EnergyDen)
131 132
   endif

Matthias Redies's avatar
Matthias Redies committed
133
   IF (mpi%irank == 0) THEN
134
      IF (banddos%dos.or.banddos%vacdos.or.input%cdinf) THEN
135 136 137 138 139
         IF (banddos%unfoldband) THEN
            eFermiPrev = 0.0
            CALL readPrevEFermi(eFermiPrev,l_error)
            CALL write_band_sc(kpts,results,eFermiPrev)
         END IF
140
#ifdef CPP_HDF
141
         CALL openBandDOSFile(banddosFile_id,input,atoms,cell,kpts,banddos)
142 143
         CALL writeBandDOSData(banddosFile_id,input,atoms,cell,kpts,results,banddos,dos,vacuum)
         CALL closeBandDOSFile(banddosFile_id)
144
#endif
145
         CALL timestart("cdngen: dos")
146
         CALL doswrite(eig_id,dimension,kpts,atoms,vacuum,input,banddos,sliceplot,noco,sym,cell,dos,mcd,results,slab,orbcomp,oneD)
Matthias Redies's avatar
Matthias Redies committed
147
         IF (banddos%dos.AND.(banddos%ndir == -3)) THEN
148
            CALL Ek_write_sl(eig_id,dimension,kpts,atoms,vacuum,input,jspmax,sym,cell,dos,slab,orbcomp,results)
149
         END IF
150
         CALL timestop("cdngen: dos")
151 152 153 154
      END IF
   END IF

   IF ((banddos%dos.OR.banddos%vacdos).AND.(banddos%ndir/=-2)) CALL juDFT_end("DOS OK",mpi%irank)
Matthias Redies's avatar
Matthias Redies committed
155
   IF (vacuum%nstm == 3) CALL juDFT_end("VACWAVE OK",mpi%irank)
156

157
   IF (mpi%irank.EQ.0) THEN
Matthias Redies's avatar
Matthias Redies committed
158
      CALL cdntot(stars,atoms,sym,vacuum,input,cell,oneD,outDen,.TRUE.,qtot,dummy)
159 160 161
      CALL closeXMLElement('valenceDensity')
   END IF ! mpi%irank = 0

162
   IF (sliceplot%slice) THEN
Matthias Redies's avatar
Matthias Redies committed
163
      IF (mpi%irank == 0) THEN
164
         CALL writeDensity(stars,noco,vacuum,atoms,cell,sphhar,input,sym,oneD,CDN_ARCHIVE_TYPE_CDN_const,CDN_INPUT_DEN_const,&
165
                           0,-1.0,0.0,.FALSE.,outDen,'cdn_slice')
166
         IF (sliceplot%iplot.EQ.1) CALL makeplots(stars, atoms, sphhar, vacuum, input, oneD, sym, cell, noco, outDen, 1, sliceplot) 
167
      END IF
168
      CALL juDFT_end("slice OK",mpi%irank)
169
   END IF
170

Uliana Alekseeva's avatar
Uliana Alekseeva committed
171
   CALL timestart("cdngen: cdncore")
172 173 174 175 176 177 178
   if(xcpot%exc_is_MetaGGA()) then
      CALL cdncore(mpi,dimension,oneD,input,vacuum,noco,sym,&
                   stars,cell,sphhar,atoms,vTot,outDen,moments,results, EnergyDen)
   else
      CALL cdncore(mpi,dimension,oneD,input,vacuum,noco,sym,&
                   stars,cell,sphhar,atoms,vTot,outDen,moments,results)
   endif
179
   call core_den%subPotDen(outDen, val_den)
Uliana Alekseeva's avatar
Uliana Alekseeva committed
180
   CALL timestop("cdngen: cdncore")
181

182
   IF(.FALSE.) CALL denMultipoleExp(input, mpi, atoms, sphhar, stars, sym, cell, oneD, outDen) ! There should be a switch in the inp file for this
183
   IF(mpi%irank.EQ.0) THEN
184
      IF(.FALSE.) CALL resMoms(input,atoms,sphhar,noco,outDen,moments%rhoLRes) ! There should be a switch in the inp file for this
185 186
   END IF

187 188
   CALL enpara%calcOutParams(input,atoms,vacuum,regCharges)

Matthias Redies's avatar
Matthias Redies committed
189
   IF (mpi%irank == 0) THEN
190
      CALL openXMLElementNoAttributes('allElectronCharges')
191
      CALL qfix(mpi,stars,atoms,sym,vacuum,sphhar,input,cell,oneD,outDen,noco%l_noco,.TRUE.,.true.,fix)
192 193
      CALL closeXMLElement('allElectronCharges')

Matthias Redies's avatar
Matthias Redies committed
194
      IF (input%jspins == 2) THEN
195 196 197
         noco_new = noco

         !Calculate and write out spin densities at the nucleus and magnetic moments in the spheres
198
         CALL magMoms(dimension,input,atoms,noco_new,vTot,moments)
199

200 201
         noco = noco_new

202 203 204
         !Generate and save the new nocoinp file if the directions of the local
         !moments are relaxed or a constraint B-field is calculated.
         IF (ANY(noco%l_relax(:atoms%ntype)).OR.noco%l_constr) THEN
205
            CALL genNewNocoInp(input,atoms,noco,noco_new)
206 207
         END IF

208
         IF (noco%l_soc) CALL orbMagMoms(input,atoms,noco,moments%clmom)
209
         
210
      END IF
Matthias Redies's avatar
Matthias Redies committed
211
   END IF ! mpi%irank == 0
212
   
213 214 215
   perform_MetaGGA = ALLOCATED(EnergyDen%mt) &
                   .AND. (xcpot%exc_is_MetaGGA() .or. xcpot%vx_is_MetaGGA())
   if(perform_MetaGGA) then
216
      call set_kinED(mpi, sphhar, atoms, sym, core_den, val_den, xcpot, &
217 218
                     input, noco, stars, cell, outDen, EnergyDen, vTot)
   endif
219
#ifdef CPP_MPI
220 221 222 223 224 225 226 227 228 229 230 231
   CALL MPI_BCAST(noco%l_ss,1,MPI_LOGICAL,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%l_mperp,1,MPI_LOGICAL,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%l_constr,1,MPI_LOGICAL,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%mix_b,1,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)

   CALL MPI_BCAST(noco%alphInit,atoms%ntype,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%alph,atoms%ntype,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%beta,atoms%ntype,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%b_con,atoms%ntype*2,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%l_relax,atoms%ntype,MPI_LOGICAL,0,mpi%mpi_comm,ierr)
   CALL MPI_BCAST(noco%qss,3,MPI_DOUBLE_PRECISION,0,mpi%mpi_comm,ierr)

232 233 234 235
   CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,outDen)
#endif

END SUBROUTINE cdngen
236

237
END MODULE m_cdngen