Commit 1069634b authored by Daniel Wortmann's avatar Daniel Wortmann

Merge remote-tracking branch 'origin/develop' into release

parents c74e6d87 3b8dcf56
......@@ -2,8 +2,8 @@ init/compileinfo.h
*~
#*
build
build.*
*.o
*.mod
*.x
*.swp
cmake_minimum_required(VERSION 2.8)
project(FLEUR LANGUAGES NONE)
project(FLEUR LANGUAGES C Fortran)
#some variables might be set in the environment
set(FLEUR_LIBRARIES ${FLEUR_LIBRARIES} $ENV{FLEUR_LIBRARIES})
set(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} $ENV{CMAKE_Fortran_FLAGS}")
if (DEFINED ENV{FLEUR_NO_SERIAL})
set(FLEUR_USE_SERIAL false)
else()
set(FLEUR_USE_SERIAL true)
endif()
include("cmake/CompilerConfig.txt")
include("cmake/Architectures.txt")
include("cmake/ReportConfig.txt")
#Here the targets and the files are defined
include("cmake/Files_and_Targets.txt")
#install(TARGETS fleur inpgen DESTINATION bin)
......
......@@ -66,30 +66,20 @@ For the compilation of FLEUR you will need:
(without you can only exploit k-point parallelism)
* HDF5 for parallel IO (if you have lots of memory, you might not need to do IO :-)
You should probably create a directory structure like this:
FLEUR now comes with a configuration script. This skript will
create a build sub-directory in the main FLEUR directory (where this file resides) and
then use cmake to generate a makefile. You should use it like
somewhere/src -this is where the source code lives, i.e. the location of this README
somewhere/build -here you can build the code
./configure.sh CONFIGURATION
change to the build directory and type:
where CONFIGURATION refers to a predefined maschine setup. Call ./configure.sh without any
argument for more information.
cmake ../src
After the configure.sh script finishes, you should do
cd build
make
This might generate the FLEUR executables. It most probable will only work on systems we know as
only for those there will be specific configurations in the cmake directory.
If you get errors, you have to configure your build system yourself.
There are two options to do that.
- Either set the environment variable FC to either ifort,pgfortran or gfortran.
Then the configuration in the cmake directory will be used from cmake.[ifort|gfortran|pgfortran].config will be used that can again work for your system or can be adjusted.
- Or you run cmake with the option -DFleur_custom_toolchain and you adjust cmake.config to your needs.
You might find inspiration in the cmake.???.config files provided for several systems.
If your build environment is recognized correctly you can obtain the following executables:
inpgen - the input generator of FLEUR
fleur - general version of FLEUR
fleur_INVS - version of FLEUR suitable for systems with inversion symmetry leading to real matrices
fleur_SOC - version of FLEUR with spin-orbit coupling
All codes might also have a _MPI attached to indicate parallel versions.
\ No newline at end of file
only for those there will be specific configurations available.
......@@ -164,11 +164,6 @@ CONTAINS
INTEGER, ALLOCATABLE :: gvac1d(:),gvac2d(:) ,kveclo(:)
INTEGER, ALLOCATABLE :: jsym(:),ksym(:)
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
COMPLEX, ALLOCATABLE :: z(:,:)
#else
REAL, ALLOCATABLE :: z(:,:)
#endif
REAL, ALLOCATABLE :: aclo(:,:,:),acnmt(:,:,:,:,:)
REAL, ALLOCATABLE :: bclo(:,:,:),bcnmt(:,:,:,:,:)
REAL, ALLOCATABLE :: cclo(:,:,:,:),ccnmt(:,:,:,:,:),we(:)
......@@ -195,7 +190,11 @@ CONTAINS
TYPE (t_orblo),ALLOCATABLE :: orblo(:,:,:,:,:)
TYPE (t_mt21), ALLOCATABLE :: mt21(:,:)
TYPE (t_lo21), ALLOCATABLE :: lo21(:,:)
TYPE (t_usdus):: usdus
TYPE (t_usdus) :: usdus
TYPE (t_zMat) :: zMat
LOGICAL :: l_real
l_real=sym%invs.AND.(.NOT.noco%l_soc).AND.(.NOT.noco%l_noco)
! ..
! ..
llpd=(atoms%lmaxd*(atoms%lmaxd+3))/2
......@@ -521,14 +520,34 @@ CONTAINS
n_end = noccbd
END IF
END IF
IF (.NOT.ALLOCATED(z)) ALLOCATE (z(dimension%nbasfcn,dimension%neigd))
z = 0
zMat%l_real = l_real
IF (l_real) THEN
IF (.NOT.ALLOCATED(zMat%z_r)) THEN
ALLOCATE (zMat%z_r(dimension%nbasfcn,dimension%neigd))
zMat%nbasfcn = dimension%nbasfcn
zMat%nbands = dimension%neigd
END IF
zMat%z_r = 0
CALL cdn_read(&
eig_id,dimension%nvd,dimension%jspd,mpi%irank,mpi%isize,&
ikpt,jspin,dimension%nbasfcn,noco%l_ss,noco%l_noco,&
noccbd,n_start,n_end,&
lapw%nmat,lapw%nv,ello,evdu,epar,kveclo,&
lapw%k1,lapw%k2,lapw%k3,bkpt,wk,nbands,eig,zMat%z_r)
ELSE
IF (.NOT.ALLOCATED(zMat%z_c)) THEN
ALLOCATE (zMat%z_c(dimension%nbasfcn,dimension%neigd))
zMat%nbasfcn = dimension%nbasfcn
zMat%nbands = dimension%neigd
END IF
zMat%z_c = 0
CALL cdn_read(&
eig_id,dimension%nvd,dimension%jspd,mpi%irank,mpi%isize,&
ikpt,jspin,dimension%nbasfcn,noco%l_ss,noco%l_noco,&
noccbd,n_start,n_end,&
lapw%nmat,lapw%nv,ello,evdu,epar,kveclo,&
lapw%k1,lapw%k2,lapw%k3,bkpt,wk,nbands,eig,z)
lapw%k1,lapw%k2,lapw%k3,bkpt,wk,nbands,eig,zMat%z_c)
endif
!IF (l_evp.AND.(isize.GT.1)) THEN
! eig(1:noccbd) = eig(n_start:n_end)
!ENDIF
......@@ -564,7 +583,11 @@ CONTAINS
nslibd = nslibd + 1
eig(nslibd) = eig(i)
we(nslibd) = we(i)
z(:,nslibd) = z(:,i)
if (l_real) THEN
zMat%z_r(:,nslibd) = zMat%z_r(:,i)
else
zMat%z_c(:,nslibd) = zMat%z_c(:,i)
endif
END IF
END DO
IF (mpi%irank==0) WRITE (16,'(a,i3)') ' eigenvalues in sliceplot%slice:',nslibd
......@@ -576,14 +599,22 @@ CONTAINS
nslibd = nslibd + 1
eig(nslibd) = eig(sliceplot%nnne)
we(nslibd) = we(sliceplot%nnne)
z(:,nslibd) = z(:,sliceplot%nnne)
if (l_real) Then
zMat%z_r(:,nslibd) = zMat%z_r(:,sliceplot%nnne)
else
zMat%z_c(:,nslibd) = zMat%z_c(:,sliceplot%nnne)
endif
ELSE
DO i = 1,nbands
IF (eig(i).GE.sliceplot%e1s .AND. eig(i).LE.sliceplot%e2s) THEN
nslibd = nslibd + 1
eig(nslibd) = eig(i)
we(nslibd) = we(i)
z(:,nslibd) = z(:,i)
if (l_real) THEN
zMat%z_r(:,nslibd) = zMat%z_r(:,i)
else
zMat%z_c(:,nslibd) = zMat%z_c(:,i)
endif
END IF
END DO
IF (mpi%irank==0) WRITE (16,FMT='(a,i3)')' eigenvalues in sliceplot%slice:',nslibd
......@@ -608,13 +639,8 @@ CONTAINS
! ----> valence density in the interstitial region
IF (.NOT.((jspin.EQ.2) .AND. noco%l_noco)) THEN
CALL timestart("cdnval: pwden")
CALL pwden(&
stars,kpts,banddos,oneD,&
input,mpi,noco,cell,atoms,sym,ikpt,&
jspin,lapw,noccbd,&
igq_fft,we,z,&
eig,bkpt,&
qpw,cdom,qis,results%force,f_b8)
CALL pwden(stars,kpts,banddos,oneD, input,mpi,noco,cell,atoms,sym,ikpt,&
jspin,lapw,noccbd,igq_fft,we, eig,bkpt,qpw,cdom,qis,results%force,f_b8,zMat,l_real)
CALL timestop("cdnval: pwden")
END IF
!+new
......@@ -623,13 +649,9 @@ CONTAINS
!
IF (banddos%dos.AND.(banddos%ndir.EQ.-3)) THEN
IF (.NOT.((jspin.EQ.2) .AND. noco%l_noco)) THEN
CALL q_int_sl(&
jspin,stars,atoms,sym,&
volsl,volintsl,&
cell,&
z,noccbd,lapw,&
nsl,zsl,nmtsl,oneD,&
qintsl(:,:))
CALL q_int_sl(jspin,stars,atoms,sym, volsl,volintsl,&
cell,noccbd,lapw, nsl,zsl,nmtsl,oneD, qintsl(:,:),zMat,l_real)
!
END IF
END IF
......@@ -638,16 +660,9 @@ CONTAINS
IF (input%film) THEN
IF (.NOT.((jspin.EQ.2) .AND. noco%l_noco)) THEN
CALL timestart("cdnval: vacden")
CALL vacden(&
vacuum,dimension,stars,oneD,&
kpts,input,&
cell,atoms,noco,banddos,&
gvac1d,gvac2d,&
we,ikpt,jspin,vz,vz0,&
noccbd,z,bkpt,lapw,&
evac,eig,&
rhtxy,rht,qvac,qvlay,&
qstars,cdomvz,cdomvxy)
CALL vacden(vacuum,dimension,stars,oneD, kpts,input, cell,atoms,noco,banddos,&
gvac1d,gvac2d, we,ikpt,jspin,vz,vz0, noccbd,bkpt,lapw, evac,eig,&
rhtxy,rht,qvac,qvlay, qstars,cdomvz,cdomvxy,zMat,l_real)
CALL timestop("cdnval: vacden")
END IF
!---> perform Brillouin zone integration and summation over the
......@@ -684,16 +699,15 @@ CONTAINS
aveccof(3,noccbd,0:atoms%lmaxd*(atoms%lmaxd+2),atoms%natd),&
bveccof(3,noccbd,0:atoms%lmaxd*(atoms%lmaxd+2),atoms%natd),&
cveccof(3,-atoms%llod:atoms%llod,noccbd,atoms%nlod,atoms%natd) )
CALL to_pulay(atoms,noccbd,sym, lapw, noco,cell,bkpt, z,noccbd,eig,usdus,&
CALL to_pulay(input,atoms,noccbd,sym, lapw, noco,cell,bkpt,noccbd,eig,usdus,&
kveclo,ispin,oneD, acof(:,0:,:,ispin),bcof(:,0:,:,ispin),&
e1cof,e2cof,aveccof,bveccof, ccof(-atoms%llod,1,1,1,ispin),acoflo,bcoflo,cveccof)
e1cof,e2cof,aveccof,bveccof, ccof(-atoms%llod,1,1,1,ispin),acoflo,bcoflo,cveccof,zMat,l_real)
CALL timestop("cdnval: to_pulay")
ELSE
CALL timestart("cdnval: abcof")
CALL abcof(atoms,noccbd,sym, cell, bkpt,lapw,noccbd,z, usdus, noco,ispin,kveclo,oneD,&
acof(:,0:,:,ispin),bcof(:,0:,:,ispin),ccof(-atoms%llod:,:,:,:,ispin))
CALL abcof(input,atoms,noccbd,sym, cell, bkpt,lapw,noccbd,usdus, noco,ispin,kveclo,oneD,&
acof(:,0:,:,ispin),bcof(:,0:,:,ispin),ccof(-atoms%llod:,:,:,:,ispin),zMat,l_real)
CALL timestop("cdnval: abcof")
END IF
......@@ -779,13 +793,12 @@ CONTAINS
IF (input%l_f) THEN
CALL timestart("cdnval: force_a12/21")
#ifndef CPP_APW
IF (.not.input%l_useapw) THEN
CALL force_a12(atoms,noccbd,sym, dimension,cell,oneD,&
we,ispin,noccbd,usdus,acof(:,0:,:,ispin),&
bcof(:,0:,:,ispin),e1cof,e2cof, acoflo,bcoflo, results,f_a12)
#endif
CALL force_a21(atoms,dimension,noccbd,sym,&
ENDIF
CALL force_a21(input,atoms,dimension,noccbd,sym,&
oneD,cell,we,ispin,epar(0:,:,ispin),noccbd,eig,usdus,acof(:,0:,:,ispin),&
bcof(:,0:,:,ispin),ccof(-atoms%llod:,:,:,:,ispin), aveccof,bveccof,cveccof,&
results,f_a21,f_b4)
......@@ -823,7 +836,7 @@ CONTAINS
cartk=matmul(bkpt,cell%bmat)
IF (banddos%ndir.GT.0) THEN
CALL sympsi(bkpt,lapw%nv(jspin),lapw%k1(:,jspin),lapw%k2(:,jspin),&
lapw%k3(:,jspin),sym,dimension,nbands,cell, z,eig,noco, ksym,jsym)
lapw%k3(:,jspin),sym,dimension,nbands,cell,eig,noco, ksym,jsym,zMat,l_real)
END IF
!
!--dw now write k-point data to tmp_dos
......@@ -837,7 +850,11 @@ CONTAINS
END IF
!---> end of loop over PE's
DEALLOCATE (z)
IF (l_real) THEN
DEALLOCATE (zMat%z_r)
ELSE
DEALLOCATE (zMat%z_c)
END IF
END IF ! --> end "IF ((mod(i_rec-1,mpi%isize).EQ.mpi%irank).OR.l_evp) THEN"
END DO !---> end of k-point loop
DEALLOCATE (we,f,g,usdus%us,usdus%dus,usdus%duds,usdus%uds,usdus%ddn)
......@@ -961,7 +978,7 @@ CONTAINS
!---> forces of equ. A8 of Yu et al.
IF ((input%l_f)) THEN
CALL timestart("cdnval: force_a8")
CALL force_a8(atoms,sphhar, ispin, vr,rho,&
CALL force_a8(input,atoms,sphhar, ispin, vr,rho,&
f_a12,f_a21,f_b4,f_b8,results%force)
CALL timestop("cdnval: force_a8")
END IF
......
......@@ -7,7 +7,7 @@
MODULE m_pwden
CONTAINS
SUBROUTINE pwden(stars,kpts,banddos,oneD, input,mpi,noco,cell,atoms,sym, &
ikpt,jspin,lapw,ne, igq_fft,we,z, eig,bkpt, qpw,cdom, qis,forces,f_b8)
ikpt,jspin,lapw,ne, igq_fft,we,eig,bkpt, qpw,cdom, qis,forces,f_b8,zMat,realdata)
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
! In this subroutine the star function expansion coefficients of
! the plane wave charge density is determined.
......@@ -89,17 +89,15 @@ CONTAINS
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_kpts),INTENT(IN) :: kpts
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_zMat),INTENT(IN) :: zMat
INTEGER, INTENT (IN) :: igq_fft(0:stars%kq1d*stars%kq2d*stars%kq3d-1)
REAL,INTENT(IN) :: we(:) !(nobd)
REAL,INTENT(IN) :: eig(:)!(dimension%neigd)
REAL,INTENT(IN) :: bkpt(3)
!-----> BASIS FUNCTION INFORMATION
INTEGER,INTENT(IN):: ne
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
COMPLEX, INTENT (IN) :: z(:,:) !(dimension%nbasfcn,dimension%neigd)
#else
REAL, INTENT (IN) :: z(:,:) !(dimension%nbasfcn,dimension%neigd)
#endif
LOGICAL,OPTIONAL,INTENT(IN)::realdata
!-----> CHARGE DENSITY INFORMATION
INTEGER,INTENT(IN) :: ikpt,jspin
COMPLEX,INTENT(INOUT) :: qpw(:,:) !(stars%n3d,dimension%jspd)
......@@ -131,13 +129,17 @@ CONTAINS
REAL, ALLOCATABLE :: ekin(:)
COMPLEX, ALLOCATABLE :: cwk(:),ecwk(:)
!
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
LOGICAL l_real
REAL CPP_BLAS_sdot
EXTERNAL CPP_BLAS_sdot
#else
COMPLEX CPP_BLAS_cdotc
EXTERNAL CPP_BLAS_cdotc
#endif
IF (PRESENT(realdata)) THEN
l_real=realdata
ELSE
l_real=zMat%l_real
ENDIF
!-------> ABBREVIATIONS
!
......@@ -189,21 +191,21 @@ CONTAINS
psi2i(0:stars%kq1d*stars%kq2d*stars%kq3d-1),&
rhomat(0:stars%kq1d*stars%kq2d*stars%kq3d-1,4) )
ELSE
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
IF (l_real) THEN
ALLOCATE ( psir(-stars%kq1d*stars%kq2d:2*stars%kq1d*stars%kq2d*(stars%kq3d+1)-1),&
psii(1),&
rhon(-stars%kq1d*stars%kq2d:stars%kq1d*stars%kq2d*(stars%kq3d+1)-1) )
IF (input%l_f) ALLOCATE ( kpsii(1),&
kpsir(-stars%kq1d*stars%kq2d:2*stars%kq1d*stars%kq2d*(stars%kq3d+1)-1),&
ekin(-stars%kq1d*stars%kq2d:2*stars%kq1d*stars%kq2d*(stars%kq3d+1)-1))
#else
ELSE
ALLOCATE ( psir(0:stars%kq1d*stars%kq2d*stars%kq3d-1),&
psii(0:stars%kq1d*stars%kq2d*stars%kq3d-1),&
rhon(0:stars%kq1d*stars%kq2d*stars%kq3d-1) )
IF (input%l_f) ALLOCATE ( kpsir(0:stars%kq1d*stars%kq2d*stars%kq3d-1),&
kpsii(0:stars%kq1d*stars%kq2d*stars%kq3d-1),&
ekin(0:stars%kq1d*stars%kq2d*stars%kq3d-1) )
#endif
ENDIF
ENDIF
!
!=======> CALCULATE CHARGE DENSITY USING FFT
......@@ -230,26 +232,24 @@ CONTAINS
IF (noco%l_noco) THEN
q0_11 = zero
q0_22 = zero
IF (.NOT.l_real ) THEN
DO nu = 1 , ne
#if ( !defined(CPP_INVERSION) )
q0_11 = q0_11 + we(nu) *&
CPP_BLAS_cdotc(lapw%nv(1),z(1,nu),1,z(1,nu),1)
q0_22 = q0_22 + we(nu) *&
CPP_BLAS_cdotc(lapw%nv(2),z(lapw%nv(1)+atoms%nlotot+1,nu),1,&
z(lapw%nv(1)+atoms%nlotot+1,nu),1)
#endif
q0_11 = q0_11 + we(nu) * CPP_BLAS_cdotc(lapw%nv(1),zMat%z_c(1,nu),1,zMat%z_c(1,nu),1)
q0_22 = q0_22 + we(nu) * CPP_BLAS_cdotc(lapw%nv(2),zMat%z_c(lapw%nv(1)+atoms%nlotot+1,nu),1, zMat%z_c(lapw%nv(1)+atoms%nlotot+1,nu),1)
ENDDO
ENDIF
q0_11 = q0_11/cell%omtil
q0_22 = q0_22/cell%omtil
ELSE
IF (l_real) THEN
DO nu = 1 , ne
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
q0=q0+we(nu)*CPP_BLAS_sdot(lapw%nv(jspin),z(1,nu),1,z(1,nu),1)
#else
q0=q0+we(nu)&
*REAL(CPP_BLAS_cdotc(lapw%nv(jspin),z(1,nu),1,z(1,nu),1))
#endif
q0=q0+we(nu)*CPP_BLAS_sdot(lapw%nv(jspin),zMat%z_r(1,nu),1,zMat%z_r(1,nu),1)
ENDDO
ELSE
DO nu = 1 , ne
q0=q0+we(nu) *REAL(CPP_BLAS_cdotc(lapw%nv(jspin),zMat%z_c(1,nu),1,zMat%z_c(1,nu),1))
ENDDO
ENDIF
q0 = q0/cell%omtil
ENDIF
!
......@@ -315,44 +315,39 @@ CONTAINS
psi2r=0.0
psi2i=0.0
!------> map WF into FFTbox
#if ( !defined(CPP_INVERSION) )
! the compiler complains about the aimag if z is real, though these
! lines are never reached in a collinear calculation
IF (noco%l_ss) THEN
DO iv = 1 , lapw%nv(1)
psi1r( iv1d(iv,1) ) = REAL( z(iv,nu) )
psi1i( iv1d(iv,1) ) = AIMAG( z(iv,nu) )
psi1r( iv1d(iv,1) ) = REAL( zMat%z_c(iv,nu) )
psi1i( iv1d(iv,1) ) = AIMAG( zMat%z_c(iv,nu) )
ENDDO
DO iv = 1 , lapw%nv(2)
psi2r( iv1d(iv,2) ) = REAL(z(lapw%nv(1)+atoms%nlotot+iv,nu))
psi2i( iv1d(iv,2) ) = AIMAG(z(lapw%nv(1)+atoms%nlotot+iv,nu))
psi2r( iv1d(iv,2) ) = REAL(zMat%z_c(lapw%nv(1)+atoms%nlotot+iv,nu))
psi2i( iv1d(iv,2) ) = AIMAG(zMat%z_c(lapw%nv(1)+atoms%nlotot+iv,nu))
ENDDO
ELSE
DO iv = 1 , lapw%nv(jspin)
psi1r( iv1d(iv,jspin) ) = REAL( z(iv,nu) )
psi1i( iv1d(iv,jspin) ) = AIMAG( z(iv,nu) )
psi2r(iv1d(iv,jspin))=REAL( z(lapw%nv(1)+atoms%nlotot+iv,nu))
psi2i(iv1d(iv,jspin))=AIMAG(z(lapw%nv(1)+atoms%nlotot+iv,nu))
psi1r( iv1d(iv,jspin) ) = REAL( zMat%z_c(iv,nu) )
psi1i( iv1d(iv,jspin) ) = AIMAG( zMat%z_c(iv,nu) )
psi2r(iv1d(iv,jspin))=REAL( zMat%z_c(lapw%nv(1)+atoms%nlotot+iv,nu))
psi2i(iv1d(iv,jspin))=AIMAG(zMat%z_c(lapw%nv(1)+atoms%nlotot+iv,nu))
ENDDO
ENDIF
#endif
ELSE
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
psir=0.0
#else
psir=0.0
psii=0.0
#endif
!------> map WF into FFTbox
IF (l_real) THEN
DO iv = 1 , lapw%nv(jspin)
psir( iv1d(iv,jspin) ) = zMat%z_r(iv,nu)
ENDDO
ELSE
DO iv = 1 , lapw%nv(jspin)
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
psir( iv1d(iv,jspin) ) = z(iv,nu)
#else
psir( iv1d(iv,jspin) ) = REAL(z(iv,nu))
psii( iv1d(iv,jspin) ) = AIMAG(z(iv,nu))
#endif
psir( iv1d(iv,jspin) ) = REAL(zMat%z_c(iv,nu))
psii( iv1d(iv,jspin) ) = AIMAG(zMat%z_c(iv,nu))
ENDDO
ENDIF
ENDIF
!
!------> do (real) inverse FFT; notice that the array psir is filled from
! 0 to ifftq3-1, but starts at -ifftq2 to give work space for rfft
......@@ -369,7 +364,7 @@ CONTAINS
CALL cfft(psi2r,psi2i,ifftq3,stars%kq3_fft,ifftq3,isn)
ELSE
isn = 1
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
IF (l_real) THEN
CALL rfft(isn,stars%kq1_fft,stars%kq2_fft,stars%kq3_fft+1,stars%kq1_fft,stars%kq2_fft,stars%kq3_fft,&
nw1,nw2,nw3,wsave,psir(ifftq3d), psir(-ifftq2))
......@@ -378,8 +373,7 @@ CONTAINS
DO in=-1,stars%kq3_fft,2
DO im=0,ifftq2-1
ir = ifftq2 * in + im
ekin(ir) = ekin(ir) - wtf(nu) * eig(nu) *&
(psir(ir)**2 + psir(ir+ifftq2)**2)
ekin(ir) = ekin(ir) - wtf(nu) * eig(nu) * (psir(ir)**2 + psir(ir+ifftq2)**2)
ENDDO
ENDDO
......@@ -394,28 +388,26 @@ CONTAINS
DO i = 1,3
s = s + xk(i)*cell%bmat(i,j)
ENDDO
kpsir( iv1d(iv,jspin) ) = s * z(iv,nu)
kpsir( iv1d(iv,jspin) ) = s * zMat%z_r(iv,nu)
ENDDO
CALL rfft(isn,stars%kq1_fft,stars%kq2_fft,stars%kq3_fft+1,stars%kq1_fft,stars%kq2_fft,stars%kq3_fft,&
nw1,nw2,nw3,wsave,kpsir(ifftq3d), kpsir(-ifftq2))
DO in=-1,stars%kq3_fft,2
DO im=0,ifftq2-1
ir = ifftq2 * in + im
ekin(ir) = ekin(ir) + wtf(nu) * 0.5 * &
(kpsir(ir)**2 + kpsir(ir+ifftq2)**2)
ekin(ir) = ekin(ir) + wtf(nu) * 0.5 * (kpsir(ir)**2 + kpsir(ir+ifftq2)**2)
ENDDO
ENDDO
ENDDO
ENDIF
#else
ELSE
CALL cfft(psir,psii,ifftq3,stars%kq1_fft,ifftq1,isn)
CALL cfft(psir,psii,ifftq3,stars%kq2_fft,ifftq2,isn)
CALL cfft(psir,psii,ifftq3,stars%kq3_fft,ifftq3,isn)
! GM forces part
IF (input%l_f) THEN
DO ir = 0,ifftq3d-1
ekin(ir) = ekin(ir) - wtf(nu)*eig(nu)*&
(psir(ir)**2+psii(ir)**2)
ekin(ir) = ekin(ir) - wtf(nu)*eig(nu)* (psir(ir)**2+psii(ir)**2)
ENDDO
DO j = 1,3
......@@ -429,8 +421,8 @@ CONTAINS
DO i = 1,3
s = s + xk(i)*cell%bmat(i,j)
ENDDO
kpsir( iv1d(iv,jspin) ) = s * REAL(z(iv,nu))
kpsii( iv1d(iv,jspin) ) = s * AIMAG(z(iv,nu))
kpsir( iv1d(iv,jspin) ) = s * REAL(zMat%z_c(iv,nu))
kpsii( iv1d(iv,jspin) ) = s * AIMAG(zMat%z_c(iv,nu))
ENDDO
CALL cfft(kpsir,kpsii,ifftq3,stars%kq1_fft,ifftq1,isn)
......@@ -438,12 +430,11 @@ CONTAINS
CALL cfft(kpsir,kpsii,ifftq3,stars%kq3_fft,ifftq3,isn)
DO ir = 0,ifftq3d-1
ekin(ir) = ekin(ir) + wtf(nu) * 0.5 *&
(kpsir(ir)**2+kpsii(ir)**2)
ekin(ir) = ekin(ir) + wtf(nu) * 0.5 * (kpsir(ir)**2+kpsii(ir)**2)
ENDDO
ENDDO
ENDIF
#endif
ENDIF
ENDIF
!----> calculate rho(r) = sum w(k)*f(nu)*conjg(psi_nu,k(r))*psi_nu,k(r)
! k,nu
......@@ -453,14 +444,10 @@ CONTAINS
!---> in the non-collinear case rho becomes a hermitian 2x2
!---> matrix (rhomat).
DO ir = 0,ifftq3d-1
rhomat(ir,1) = rhomat(ir,1) &
+ wtf(nu)*( psi1r(ir)**2 + psi1i(ir)**2 )
rhomat(ir,2) = rhomat(ir,2) &
+ wtf(nu)*( psi2r(ir)**2 + psi2i(ir)**2 )
rhomat(ir,3) = rhomat(ir,3) + wtf(nu)*&
(psi2r(ir)*psi1r(ir)+psi2i(ir)*psi1i(ir))
rhomat(ir,4) = rhomat(ir,4) + wtf(nu)*&
(psi2r(ir)*psi1i(ir)-psi2i(ir)*psi1r(ir))
rhomat(ir,1) = rhomat(ir,1) + wtf(nu)*( psi1r(ir)**2 + psi1i(ir)**2 )
rhomat(ir,2) = rhomat(ir,2) + wtf(nu)*( psi2r(ir)**2 + psi2i(ir)**2 )
rhomat(ir,3) = rhomat(ir,3) + wtf(nu)* (psi2r(ir)*psi1r(ir)+psi2i(ir)*psi1i(ir))
rhomat(ir,4) = rhomat(ir,4) + wtf(nu)* (psi2r(ir)*psi1i(ir)-psi2i(ir)*psi1r(ir))
ENDDO
!---> in a non-collinear calculation the interstitial charge
!---> cannot be calculated by a simple substraction if the
......@@ -488,42 +475,32 @@ CONTAINS
CMPLX(psi1r(igq_fft(ik)),psi1i(igq_fft(ik)))
ENDDO
DO istr = 1,stars%ng3_fft
CALL pwint(&
stars,atoms,sym,&
oneD,cell,stars%kv3(1,istr),x)
qis(nu,ikpt,1) = qis(nu,ikpt,1)&
+ REAL(cwk(istr)*x)/cell%omtil/REAL(ifftq3)
CALL pwint(stars,atoms,sym, oneD,cell,stars%kv3(1,istr),x)
qis(nu,ikpt,1) = qis(nu,ikpt,1) + REAL(cwk(istr)*x)/cell%omtil/REAL(ifftq3)
ENDDO
cwk=0.0
DO ik = 0 , stars%kmxq_fft - 1
cwk(stars%igfft(ik,1))=cwk(stars%igfft(ik,1))+CONJG(stars%pgfft(ik))*&
CMPLX(psi2r(igq_fft(ik)),psi2i(igq_fft(ik)))
cwk(stars%igfft(ik,1))=cwk(stars%igfft(ik,1))+CONJG(stars%pgfft(ik))* CMPLX(psi2r(igq_fft(ik)),psi2i(igq_fft(ik)))
ENDDO
DO istr = 1,stars%ng3_fft
CALL pwint(&
stars,atoms,sym,&
oneD,cell,&
stars%kv3(1,istr),&
x)
qis(nu,ikpt,input%jspins) = qis(nu,ikpt,input%jspins)&
+ REAL(cwk(istr)*x)/cell%omtil/REAL(ifftq3)
CALL pwint(stars,atoms,sym, oneD,cell, stars%kv3(1,istr), x)
qis(nu,ikpt,input%jspins) = qis(nu,ikpt,input%jspins) + REAL(cwk(istr)*x)/cell%omtil/REAL(ifftq3)
ENDDO
ENDIF
ELSE
#if ( defined(CPP_INVERSION) && !defined(CPP_SOC) )
IF (l_real) THEN
DO in=-1,stars%kq3_fft,2
DO im=0,ifftq2-1
ir = ifftq2 * in + im
rhon(ir) = rhon(ir) + wtf(nu) * ( psir(ir)**2 +&
psir(ir+ifftq2)**2 )
rhon(ir) = rhon(ir) + wtf(nu) * ( psir(ir)**2 + psir(ir+ifftq2)**2 )
ENDDO
ENDDO
#else
ELSE
DO ir = 0