Commit 6cfa7e17 authored by Matthias Redies's avatar Matthias Redies

this should finally fix the branch

parent fd7799fa
......@@ -49,9 +49,6 @@ CONTAINS
REAL, ALLOCATABLE :: eig_irr(:, :)
CALL timestart("hybrid code")
call hybmpi%copy_mpi(mpi)
call work_pack%init(fi, hybmpi%rank, hybmpi%size)
INQUIRE (file="v_x.1", exist=hybdat%l_addhf)
......@@ -67,6 +64,7 @@ CONTAINS
IF (.NOT. hybdat%l_calhf) THEN
hybdat%l_subvxc = hybdat%l_subvxc .AND. hybdat%l_addhf
else
call hybmpi%copy_mpi(mpi)
results%te_hfex%core = 0
!Check if we are converged well enough to calculate a new potential
......@@ -77,7 +75,6 @@ CONTAINS
call first_iteration_alloc(fi, hybdat)
init_vex = .FALSE.
END IF
hybdat%l_subvxc = (hybdat%l_subvxc .AND. hybdat%l_addhf)
IF (.NOT. ALLOCATED(results%w_iks)) allocate(results%w_iks(fi%input%neig, fi%kpts%nkpt, fi%input%jspins))
......@@ -104,7 +101,12 @@ CONTAINS
call hybdat%coul(i)%alloc(fi, mpdata%num_radbasfn, mpdata%n_g, i)
enddo
! use jsp=1 for coulomb work-planning
call hybdat%set_states(fi, results, 1)
call work_pack%init(fi, hybdat, 1, hybmpi%rank, hybmpi%size)
CALL coulombmatrix(mpi, fi, mpdata, hybdat, xcpot, work_pack)
call work_pack%free()
do i =1,fi%kpts%nkpt
call hybdat%coul(i)%mpi_ibc(fi, hybmpi, work_pack%owner_nk(i))
......@@ -116,16 +118,15 @@ CONTAINS
CALL timestart("Calculation of non-local HF potential")
DO jsp = 1, fi%input%jspins
call timestart("HF_setup")
CALL HF_setup(mpdata,fi%hybinp, fi%input, fi%sym, fi%kpts, fi%atoms, &
mpi, fi%noco, nococonv,fi%cell, fi%oneD, results, jsp, enpara, &
hybdat, fi%sym%invs, v%mt(:, 0, :, :), eig_irr)
CALL HF_setup(mpdata,fi, mpi, nococonv, results, jsp, enpara, &
hybdat, v%mt(:, 0, :, :), eig_irr)
call work_pack%init(fi, hybdat, jsp, hybmpi%rank, hybmpi%size)
call timestop("HF_setup")
DO i = 1,work_pack%k_packs(1)%size
nk = work_pack%k_packs(i)%nk
CALL lapw%init(fi%input, fi%noco, nococonv,fi%kpts, fi%atoms, fi%sym, nk, fi%cell, l_zref)
CALL hsfock(fi,nk, mpdata, lapw, jsp, hybdat, eig_irr, &
CALL hsfock(fi, work_pack%k_packs(i), mpdata, lapw, jsp, hybdat, eig_irr, &
nococonv, stars, results, xcpot, mpi)
END DO
END DO
......
......@@ -59,7 +59,7 @@ MODULE m_exchange_valence_hf
INTEGER, PARAMETER:: maxmem = 600
CONTAINS
SUBROUTINE exchange_valence_hf(ik, fi, z_k, c_phase_k, mpdata, jsp, hybdat, lapw, eig_irr, results, &
SUBROUTINE exchange_valence_hf(k_pack, fi, z_k, c_phase_k, mpdata, jsp, hybdat, lapw, eig_irr, results, &
n_q, wl_iks, xcpot, nococonv, stars, nsest, indx_sest, mpi_var, mat_ex)
USE m_wrapper
......@@ -71,12 +71,14 @@ CONTAINS
USE m_io_hybinp
USE m_kp_perturbation
use m_spmm
use m_work_package
#ifdef CPP_MPI
use mpi
#endif
IMPLICIT NONE
type(t_fleurinput), intent(in) :: fi
type(t_k_package), intent(in) :: k_pack
type(t_mat), intent(in) :: z_k
TYPE(t_results), INTENT(IN) :: results
TYPE(t_xcpot_inbuild), INTENT(IN) :: xcpot
......@@ -94,7 +96,6 @@ CONTAINS
! scalars
INTEGER, INTENT(IN) :: jsp
INTEGER, INTENT(IN) :: ik
! arrays
INTEGER, INTENT(IN) :: n_q(:)
......@@ -103,16 +104,15 @@ CONTAINS
REAL, INTENT(IN) :: eig_irr(:, :)
REAL, INTENT(IN) :: wl_iks(:, :)
complex, intent(in) :: c_phase_k(hybdat%nbands(ik))
complex, intent(in) :: c_phase_k(hybdat%nbands(k_pack%nk))
! local scalars
INTEGER :: iband, iband1, jq, iq
INTEGER :: i, ierr
INTEGER :: i, ierr, ik
INTEGER :: j, iq_p
INTEGER :: n1, n2, nn2, cnt_read_z
INTEGER :: ikqpt, iob, m,n,k,lda,ldb,ldc
INTEGER :: ok, psize, n_parts, ipart, ibando
integer, allocatable :: start_idx(:), psizes(:)
REAL, SAVE :: divergence
......@@ -123,18 +123,19 @@ CONTAINS
! local arrays
COMPLEX :: exchcorrect(fi%kpts%nkptf)
COMPLEX :: dcprod(hybdat%nbands(ik), hybdat%nbands(ik), 3)
COMPLEX :: exch_vv(hybdat%nbands(ik), hybdat%nbands(ik))
COMPLEX :: dcprod(hybdat%nbands(k_pack%nk), hybdat%nbands(k_pack%nk), 3)
COMPLEX :: exch_vv(hybdat%nbands(k_pack%nk), hybdat%nbands(k_pack%nk))
COMPLEX :: hessian(3, 3), ctmp
COMPLEX :: proj_ibsc(3, MAXVAL(hybdat%nobd(:, jsp)), hybdat%nbands(ik))
COMPLEX :: proj_ibsc(3, MAXVAL(hybdat%nobd(:, jsp)), hybdat%nbands(k_pack%nk))
COMPLEX :: olap_ibsc(3, 3, MAXVAL(hybdat%nobd(:, jsp)), MAXVAL(hybdat%nobd(:, jsp)))
COMPLEX, ALLOCATABLE :: phase_vv(:, :)
REAL :: kqpt(3), kqpthlp(3), target_psize, rtmp
REAL :: kqpt(3), kqpthlp(3), rtmp
LOGICAL :: occup(fi%input%neig), conjg_mtir
type(t_mat) :: carr1_v, cprod_vv, carr3_vv, dot_result
character(len=300) :: errmsg
CALL timestart("valence exchange calculation")
ik = k_pack%nk
IF (initialize) THEN !it .eq. 1 .and. ik .eq. 1) THEN
call calc_divergence(fi%cell, fi%kpts, divergence)
......@@ -156,23 +157,16 @@ CONTAINS
cnt_read_z = predict_max_read_z(fi, hybdat, jsp)
#endif
DO jq = 1,fi%kpts%EIBZ(ik)%nkpt
iq = fi%kpts%EIBZ(ik)%pointer(jq)
iq = k_pack%q_packs(jq)%ptr
iq_p = fi%kpts%bkp(iq)
ikqpt = fi%kpts%get_nk(fi%kpts%to_first_bz(fi%kpts%bkf(:,ik) + fi%kpts%bkf(:,iq)))
! arrays should be less than 5 gb
if(mat_ex%l_real) then
target_psize = 5e9/( 8.0 * maxval(hybdat%nbasm) * hybdat%nbands(ik))
else
target_psize = 5e9/(16.0 * maxval(hybdat%nbasm) * hybdat%nbands(ik))
endif
n_parts = ceiling(hybdat%nobd(ikqpt, jsp)/target_psize)
call split_iob_loop(hybdat, hybdat%nobd(ikqpt, jsp), n_parts, start_idx, psizes)
n_parts = size(k_pack%q_packs(jq)%band_packs)
do ipart = 1, n_parts
if(n_parts > 1) write (*,*) "Part (" // int2str(ipart) //"/"// int2str(n_parts) // ")"
psize = psizes(ipart)
ibando = start_idx(ipart)
psize = k_pack%q_packs(jq)%band_packs(ipart)%psize
ibando = k_pack%q_packs(jq)%band_packs(ipart)%start_idx
call cprod_vv%alloc(mat_ex%l_real, hybdat%nbasm(iq), psize * hybdat%nbands(ik))
IF (mat_ex%l_real) THEN
......@@ -511,36 +505,6 @@ CONTAINS
enddo
end function calc_divergence2
subroutine split_iob_loop(hybdat, n_total, n_parts, start_idx, psize)
use m_types
implicit none
type(t_hybdat), intent(inout) :: hybdat
integer, intent(in) :: n_total, n_parts
integer, allocatable, intent(inout) :: start_idx(:), psize(:)
integer :: n_loops, i, big_size, small_size, end_idx
if(allocated(start_idx)) deallocate(start_idx)
if(allocated(psize)) deallocate(psize)
allocate(start_idx(n_parts), psize(n_parts))
small_size = floor((1.0*n_total)/n_parts)
big_size = small_size +1
end_idx = 0
do i = 1,n_parts
psize(i) = merge(big_size, small_size,i <= mod(n_total, n_parts))
start_idx(i) = end_idx + 1
end_idx = start_idx(i) + psize(i) - 1
enddo
if(hybdat%l_print_iob_splitting) then
write (*,*) "Split iob loop into " // int2str(n_parts) // " parts"
write (*,*) "sizes: ", psize(1), psize(n_parts)
hybdat%l_print_iob_splitting = .False.
endif
end subroutine split_iob_loop
subroutine recombine_parts(in_part, ipart, psizes, out_total)
use m_types
type(t_mat), intent(in) :: in_part
......
......@@ -10,7 +10,7 @@ MODULE m_hf_setup
CONTAINS
SUBROUTINE hf_setup(mpdata, fi, mpi,nococonv, results, jsp, enpara, &
hybdat, l_real, vr0, eig_irr)
hybdat, vr0, eig_irr)
USE m_types
USE m_constants
USE m_eig66_io
......@@ -33,7 +33,6 @@ CONTAINS
INTEGER, INTENT(IN) :: jsp
REAL, INTENT(IN) :: vr0(:, :, :)
LOGICAL, INTENT(IN) :: l_real
REAL, ALLOCATABLE, INTENT(INOUT) :: eig_irr(:, :)
......@@ -53,7 +52,6 @@ CONTAINS
REAL :: zDebug_r(lapw_dim_nbasfcn,fi%input%neig)
COMPLEX :: zDebug_c(lapw_dim_nbasfcn,fi%input%neig)
call hybdat%set_states(fi, results, jsp)
IF (hybdat%l_calhf) THEN
! Preparations for HF and hybinp functional calculation
......
......@@ -3,12 +3,17 @@ module m_work_package
implicit none
type t_band_package
integer :: start_idx, end_idx
integer :: start_idx, psize, rank, size
contains
procedure :: init => t_band_package_init
end type t_band_package
type t_q_package
integer :: n_q
integer :: comm_q_group
integer :: rank, size, ptr
type(t_band_package), allocatable :: band_packs(:)
contains
procedure :: init => t_q_package_init
procedure :: free => t_q_package_free
end type t_q_package
type t_k_package
......@@ -17,6 +22,7 @@ module m_work_package
contains
procedure :: init => t_k_package_init
procedure :: print => t_k_package_print
procedure :: free => t_k_package_free
end type t_k_package
type t_work_package
......@@ -27,30 +33,117 @@ module m_work_package
procedure :: print => t_work_package_print
procedure :: owner_nk => t_work_package_owner_nk
procedure :: has_nk => t_work_package_has_nk
procedure :: free => t_work_package_free
end type t_work_package
contains
subroutine t_work_package_init(work_pack, fi, rank, size)
subroutine t_work_package_free(work_pack)
implicit none
class(t_work_package), intent(inout) :: work_pack
integer :: i
if(allocated(work_pack%k_packs)) then
do i = 1, size(work_pack%k_packs)
call work_pack%k_packs(i)%free()
enddo
deallocate(work_pack%k_packs)
endif
end subroutine t_work_package_free
subroutine t_k_package_free(k_pack)
implicit none
class(t_k_package), intent(inout) :: k_pack
integer :: i
if(allocated(k_pack%q_packs)) then
do i = 1, size(k_pack%q_packs)
call k_pack%q_packs(i)%free()
enddo
deallocate(k_pack%q_packs)
endif
end subroutine t_k_package_free
subroutine t_q_package_free(q_pack)
implicit none
class(t_q_package), intent(inout) :: q_pack
if(allocated(q_pack%band_packs)) deallocate(q_pack%band_packs)
end subroutine t_q_package_free
subroutine t_work_package_init(work_pack, fi, hybdat, jsp, rank, size)
implicit none
class(t_work_package), intent(inout) :: work_pack
type(t_fleurinput), intent(in) :: fi
integer, intent(in) :: rank, size
type(t_hybdat), intent(in) :: hybdat
integer, intent(in) :: rank, size, jsp
work_pack%rank = rank
work_pack%size = size
call split_into_work_packages(work_pack, fi)
call split_into_work_packages(work_pack, fi, hybdat, jsp)
end subroutine
subroutine t_k_package_init(k_pack, fi, nk)
subroutine t_k_package_init(k_pack, fi, hybdat, jsp, nk)
implicit none
class(t_k_package), intent(inout) :: k_pack
type(t_fleurinput), intent(in) :: fi
integer, intent(in) :: nk
integer :: nkpt_eibz
type(t_fleurinput), intent(in) :: fi
type(t_hybdat), intent(in) :: hybdat
integer, intent(in) :: nk, jsp
integer :: iq, jq
k_pack%nk = nk
allocate(k_pack%q_packs(fi%kpts%EIBZ(nk)%nkpt))
end subroutine
do iq = 1,fi%kpts%EIBZ(nk)%nkpt
jq = fi%kpts%EIBZ(nk)%pointer(iq)
call k_pack%q_packs(iq)%init(fi, hybdat, jsp, nk, iq, jq)
enddo
end subroutine t_k_package_init
subroutine t_q_package_init(q_pack, fi, hybdat, jsp, nk, rank, ptr)
implicit none
class(t_q_package), intent(inout) :: q_pack
type(t_fleurinput), intent(in) :: fi
type(t_hybdat), intent(in) :: hybdat
integer, intent(in) :: rank, ptr, jsp, nk
real :: target_psize
integer :: n_parts, ikqpt, i
integer, allocatable :: start_idx(:), psize(:)
q_pack%rank = rank
q_pack%size = fi%kpts%EIBZ(nk)%nkpt
q_pack%ptr = ptr
! arrays should be less than 5 gb
if(fi%sym%invs) then
target_psize = 5e9/( 8.0 * maxval(hybdat%nbasm) * MIN(fi%hybinp%bands1, fi%input%neig))
else
target_psize = 5e9/(16.0 * maxval(hybdat%nbasm) * MIN(fi%hybinp%bands1, fi%input%neig))
endif
ikqpt = fi%kpts%get_nk(fi%kpts%to_first_bz(fi%kpts%bkf(:,nk) + fi%kpts%bkf(:,ptr)))
n_parts = ceiling(hybdat%nobd(ikqpt, jsp)/target_psize)
allocate(start_idx(n_parts), psize(n_parts))
allocate(q_pack%band_packs(n_parts))
call split_band_loop(hybdat%nobd(ikqpt, jsp), n_parts, start_idx, psize)
do i = 1, n_parts
call q_pack%band_packs(i)%init(start_idx(i), psize(i), i, n_parts)
enddo
end subroutine t_q_package_init
subroutine t_band_package_init(band_pack, start_idx, psize, rank, size)
implicit none
class(t_band_package), intent(inout) :: band_pack
integer, intent(in) :: rank, size, start_idx, psize
band_pack%start_idx = start_idx
band_pack%psize = psize
band_pack%rank = rank
band_pack%size = size
end subroutine t_band_package_init
subroutine t_work_package_print(work_pack)
implicit none
......@@ -71,10 +164,12 @@ contains
write (*,*) "nk = ", k_pack%nk
end subroutine t_k_package_print
subroutine split_into_work_packages(work_pack, fi)
subroutine split_into_work_packages(work_pack, fi, hybdat, jsp)
implicit none
class(t_work_package), intent(inout) :: work_pack
type(t_fleurinput), intent(in) :: fi
type(t_hybdat), intent(in) :: hybdat
integer, intent(in) :: jsp
integer :: my_num_ks, k_cnt, i
if(fi%kpts%nkpt < work_pack%size) call judft_error("not enough k-points for work_packages")
......@@ -94,7 +189,7 @@ contains
work_pack%k_packs(k_cnt)%rank = k_cnt -1
work_pack%k_packs(k_cnt)%size = my_num_ks
call work_pack%k_packs(k_cnt)%init(fi, i)
call work_pack%k_packs(k_cnt)%init(fi, hybdat, jsp, i)
k_cnt = k_cnt + 1
enddo
end subroutine split_into_work_packages
......@@ -125,4 +220,28 @@ contains
endif
enddo
end function t_work_package_has_nk
subroutine split_band_loop(n_total, n_parts, start_idx, psize)
use m_types
implicit none
integer, intent(in) :: n_total, n_parts
integer, allocatable, intent(inout) :: start_idx(:), psize(:)
integer :: i, big_size, small_size, end_idx
if(allocated(start_idx)) deallocate(start_idx)
if(allocated(psize)) deallocate(psize)
allocate(start_idx(n_parts), psize(n_parts))
small_size = floor((1.0*n_total)/n_parts)
big_size = small_size +1
end_idx = 0
do i = 1,n_parts
psize(i) = merge(big_size, small_size,i <= mod(n_total, n_parts))
start_idx(i) = end_idx + 1
end_idx = start_idx(i) + psize(i) - 1
enddo
end subroutine split_band_loop
end module m_work_package
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment