Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_cdncore
8 :
9 : CONTAINS
10 :
11 670 : SUBROUTINE cdncore(fmpi ,input,vacuum,noco,nococonv,sym,&
12 : stars,cell,sphhar,atoms,vTot,outDen,moments,results, EnergyDen)
13 :
14 : USE m_constants
15 : USE m_judft
16 : USE m_cdn_io
17 : USE m_cdnovlp
18 : USE m_cored
19 : USE m_coredr
20 : USE m_types
21 : USE m_xmlOutput
22 :
23 : #ifdef CPP_MPI
24 : USE m_mpi_bc_coreden
25 : #endif
26 :
27 : IMPLICIT NONE
28 :
29 :
30 : TYPE(t_mpi), INTENT(IN) :: fmpi
31 :
32 :
33 : TYPE(t_input), INTENT(IN) :: input
34 : TYPE(t_vacuum), INTENT(IN) :: vacuum
35 : TYPE(t_noco), INTENT(IN) :: noco
36 : TYPE(t_nococonv), INTENT(IN) :: nococonv
37 : TYPE(t_sym), INTENT(IN) :: sym
38 : TYPE(t_stars), INTENT(IN) :: stars
39 : TYPE(t_cell), INTENT(IN) :: cell
40 : TYPE(t_sphhar), INTENT(IN) :: sphhar
41 : TYPE(t_atoms), INTENT(IN) :: atoms
42 : TYPE(t_potden), INTENT(IN) :: vTot
43 : TYPE(t_potden), INTENT(INOUT) :: outDen
44 : TYPE(t_moments), INTENT(INOUT) :: moments
45 : TYPE(t_results), INTENT(INOUT) :: results
46 : TYPE(t_potden), INTENT(INOUT), OPTIONAL :: EnergyDen
47 :
48 : INTEGER :: jspin, n, iType, ierr
49 : REAL :: seig, rhoint, momint,rho11,rho22
50 : COMPLEX :: rho21
51 : LOGICAL, PARAMETER :: l_st=.FALSE.
52 : LOGICAL :: l_coreDenPresent
53 :
54 670 : REAL :: rh(atoms%msh,atoms%ntype,input%jspins)
55 670 : REAL :: qint(atoms%ntype,input%jspins)
56 670 : REAL :: tec(atoms%ntype,input%jspins)
57 670 : REAL :: rhTemp(atoms%msh,atoms%ntype,input%jspins)
58 670 : REAL ,ALLOCATABLE :: vr0(:,:,:)
59 :
60 670 : results%seigc = 0.0
61 670 : IF (fmpi%irank==0) THEN
62 866 : DO jspin = 1,input%jspins
63 1784 : DO n = 1,atoms%ntype
64 1449 : moments%svdn(n,jspin) = outDen%mt(1,0,n,jspin) / (sfp_const*atoms%rmsh(1,n)*atoms%rmsh(1,n))
65 : END DO
66 : END DO
67 : END IF
68 :
69 670 : l_CoreDenPresent = .FALSE.
70 670 : IF (input%kcrel==0) THEN
71 : ! Generate input file ecore for subsequent GW calculation
72 : ! 11.2.2004 Arno Schindlmayr
73 666 : IF ((input%gw==1 .or. input%gw==3).AND.(fmpi%irank==0)) THEN
74 0 : OPEN (15,file='ecore',status='unknown', action='write',form='unformatted')
75 : END IF
76 :
77 1553260 : rh = 0.0
78 3548 : tec = 0.0
79 3548 : qint = 0.0
80 666 : IF (input%frcor) THEN
81 0 : IF (fmpi%irank==0) THEN
82 0 : IF(isCoreDensityPresent()) THEN
83 0 : CALL readCoreDensity(input,atoms,rh,tec,qint)
84 0 : l_coreDenPresent = .TRUE.
85 : END IF
86 : END IF
87 : #ifdef CPP_MPI
88 0 : CALL MPI_BCAST(l_CoreDenPresent,1,MPI_LOGICAL,0,fmpi%mpi_comm,ierr)
89 0 : CALL mpi_bc_coreDen(fmpi,atoms,input,rh,tec,qint)
90 : #endif
91 : END IF
92 : END IF
93 :
94 1423500 : vr0=vtot%mt(:,0,:,:)
95 : IF (.false.) THEN !there should be a imput switch here!!
96 : vr0(:,:,1)=0.5*(vr0(:,:,1)+vr0(:,:,2))
97 : vr0(:,:,2)=vr0(:,:,1)
98 : ENDIF
99 :
100 : !add in core density
101 670 : IF (fmpi%irank==0) THEN
102 335 : IF (input%kcrel==0) THEN
103 860 : DO jspin = 1,input%jspins
104 527 : IF(PRESENT(EnergyDen)) THEN
105 0 : CALL cored(input,jspin,atoms,outDen%mt,sphhar,l_CoreDenPresent,vr0(:,:,jspin), qint,rh ,tec,seig, EnergyDen%mt)
106 : ELSE
107 527 : CALL cored(input,jspin,atoms,outDen%mt,sphhar,l_CoreDenPresent,vr0(:,:,jspin), qint,rh ,tec,seig)
108 : ENDIF
109 :
110 776297 : rhTemp(:,:,jspin) = rh(:,:,jspin)
111 860 : results%seigc = results%seigc + seig
112 : END DO
113 : ELSE
114 2 : IF(PRESENT(EnergyDen)) call juDFT_error("Energyden not implemented for relativistic core calculations")
115 2 : CALL coredr(input,atoms,seig, outDen%mt,sphhar,vr0,qint,rh)
116 2 : results%seigc = results%seigc + seig
117 : END IF
118 : END IF
119 1732 : DO jspin = 1,input%jspins
120 1062 : IF (fmpi%irank==0) THEN
121 1449 : DO n = 1,atoms%ntype
122 1449 : moments%stdn(n,jspin) = outDen%mt(1,0,n,jspin) / (sfp_const*atoms%rmsh(1,n)*atoms%rmsh(1,n))
123 : END DO
124 : END IF
125 1732 : IF ((noco%l_noco.and..not.input%ctail).AND.(fmpi%irank==0)) THEN
126 190 : IF (jspin==2) THEN
127 :
128 95 : IF(PRESENT(EnergyDen)) call juDFT_error("Energyden not implemented for noco")
129 : !pk non-collinear (start)
130 : !add the coretail-charge to the constant interstitial
131 : !charge (star 0), taking into account the direction of
132 : !magnetisation of this atom
133 245 : DO iType = 1,atoms%ntype
134 : !rho11=qint(iType,1)/(cell%volint * input%jspins)
135 : !rho22=qint(iType,2)/(cell%volint * input%jspins)
136 : !rho21=0.0
137 : !call nococonv%rotdenmat(itype,rho11,rho22,rho21,toGlobal=.true.) !Todo: shouldn'1 this be toGlobal=.true.???
138 : !outDen%pw(1,1) = outDen%pw(1,1)+rho11
139 : !outDen%pw(1,2) = outDen%pw(1,2)+rho22
140 : !outDen%pw(1,3) = outDen%pw(1,3)+rho21
141 150 : rho11=(qint(iType,1)/(cell%volint * input%jspins)+qint(iType,2)/(cell%volint * input%jspins))/2
142 150 : outDen%pw(1,1) = outDen%pw(1,1)+rho11
143 245 : outDen%pw(1,2) = outDen%pw(1,2)+rho11
144 : END DO
145 : !pk non-collinear (end)
146 : END IF
147 : END IF
148 : END DO
149 1732 : DO jspin = 1,input%jspins
150 1732 : IF (input%ctail) THEN
151 404 : IF (noco%l_noco.and.jspin==1) THEN
152 0 : rh(:,:,1)=(rh(:,:,1)+rh(:,:,2))/2.
153 0 : rh(:,:,2)=rh(:,:,1)
154 : END IF
155 404 : IF(PRESENT(EnergyDen)) call juDFT_error("Energyden not implemented for ctail")
156 : !+gu hope this works as well
157 : CALL cdnovlp(fmpi,sphhar,stars,atoms,sym,vacuum,&
158 : cell,input ,l_st,jspin,rh(:,:,jspin),&
159 404 : outDen%pw,outDen%mt,outDen%vac,vTot%pw_w,vTot%mt)
160 658 : ELSE IF ((fmpi%irank==0).AND.(.NOT.noco%l_noco)) THEN
161 344 : DO iType = 1,atoms%ntype
162 344 : outDen%pw(1,jspin) = outDen%pw(1,jspin) + qint(iType,jspin) / (input%jspins * cell%volint)
163 : END DO
164 : END IF
165 : END DO
166 :
167 670 : IF (input%kcrel==0) THEN
168 666 : IF (fmpi%irank==0) THEN
169 333 : CALL writeCoreDensity(input,atoms,rhTemp,tec,qint)
170 776963 : outDen%mtCore = rhTemp
171 2107 : outDen%tec = tec
172 2107 : outDen%qint = qint
173 : END IF
174 666 : IF ((input%gw==1 .or. input%gw==3).AND.(fmpi%irank==0)) CLOSE(15)
175 : END IF
176 :
177 670 : END SUBROUTINE cdncore
178 :
179 : END MODULE m_cdncore
|