diff --git a/eigen/eigen.F90 b/eigen/eigen.F90 index 2fa9213f27a3ef87fa4305fa58749b6213c5fdef..0ec60beb204f62f10f5aab7c531b2ccb48fe8c5c 100644 --- a/eigen/eigen.F90 +++ b/eigen/eigen.F90 @@ -7,267 +7,266 @@ !The eigenvalue problem is constructed and solved in this routine. The following steps are performed: - !> 1. Preparation: generate energy parameters, open eig-file - !> 2. CALL to mt_setup() : this constructs the local Hamiltonian (i.e. the Hamiltonian in the \f$u,\dot u, u_{lo} \f$ basis) LDA+U is also added here - !> 3. within the (collinear)spin and k-point loop: CALL to eigen_hssetup() to generate the matrices, CALL to eigen_diag() to perform diagonalization - !> 4. writing (saving) of eigenvectors - !> - !> The matrices generated and diagonalized here are of type m_mat as defined in m_types_mat. - !>@author D. Wortmann - SUBROUTINE eigen(mpi,stars,sphhar,atoms,obsolete,xcpot,sym,kpts,DIMENSION,vacuum,input,& - cell,enpara,banddos,noco,oneD,hybrid,iter,eig_id,results,inden,v,vx) + !>The eigenvalue problem is constructed and solved in this routine. The following steps are performed: + !> 1. Preparation: generate energy parameters, open eig-file + !> 2. CALL to mt_setup() : this constructs the local Hamiltonian (i.e. the Hamiltonian in the \f$u,\dot u, u_{lo} \f$ basis) LDA+U is also added here + !> 3. within the (collinear)spin and k-point loop: CALL to eigen_hssetup() to generate the matrices, CALL to eigen_diag() to perform diagonalization + !> 4. writing (saving) of eigenvectors + !> + !> The matrices generated and diagonalized here are of type m_mat as defined in m_types_mat. + !>@author D. Wortmann + SUBROUTINE eigen(mpi,stars,sphhar,atoms,obsolete,xcpot,sym,kpts,DIMENSION,vacuum,input,& + cell,enpara,banddos,noco,oneD,hybrid,iter,eig_id,results,inden,v,vx) #include"cpp_double.h" - USE m_constants, ONLY : pi_const,sfp_const - USE m_types - USE m_apws - USE m_eigen_hssetup - USE m_pot_io - USE m_eigen_diag - USE m_add_vnonlocal - USE m_subvxc - !USE m_hsefunctional - USE m_mt_setup - USE m_util - USE m_io_hybrid - !USE m_icorrkeys - USE m_eig66_io, ONLY : open_eig, write_eig, close_eig,read_eig - USE m_xmlOutput + USE m_constants, ONLY : pi_const,sfp_const + USE m_types + USE m_apws + USE m_eigen_hssetup + USE m_pot_io + USE m_eigen_diag + USE m_add_vnonlocal + USE m_subvxc + !USE m_hsefunctional + USE m_mt_setup + USE m_util + USE m_io_hybrid + !USE m_icorrkeys + USE m_eig66_io, ONLY : open_eig, write_eig, close_eig,read_eig + USE m_xmlOutput #ifdef CPP_MPI - USE m_mpi_bc_potden + USE m_mpi_bc_potden #endif - USE m_symmetrize_matrix - USE m_unfold_band_kpts !used for unfolding bands - USE m_types_mpimat - - - IMPLICIT NONE - TYPE(t_results),INTENT(INOUT):: results - CLASS(t_xcpot),INTENT(IN) :: xcpot - TYPE(t_mpi),INTENT(IN) :: mpi - TYPE(t_dimension),INTENT(IN) :: DIMENSION - TYPE(t_oneD),INTENT(IN) :: oneD - TYPE(t_hybrid),INTENT(INOUT) :: hybrid - TYPE(t_enpara),INTENT(INOUT) :: enpara - TYPE(t_obsolete),INTENT(IN) :: obsolete - TYPE(t_input),INTENT(IN) :: input - TYPE(t_vacuum),INTENT(IN) :: vacuum - TYPE(t_noco),INTENT(IN) :: noco - TYPE(t_banddos),INTENT(IN) :: banddos - TYPE(t_sym),INTENT(IN) :: sym - TYPE(t_stars),INTENT(IN) :: stars - TYPE(t_cell),INTENT(IN) :: cell - TYPE(t_kpts),INTENT(INOUT) :: kpts - TYPE(t_sphhar),INTENT(IN) :: sphhar - TYPE(t_atoms),INTENT(IN) :: atoms - TYPE(t_potden),INTENT(IN) :: inden,vx - TYPE(t_potden),INTENT(INOUT) :: v !u_setup will modify the potential matrix - + USE m_symmetrize_matrix + USE m_unfold_band_kpts !used for unfolding bands + USE m_types_mpimat + + IMPLICIT NONE + TYPE(t_results),INTENT(INOUT):: results + CLASS(t_xcpot),INTENT(IN) :: xcpot + TYPE(t_mpi),INTENT(IN) :: mpi + TYPE(t_dimension),INTENT(IN) :: DIMENSION + TYPE(t_oneD),INTENT(IN) :: oneD + TYPE(t_hybrid),INTENT(INOUT) :: hybrid + TYPE(t_enpara),INTENT(INOUT) :: enpara + TYPE(t_obsolete),INTENT(IN) :: obsolete + TYPE(t_input),INTENT(IN) :: input + TYPE(t_vacuum),INTENT(IN) :: vacuum + TYPE(t_noco),INTENT(IN) :: noco + TYPE(t_banddos),INTENT(IN) :: banddos + TYPE(t_sym),INTENT(IN) :: sym + TYPE(t_stars),INTENT(IN) :: stars + TYPE(t_cell),INTENT(IN) :: cell + TYPE(t_kpts),INTENT(INOUT) :: kpts + TYPE(t_sphhar),INTENT(IN) :: sphhar + TYPE(t_atoms),INTENT(IN) :: atoms + TYPE(t_potden),INTENT(IN) :: inden,vx + TYPE(t_potden),INTENT(INOUT) :: v !u_setup will modify the potential matrix + #ifdef CPP_MPI - INCLUDE 'mpif.h' + INCLUDE 'mpif.h' #endif ! EXTERNAL MPI_BCAST !only used by band_unfolding to broadcast the gvec - ! Scalar Arguments - INTEGER,INTENT(IN) :: iter - INTEGER,INTENT(IN) :: eig_id - - ! Local Scalars - INTEGER jsp,nk,nred,ne_all,ne_found - INTEGER ne,lh0 - INTEGER isp,i,j,err - LOGICAL l_wu,l_file,l_real,l_zref - - ! Local Arrays - INTEGER :: ierr(3) - INTEGER :: neigBuffer(kpts%nkpt,input%jspins) - - COMPLEX :: unfoldingBuffer(SIZE(results%unfolding_weights,1),kpts%nkpt,input%jspins) ! needed for unfolding bandstructure mpi case - - INTEGER, PARAMETER :: lmaxb = 3 - REAL, ALLOCATABLE :: bkpt(:) - REAL, ALLOCATABLE :: eig(:) - COMPLEX, ALLOCATABLE :: vs_mmp(:,:,:,:) - - INTEGER :: jsp_m, i_kpt_m, i_m - - TYPE(t_tlmplm) :: td - TYPE(t_usdus) :: ud - TYPE(t_lapw) :: lapw - CLASS(t_mat), ALLOCATABLE :: zMat - CLASS(t_mat), ALLOCATABLE :: hmat,smat - CLASS(t_mat), ALLOCATABLE :: smat_unfold !used for unfolding bandstructure - - ! Variables for HF or hybrid functional calculation - INTEGER :: comm(kpts%nkpt),irank2(kpts%nkpt),isize2(kpts%nkpt), dealloc_stat - character(len=300) :: errmsg - - call ud%init(atoms,DIMENSION%jspd) - ALLOCATE (eig(DIMENSION%neigd),bkpt(3)) - - l_real=sym%invs.AND..NOT.noco%l_noco - - ! check if z-reflection trick can be used - l_zref=(sym%zrfs.AND.(SUM(ABS(kpts%bk(3,:kpts%nkpt))).LT.1e-9).AND..NOT.noco%l_noco) - IF (mpi%n_size > 1) l_zref = .FALSE. - + ! Scalar Arguments + INTEGER,INTENT(IN) :: iter + INTEGER,INTENT(IN) :: eig_id + + ! Local Scalars + INTEGER jsp,nk,nred,ne_all,ne_found + INTEGER ne,lh0 + INTEGER isp,i,j,err + LOGICAL l_wu,l_file,l_real,l_zref + + ! Local Arrays + INTEGER :: ierr(3) + INTEGER :: neigBuffer(kpts%nkpt,input%jspins) + + COMPLEX :: unfoldingBuffer(SIZE(results%unfolding_weights,1),kpts%nkpt,input%jspins) ! needed for unfolding bandstructure mpi case + + INTEGER, PARAMETER :: lmaxb = 3 + REAL, ALLOCATABLE :: bkpt(:) + REAL, ALLOCATABLE :: eig(:) + COMPLEX, ALLOCATABLE :: vs_mmp(:,:,:,:) + + INTEGER :: jsp_m, i_kpt_m, i_m + + TYPE(t_tlmplm) :: td + TYPE(t_usdus) :: ud + TYPE(t_lapw) :: lapw + CLASS(t_mat), ALLOCATABLE :: zMat + CLASS(t_mat), ALLOCATABLE :: hmat,smat + CLASS(t_mat), ALLOCATABLE :: smat_unfold !used for unfolding bandstructure + + ! Variables for HF or hybrid functional calculation + INTEGER :: comm(kpts%nkpt),irank2(kpts%nkpt),isize2(kpts%nkpt), dealloc_stat + character(len=300) :: errmsg + + call ud%init(atoms,DIMENSION%jspd) + ALLOCATE (eig(DIMENSION%neigd),bkpt(3)) + + l_real=sym%invs.AND..NOT.noco%l_noco + + ! check if z-reflection trick can be used + l_zref=(sym%zrfs.AND.(SUM(ABS(kpts%bk(3,:kpts%nkpt))).LT.1e-9).AND..NOT.noco%l_noco) + IF (mpi%n_size > 1) l_zref = .FALSE. + #ifdef CPP_MPI - CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,v) + CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,v) #endif - !IF (mpi%irank.EQ.0) CALL openXMLElementFormPoly('iteration',(/'numberForCurrentRun','overallNumber '/),(/iter,v%iter/),& - ! RESHAPE((/19,13,5,5/),(/2,2/))) - - ! Set up and solve the eigenvalue problem - ! loop over spins - ! set up k-point independent t(l'm',lm) matrices - CALL mt_setup(atoms,sym,sphhar,input,noco,enpara,inden,v,mpi,results,DIMENSION,td,ud) - - neigBuffer = 0 - results%neig = 0 - results%eig = 1.0e300 - unfoldingBuffer = CMPLX(0.0,0.0) - - DO jsp = 1,MERGE(1,input%jspins,noco%l_noco) - k_loop:DO nk = mpi%n_start,kpts%nkpt,mpi%n_stride - - ! Set up lapw list - CALL lapw%init(input,noco, kpts,atoms,sym,nk,cell,l_zref, mpi) - call timestart("Setup of H&S matrices") - CALL eigen_hssetup(jsp,mpi,DIMENSION,hybrid,enpara,input,vacuum,noco,sym,& - stars,cell,sphhar,atoms,ud,td,v,lapw,l_real,smat,hmat) - CALL timestop("Setup of H&S matrices") - - IF(hybrid%l_hybrid) THEN - - DO i = 1, hmat%matsize1 - DO j = 1, i - IF (hmat%l_real) THEN - IF ((i.LE.5).AND.(j.LE.5)) THEN - WRITE(1233,'(2i7,2f15.8)') i, j, hmat%data_r(i,j), hmat%data_r(j,i) - END IF - ELSE - ENDIF - END DO - END DO + !IF (mpi%irank.EQ.0) CALL openXMLElementFormPoly('iteration',(/'numberForCurrentRun','overallNumber '/),(/iter,v%iter/),& + ! RESHAPE((/19,13,5,5/),(/2,2/))) + + ! Set up and solve the eigenvalue problem + ! loop over spins + ! set up k-point independent t(l'm',lm) matrices + CALL mt_setup(atoms,sym,sphhar,input,noco,enpara,inden,v,mpi,results,DIMENSION,td,ud) - ! Write overlap matrix smat to direct access file olap - print *,"Wrong overlap matrix used, fix this later" - CALL write_olap(smat,(jsp-1)*kpts%nkpt+nk) ! Note: At this moment this only works without MPI parallelization - PRINT *,"TODO" + neigBuffer = 0 + results%neig = 0 + results%eig = 1.0e300 + unfoldingBuffer = CMPLX(0.0,0.0) + + DO jsp = 1,MERGE(1,input%jspins,noco%l_noco) + k_loop:DO nk = mpi%n_start,kpts%nkpt,mpi%n_stride + + ! Set up lapw list + CALL lapw%init(input,noco, kpts,atoms,sym,nk,cell,l_zref, mpi) + call timestart("Setup of H&S matrices") + CALL eigen_hssetup(jsp,mpi,DIMENSION,hybrid,enpara,input,vacuum,noco,sym,& + stars,cell,sphhar,atoms,ud,td,v,lapw,l_real,smat,hmat) + CALL timestop("Setup of H&S matrices") + + IF(hybrid%l_hybrid) THEN + + DO i = 1, hmat%matsize1 + DO j = 1, i + IF (hmat%l_real) THEN + IF ((i.LE.5).AND.(j.LE.5)) THEN + WRITE(1233,'(2i7,2f15.8)') i, j, hmat%data_r(i,j), hmat%data_r(j,i) + END IF + ELSE + ENDIF + END DO + END DO + + ! Write overlap matrix smat to direct access file olap + print *,"Wrong overlap matrix used, fix this later" + CALL write_olap(smat,(jsp-1)*kpts%nkpt+nk) ! Note: At this moment this only works without MPI parallelization + PRINT *,"TODO" ! STOP "TODO" - PRINT *,"BASIS:", lapw%nv(jsp), atoms%nlotot - IF (hybrid%l_addhf) CALL add_Vnonlocal(nk,lapw,atoms,hybrid,dimension,kpts,jsp,results,xcpot,noco,hmat) - - IF(hybrid%l_subvxc) THEN - CALL subvxc(lapw,kpts%bk(:,nk),DIMENSION,input,jsp,v%mt(:,0,:,:),atoms,ud,hybrid,enpara%el0,enpara%ello0,& - sym,cell,sphhar,stars,xcpot,mpi,oneD,hmat,vx) - END IF - END IF ! hybrid%l_hybrid - - l_wu=.FALSE. - ne_all=DIMENSION%neigd - if (allocated(zmat)) then - deallocate(zmat, stat=dealloc_stat, errmsg=errmsg) - if(dealloc_stat /= 0) call juDFT_error("deallocate failed for zmat",& - hint=errmsg, calledby="eigen.F90") - endif - - !Try to symmetrize matrix - CALL symmetrize_matrix(mpi,noco,kpts,nk,hmat,smat) - - IF (banddos%unfoldband) THEN - select type(smat) - type is (t_mat) - allocate(t_mat::smat_unfold) - select type(smat_unfold) - type is (t_mat) - smat_unfold=smat - end select - type is (t_mpimat) - allocate(t_mpimat::smat_unfold) - select type(smat_unfold) - type is (t_mpimat) - smat_unfold=smat - end select - end select - END IF - - CALL eigen_diag(mpi,hmat,smat,nk,jsp,iter,ne_all,eig,zMat) - CALL smat%free() - DEALLOCATE(hmat,smat, stat=dealloc_stat, errmsg=errmsg) - if(dealloc_stat /= 0) call juDFT_error("deallocate failed for hmat or smat",& - hint=errmsg, calledby="eigen.F90") - - ! Output results - CALL timestart("EV output") + PRINT *,"BASIS:", lapw%nv(jsp), atoms%nlotot + IF (hybrid%l_addhf) CALL add_Vnonlocal(nk,lapw,atoms,hybrid,dimension,kpts,jsp,results,xcpot,noco,hmat) + + IF(hybrid%l_subvxc) THEN + CALL subvxc(lapw,kpts%bk(:,nk),DIMENSION,input,jsp,v%mt(:,0,:,:),atoms,ud,hybrid,enpara%el0,enpara%ello0,& + sym,cell,sphhar,stars,xcpot,mpi,oneD,hmat,vx) + END IF + END IF ! hybrid%l_hybrid + + l_wu=.FALSE. + ne_all=DIMENSION%neigd + + !Try to symmetrize matrix + CALL symmetrize_matrix(mpi,noco,kpts,nk,hmat,smat) + + IF (banddos%unfoldband) THEN + select type(smat) + type is (t_mat) + allocate(t_mat::smat_unfold) + select type(smat_unfold) + type is (t_mat) + smat_unfold=smat + end select + type is (t_mpimat) + allocate(t_mpimat::smat_unfold) + select type(smat_unfold) + type is (t_mpimat) + smat_unfold=smat + end select + end select + END IF + + CALL eigen_diag(mpi,hmat,smat,nk,jsp,iter,ne_all,eig,zMat) + + CALL smat%free() + CALL hmat%free() + DEALLOCATE(hmat,smat, stat=dealloc_stat, errmsg=errmsg) + if(dealloc_stat /= 0) call juDFT_error("deallocate failed for hmat or smat",& + hint=errmsg, calledby="eigen.F90") + + ! Output results + CALL timestart("EV output") #if defined(CPP_MPI) - ! Collect number of all eigenvalues - ne_found=ne_all - CALL MPI_ALLREDUCE(ne_found,ne_all,1,MPI_INTEGER,MPI_SUM, mpi%sub_comm,ierr) - ne_all=MIN(DIMENSION%neigd,ne_all) + ! Collect number of all eigenvalues + ne_found=ne_all + CALL MPI_ALLREDUCE(ne_found,ne_all,1,MPI_INTEGER,MPI_SUM, mpi%sub_comm,ierr) + ne_all=MIN(DIMENSION%neigd,ne_all) #else - ne_found=ne_all -#endif - IF (.NOT.zmat%l_real) THEN - zMat%data_c(:lapw%nmat,:ne_found) = CONJG(zMat%data_c(:lapw%nmat,:ne_found)) - END IF - CALL write_eig(eig_id, nk,jsp,ne_found,ne_all,& - eig(:ne_found),n_start=mpi%n_size,n_end=mpi%n_rank,zmat=zMat) - neigBuffer(nk,jsp) = ne_found + ne_found=ne_all +#endif + IF (.NOT.zMat%l_real) THEN + zMat%data_c(:lapw%nmat,:ne_found) = CONJG(zMat%data_c(:lapw%nmat,:ne_found)) + END IF + CALL write_eig(eig_id, nk,jsp,ne_found,ne_all,& + eig(:ne_found),n_start=mpi%n_size,n_end=mpi%n_rank,zMat=zMat) + neigBuffer(nk,jsp) = ne_found #if defined(CPP_MPI) - ! RMA synchronization - CALL MPI_BARRIER(mpi%MPI_COMM,ierr) + ! RMA synchronization + CALL MPI_BARRIER(mpi%MPI_COMM,ierr) #endif - CALL timestop("EV output") + CALL timestop("EV output") - IF (banddos%unfoldband) THEN + IF (banddos%unfoldband) THEN CALL calculate_plot_w_n(banddos,cell,kpts,smat_unfold,zMat,lapw,nk,jsp,eig,results,input,atoms,unfoldingBuffer,mpi) - DEALLOCATE(smat_unfold, stat=dealloc_stat, errmsg=errmsg) - if(dealloc_stat /= 0) call juDFT_error("deallocate failed for smat_unfold",& - hint=errmsg, calledby="eigen.F90") - END IF + CALL smat_unfold%free() + DEALLOCATE(smat_unfold, stat=dealloc_stat, errmsg=errmsg) + if(dealloc_stat /= 0) call juDFT_error("deallocate failed for smat_unfold",& + hint=errmsg, calledby="eigen.F90") + END IF - END DO k_loop - END DO ! spin loop ends + call zMat%free() + deallocate(zMat) + END DO k_loop + END DO ! spin loop ends #ifdef CPP_MPI - IF (banddos%unfoldband) THEN - results%unfolding_weights = CMPLX(0.0,0.0) + IF (banddos%unfoldband) THEN + results%unfolding_weights = CMPLX(0.0,0.0) CALL MPI_ALLREDUCE(unfoldingBuffer,results%unfolding_weights,SIZE(results%unfolding_weights,1)*SIZE(results%unfolding_weights,2)*SIZE(results%unfolding_weights,3),CPP_MPI_COMPLEX,MPI_SUM,mpi%mpi_comm,ierr) - END IF - CALL MPI_ALLREDUCE(neigBuffer,results%neig,kpts%nkpt*input%jspins,MPI_INTEGER,MPI_SUM,mpi%sub_comm,ierr) - CALL MPI_BARRIER(mpi%MPI_COMM,ierr) + END IF + CALL MPI_ALLREDUCE(neigBuffer,results%neig,kpts%nkpt*input%jspins,MPI_INTEGER,MPI_SUM,mpi%mpi_comm,ierr) + CALL MPI_BARRIER(mpi%MPI_COMM,ierr) #else - results%neig(:,:) = neigBuffer(:,:) - results%unfolding_weights(:,:,:) = unfoldingBuffer(:,:,:) + results%neig(:,:) = neigBuffer(:,:) + results%unfolding_weights(:,:,:) = unfoldingBuffer(:,:,:) #endif - ! Sorry for the following strange workaround to fill the results%eig array. - ! At some point someone should have a closer look at how the eigenvalues are - ! distributed and fill the array without using the eigenvalue-IO. - DO jsp = 1,MERGE(1,input%jspins,noco%l_noco) - DO nk = 1,kpts%nkpt - CALL read_eig(eig_id,nk,jsp,results%neig(nk,jsp),results%eig(:,nk,jsp)) + ! Sorry for the following strange workaround to fill the results%eig array. + ! At some point someone should have a closer look at how the eigenvalues are + ! distributed and fill the array without using the eigenvalue-IO. + DO jsp = 1,MERGE(1,input%jspins,noco%l_noco) + DO nk = 1,kpts%nkpt + CALL read_eig(eig_id,nk,jsp,results%neig(nk,jsp),results%eig(:,nk,jsp)) #ifdef CPP_MPI - CALL MPI_BARRIER(mpi%MPI_COMM,ierr) + CALL MPI_BARRIER(mpi%MPI_COMM,ierr) #endif - END DO - END DO - - !IF (hybrid%l_hybrid.OR.hybrid%l_calhf) CALL close_eig(eig_id) - - IF( input%jspins .EQ. 1 .AND. hybrid%l_hybrid ) THEN - results%te_hfex%valence = 2*results%te_hfex%valence - IF(hybrid%l_calhf) results%te_hfex%core = 2*results%te_hfex%core - END IF - enpara%epara_min = MINVAL(enpara%el0) - enpara%epara_min = MIN(MINVAL(enpara%ello0),enpara%epara_min) - END SUBROUTINE eigen + END DO + END DO + + !IF (hybrid%l_hybrid.OR.hybrid%l_calhf) CALL close_eig(eig_id) + + IF( input%jspins .EQ. 1 .AND. hybrid%l_hybrid ) THEN + results%te_hfex%valence = 2*results%te_hfex%valence + IF(hybrid%l_calhf) results%te_hfex%core = 2*results%te_hfex%core + END IF + enpara%epara_min = MINVAL(enpara%el0) + enpara%epara_min = MIN(MINVAL(enpara%ello0),enpara%epara_min) + END SUBROUTINE eigen END MODULE m_eigen diff --git a/eigen/eigen_hssetup.F90 b/eigen/eigen_hssetup.F90 index 87206a91dbf264690553557a2abeff178b1e4195..2b8a0b76a518c251873973d696360193ea5eb2bb 100644 --- a/eigen/eigen_hssetup.F90 +++ b/eigen/eigen_hssetup.F90 @@ -97,7 +97,7 @@ CONTAINS ALLOCATE(hmat_final,mold=smat(1,1)) CALL eigen_redist_matrix(mpi,lapw,atoms,smat,smat_final) CALL eigen_redist_matrix(mpi,lapw,atoms,hmat,hmat_final,smat_final) - + END SUBROUTINE eigen_hssetup END MODULE m_eigen_hssetup diff --git a/io/banddos_io.F90 b/io/banddos_io.F90 index 77405a5e0f1d998bf2a15a07b750e5f36801a41a..4d2646b7d392b1351d18f83082924c1fec8e7fd3 100644 --- a/io/banddos_io.F90 +++ b/io/banddos_io.F90 @@ -27,16 +27,17 @@ MODULE m_banddos_io CONTAINS - SUBROUTINE openBandDOSFile(fileID, input, atoms, cell, kpts) + SUBROUTINE openBandDOSFile(fileID, input, atoms, cell, kpts, banddos) USE m_types USE hdf5 USE m_cdn_io - TYPE(t_input), INTENT(IN) :: input - TYPE(t_atoms), INTENT(IN) :: atoms - TYPE(t_cell), INTENT(IN) :: cell - TYPE(t_kpts), INTENT(IN) :: kpts + TYPE(t_input), INTENT(IN) :: input + TYPE(t_atoms), INTENT(IN) :: atoms + TYPE(t_cell), INTENT(IN) :: cell + TYPE(t_kpts), INTENT(IN) :: kpts + TYPE(t_banddos), INTENT(IN) :: banddos INTEGER(HID_T), INTENT(OUT) :: fileID @@ -67,6 +68,7 @@ MODULE m_banddos_io INTEGER :: hdfError, dimsInt(7) INTEGER :: version + INTEGER :: fakeLogical REAL :: eFermiPrev LOGICAL :: l_error @@ -98,6 +100,9 @@ MODULE m_banddos_io CALL h5gcreate_f(fileID, '/general', generalGroupID, hdfError) CALL io_write_attint0(generalGroupID,'spins',input%jspins) CALL io_write_attreal0(generalGroupID,'lastFermiEnergy',eFermiPrev) + fakeLogical = 0 + IF (banddos%unfoldband) fakeLogical = 1 + CALL io_write_attint0(generalGroupID,'bandUnfolding',fakeLogical) CALL h5gclose_f(generalGroupID, hdfError) CALL h5gcreate_f(fileID, '/cell', cellGroupID, hdfError) @@ -236,12 +241,14 @@ MODULE m_banddos_io INTEGER(HID_T), INTENT(IN) :: fileID INTEGER(HID_T) :: eigenvaluesGroupID + INTEGER(HID_T) :: bandUnfoldingGroupID INTEGER(HID_T) :: eigenvaluesSpaceID, eigenvaluesSetID INTEGER(HID_T) :: numFoundEigsSpaceID, numFoundEigsSetID INTEGER(HID_T) :: lLikeChargeSpaceID, lLikeChargeSetID INTEGER(HID_T) :: jsymSpaceID, jsymSetID INTEGER(HID_T) :: ksymSpaceID, ksymSetID + INTEGER(HID_T) :: bUWeightsSpaceID, bUWeightsSetID INTEGER :: hdfError, dimsInt(7) @@ -296,6 +303,20 @@ MODULE m_banddos_io CALL h5gclose_f(eigenvaluesGroupID, hdfError) + IF (banddos%unfoldband) THEN + CALL h5gcreate_f(fileID, '/bandUnfolding', bandUnfoldingGroupID, hdfError) + + dims(:3)=(/neigd,kpts%nkpt,input%jspins/) + dimsInt = dims + CALL h5screate_simple_f(3,dims(:3),bUWeightsSpaceID,hdfError) + CALL h5dcreate_f(bandUnfoldingGroupID, "weights", H5T_NATIVE_DOUBLE, bUWeightsSpaceID, buWeightsSetID, hdfError) + CALL h5sclose_f(bUWeightsSpaceID,hdfError) + CALL io_write_real3(bUWeightsSetID,(/1,1,1/),dimsInt(:3), REAL(results%unfolding_weights(:neigd,:,:))) + CALL h5dclose_f(bUWeightsSetID, hdfError) + + CALL h5gclose_f(bandUnfoldingGroupID, hdfError) + END IF + END SUBROUTINE SUBROUTINE io_write_string1(datasetID,dims,stringLength,dataArray) diff --git a/kpoints/unfoldBandKPTS.f90 b/kpoints/unfoldBandKPTS.f90 index 16f8d69f7ddc27b399e6929c10aa5d260686b3d5..baa6750db4bbfff98b576e870d1621dfc9fc4314 100644 --- a/kpoints/unfoldBandKPTS.f90 +++ b/kpoints/unfoldBandKPTS.f90 @@ -209,17 +209,20 @@ CONTAINS LOGICAL :: write_to_file=.false. CLASS(t_mat), ALLOCATABLE :: zMat_s - CALL build_primitive_cell(banddos,p_cell,cell) +! method_rubel=.true. !this switch is to switch between overlap matrix and rubel method (without overlap matrix) - DO j = 1, lapw%nv(jsp) - DO i = 1, j-1 - IF(smat_unfold%l_real) THEN - smat_unfold%data_r(j,i) = smat_unfold%data_r(i,j) - ELSE - smat_unfold%data_c(j,i) = CONJG(smat_unfold%data_c(i,j)) - END IF - END DO - END DO + CALL build_primitive_cell(banddos,p_cell,cell) + IF (.not. method_rubel) THEN + DO j = 1, lapw%nv(jsp) + DO i = 1, j-1 + IF(smat_unfold%l_real) THEN + smat_unfold%data_r(j,i) = smat_unfold%data_r(i,j) + ELSE + smat_unfold%data_c(j,i) = CONJG(smat_unfold%data_c(i,j)) + END IF + END DO + END DO + END IF ! write_to_file=.true. IF (write_to_file) THEN IF (i_kpt==1) THEN @@ -233,7 +236,6 @@ CONTAINS ! write(222,'(234f15.8)') zMat%data_r ! write(223,'(234f15.8)') smat_unfold%data_r -! method_rubel=.true. !this switch is to switch between overlap matrix and rubel method (without overlap matrix) IF (zmat%l_real) THEN ALLOCATE(w_n(zMat%matsize2)) diff --git a/main/cdngen.F90 b/main/cdngen.F90 index 6aa753f28383568aed1419eab5ba498ab6f05499..d320a06a764469d722abd45755deb409645f75ac 100644 --- a/main/cdngen.F90 +++ b/main/cdngen.F90 @@ -127,7 +127,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& CALL write_band_sc(kpts,results,eFermiPrev) END IF #ifdef CPP_HDF - CALL openBandDOSFile(banddosFile_id,input,atoms,cell,kpts) + CALL openBandDOSFile(banddosFile_id,input,atoms,cell,kpts,banddos) CALL writeBandDOSData(banddosFile_id,input,atoms,cell,kpts,results,banddos,dos,vacuum) CALL closeBandDOSFile(banddosFile_id) #endif diff --git a/types/types_mpimat.F90 b/types/types_mpimat.F90 index e4d24607fb9bf14020ce49c1038ae9aafd62e5c3..fa274e6f929c2e0ba274c38a98d332393b938821 100644 --- a/types/types_mpimat.F90 +++ b/types/types_mpimat.F90 @@ -46,7 +46,7 @@ MODULE m_types_mpimat PROCEDURE,PASS :: generate_full_matrix ! construct full matrix if only upper triangle of hermitian matrix is given PROCEDURE,PASS :: print_matrix PROCEDURE,PASS :: from_non_dist - FINAL :: finalize + FINAL :: finalize, finalize_1d, finalize_2d, finalize_3d END TYPE t_mpimat PUBLIC t_mpimat @@ -295,6 +295,44 @@ CONTAINS CALL mpimat_free(mat) END SUBROUTINE finalize + SUBROUTINE finalize_1d(mat) + IMPLICIT NONE + + TYPE(t_mpimat),INTENT(INOUT) :: mat(:) + INTEGER :: i + DO i = 1,size(mat) + CALL mpimat_free(mat(i)) + ENDDO + END SUBROUTINE finalize_1d + + SUBROUTINE finalize_2d(mat) + IMPLICIT NONE + + TYPE(t_mpimat),INTENT(INOUT) :: mat(:,:) + INTEGER :: i,j + + DO i = 1,size(mat, dim=1) + DO j = 1,size(mat, dim=2) + CALL mpimat_free(mat(i,j)) + ENDDO + ENDDO + END SUBROUTINE finalize_2d + + SUBROUTINE finalize_3d(mat) + IMPLICIT NONE + + TYPE(t_mpimat),INTENT(INOUT) :: mat(:,:,:) + INTEGER :: i,j,k + + DO i = 1,size(mat, dim=1) + DO j = 1,size(mat, dim=2) + DO k = 1,size(mat, dim=3) + CALL mpimat_free(mat(i,j,k)) + ENDDO + ENDDO + ENDDO + END SUBROUTINE finalize_3d + SUBROUTINE mpimat_free(mat) IMPLICIT NONE CLASS(t_mpimat),INTENT(INOUT) :: mat @@ -370,8 +408,10 @@ CONTAINS mat%blacsdata%blacs_desc(4)=global_size2 mat%global_size1=global_size1 mat%global_size2=global_size2 +#ifdef CPP_SCALAPACK mat%matsize1=NUMROC( global_size1,mat%blacsdata%blacs_desc(5), mat%blacsdata%myrow, mat%blacsdata%blacs_desc(7), mat%blacsdata%nprow ) mat%matsize1=NUMROC( global_size2,mat%blacsdata%blacs_desc(6), mat%blacsdata%mycol, mat%blacsdata%blacs_desc(8), mat%blacsdata%npcol ) +#endif ELSE mat%matsize1=templ%matsize1 mat%matsize2=templ%matsize2 @@ -467,12 +507,12 @@ CONTAINS k = k + 1 ENDDO ENDDO -#ifdef CPP_BLACSDEFAULT +!#ifdef CPP_BLACSDEFAULT !Get the Blacs default context CALL BLACS_GET(0,0,ictextblacs) -#else - ictextblacs=mpi_subcom -#endif +!#else +! ictextblacs=mpi_subcom +!#endif ! Create the Grid CALL BLACS_GRIDMAP(ictextblacs,iusermap,size(iusermap,1),blacsdata%nprow,blacsdata%npcol) ! Now control, whether the BLACS grid is the one we wanted