LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-07-03 09:55:59 +02:00
Code Issues Releases Wiki Activity

Fixed bug of dimensions in davidson

Browse Source
This commit is contained in:
Anthony Scemama 2016-10-07 12:22:43 +02:00
parent f7a2710f5c
commit b75e7dd2c2
2 changed files with 77 additions and 49 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

71
src/Davidson/davidson_parallel.irp.f
View File

@ -64,10 +64,10 @@ subroutine davidson_process(blockb, blocke, N, idx, vt, st)
st (:,org_j) = 0d0 st (:,org_j) = 0d0
end if end if
do istate=1,N_states_diag do istate=1,N_states_diag
vt (istate,org_i) += hij*dav_ut(istate,org_j) vt(istate,org_i) = vt(istate,org_i) + hij*dav_ut(istate,org_j)
st (istate,org_i) += s2*dav_ut(istate,org_j) st(istate,org_i) = st(istate,org_i) + s2 *dav_ut(istate,org_j)
vt (istate,org_j) += hij*dav_ut(istate,org_i) vt(istate,org_j) = vt(istate,org_j) + hij*dav_ut(istate,org_i)
st (istate,org_j) += s2*dav_ut(istate,org_i) st(istate,org_j) = st(istate,org_j) + s2 *dav_ut(istate,org_i)
enddo enddo
endif endif
enddo enddo
@ -117,13 +117,33 @@ subroutine davidson_collect(blockb, blocke, N, idx, vt, st , v0t, s0t)
end subroutine end subroutine
subroutine davidson_init(zmq_to_qp_run_socket) subroutine davidson_init(zmq_to_qp_run_socket,n,n_st_8,ut)
use f77_zmq use f77_zmq
implicit none implicit none
integer(ZMQ_PTR), intent(out) :: zmq_to_qp_run_socket integer(ZMQ_PTR), intent(out) :: zmq_to_qp_run_socket
integer, intent(in) :: n, n_st_8
double precision, intent(in) :: ut(n_st_8,n)
integer :: i,k
dav_size = n
touch dav_size
do i=1,n
do k=1,N_int
dav_det(k,1,i) = psi_det(k,1,i)
dav_det(k,2,i) = psi_det(k,2,i)
enddo
enddo
do i=1,n
do k=1,N_states_diag
dav_ut(k,i) = ut(k,i)
enddo
enddo
touch dav_det dav_ut
touch dav_size dav_det dav_ut
call new_parallel_job(zmq_to_qp_run_socket,"davidson") call new_parallel_job(zmq_to_qp_run_socket,"davidson")
end subroutine end subroutine
@ -309,15 +329,16 @@ end subroutine
subroutine davidson_collector(zmq_to_qp_run_socket, zmq_socket_pull , v0, s0) subroutine davidson_collector(zmq_to_qp_run_socket, zmq_socket_pull , v0, s0, LDA)
use f77_zmq use f77_zmq
implicit none implicit none
integer :: LDA
integer(ZMQ_PTR), intent(in) :: zmq_to_qp_run_socket integer(ZMQ_PTR), intent(in) :: zmq_to_qp_run_socket
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
double precision ,intent(inout) :: v0(dav_size, N_states_diag) double precision ,intent(inout) :: v0(LDA, N_states_diag)
double precision ,intent(inout) :: s0(dav_size, N_states_diag) double precision ,intent(inout) :: s0(LDA, N_states_diag)
integer :: more, task_id integer :: more, task_id
@ -353,10 +374,11 @@ subroutine davidson_collector(zmq_to_qp_run_socket, zmq_socket_pull , v0, s0)
end subroutine end subroutine
subroutine davidson_run(zmq_to_qp_run_socket , v0, s0) subroutine davidson_run(zmq_to_qp_run_socket , v0, s0, LDA)
use f77_zmq use f77_zmq
implicit none implicit none
integer :: LDA
integer(ZMQ_PTR), intent(in) :: zmq_to_qp_run_socket integer(ZMQ_PTR), intent(in) :: zmq_to_qp_run_socket
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
integer(ZMQ_PTR) :: zmq_collector integer(ZMQ_PTR) :: zmq_collector
@ -366,8 +388,8 @@ subroutine davidson_run(zmq_to_qp_run_socket , v0, s0)
integer :: i integer :: i
integer, external :: omp_get_thread_num integer, external :: omp_get_thread_num
double precision , intent(inout) :: v0(dav_size, N_states_diag) double precision , intent(inout) :: v0(LDA, N_states_diag)
double precision , intent(inout) :: s0(dav_size, N_states_diag) double precision , intent(inout) :: s0(LDA, N_states_diag)
call zmq_set_running(zmq_to_qp_run_socket) call zmq_set_running(zmq_to_qp_run_socket)
@ -381,7 +403,7 @@ subroutine davidson_run(zmq_to_qp_run_socket , v0, s0)
!$OMP PARALLEL DEFAULT(shared) private(i) num_threads(nproc+2) !$OMP PARALLEL DEFAULT(shared) private(i) num_threads(nproc+2)
i = omp_get_thread_num() i = omp_get_thread_num()
if (i==0) then if (i==0) then
call davidson_collector(zmq_collector, zmq_socket_pull , v0, s0) call davidson_collector(zmq_collector, zmq_socket_pull , v0, s0, size(v0,1))
call end_zmq_to_qp_run_socket(zmq_collector) call end_zmq_to_qp_run_socket(zmq_collector)
call end_zmq_pull_socket(zmq_socket_pull) call end_zmq_pull_socket(zmq_socket_pull)
call davidson_miniserver_end() call davidson_miniserver_end()
@ -471,15 +493,28 @@ end subroutine
BEGIN_PROVIDER [ integer(bit_kind), dav_det, (N_int, 2, dav_size) ] BEGIN_PROVIDER [ integer(bit_kind), dav_det, (N_int, 2, dav_size) ]
END_PROVIDER &BEGIN_PROVIDER [ double precision, dav_ut, (N_states_diag, dav_size) ]
use bitmasks
implicit none
BEGIN_PROVIDER [ double precision, dav_ut, (N_states_diag, dav_size) ] BEGIN_DOC
! Temporary arrays for parallel davidson
!
! Touched in davidson_miniserver_get
END_DOC
dav_det = 0_bit_kind
dav_ut = -huge(1.d0)
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ integer, dav_size ] BEGIN_PROVIDER [ integer, dav_size ]
implicit none
BEGIN_DOC
! Size of the arrays for Davidson
!
! Touched in davidson_miniserver_get
END_DOC
dav_size = 1
END_PROVIDER END_PROVIDER

27
src/Davidson/u0Hu0.irp.f
View File

@ -56,10 +56,10 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
integer :: N_st_8 integer :: N_st_8
integer, external :: align_double integer, external :: align_double
!!!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: vt, ut !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: vt, ut
if(N_st /= N_states_diag) stop "H_u_0_nstates N_st /= N_states_diag" if(N_st /= N_states_diag) stop "H_u_0_nstates N_st /= N_states_diag"
N_st_8 = N_states_diag ! align_double(N_st) N_st_8 = align_double(N_st)
ASSERT (Nint > 0) ASSERT (Nint > 0)
ASSERT (Nint == N_int) ASSERT (Nint == N_int)
@ -165,7 +165,7 @@ subroutine H_u_0_nstates(v_0,u_0,H_jj,n,keys_tmp,Nint,N_st,sze_8)
v_0(i,istate) += H_jj(i) * u_0(i,istate) v_0(i,istate) += H_jj(i) * u_0(i,istate)
enddo enddo
enddo enddo
!deallocate (shortcut, sort_idx, sorted, version, ut) deallocate (shortcut, sort_idx, sorted, version, ut)
end end
BEGIN_PROVIDER [ double precision, psi_energy, (N_states) ] BEGIN_PROVIDER [ double precision, psi_energy, (N_states) ]
@ -210,12 +210,12 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
integer, external :: align_double integer, external :: align_double
integer :: workload, blockb, blocke integer :: workload, blockb, blocke
!!!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: vt, ut ! !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: vt, ut
integer(ZMQ_PTR) :: handler integer(ZMQ_PTR) :: handler
if(N_st /= N_states_diag .or. sze_8 < N_det) stop "assert fail in H_S2_u_0_nstates" if(N_st /= N_states_diag .or. sze_8 < N_det) stop "assert fail in H_S2_u_0_nstates"
N_st_8 = N_st !! align_double(N_st) N_st_8 = N_st ! align_double(N_st)
ASSERT (Nint > 0) ASSERT (Nint > 0)
ASSERT (Nint == N_int) ASSERT (Nint == N_int)
@ -228,8 +228,6 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
v_0 = 0.d0 v_0 = 0.d0
s_0 = 0.d0 s_0 = 0.d0
if(n /= N_det) stop "n /= N_det"
do i=1,n do i=1,n
do istate=1,N_st do istate=1,N_st
ut(istate,i) = u_0(i,istate) ut(istate,i) = u_0(i,istate)
@ -238,15 +236,10 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
call sort_dets_ab_v(keys_tmp, sorted(1,1,1), sort_idx(1,1), shortcut(0,1), version(1,1,1), n, Nint) call sort_dets_ab_v(keys_tmp, sorted(1,1,1), sort_idx(1,1), shortcut(0,1), version(1,1,1), n, Nint)
call sort_dets_ba_v(keys_tmp, sorted(1,1,2), sort_idx(1,2), shortcut(0,2), version(1,1,2), n, Nint) call sort_dets_ba_v(keys_tmp, sorted(1,1,2), sort_idx(1,2), shortcut(0,2), version(1,1,2), n, Nint)
dav_size = n
touch dav_size
dav_det = psi_det
dav_ut = ut
workload = 0 workload = 0
blockb = shortcut(0,1) blockb = shortcut(0,1)
blocke = blockb blocke = blockb
call davidson_init(handler) call davidson_init(handler,n,N_st_8,ut)
do sh=shortcut(0,1),1,-1 do sh=shortcut(0,1),1,-1
workload += (shortcut(sh+1,1) - shortcut(sh,1))**2 workload += (shortcut(sh+1,1) - shortcut(sh,1))**2
if(workload > 100000) then if(workload > 100000) then
@ -258,7 +251,7 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
enddo enddo
if(blockb > 0) call davidson_add_task(handler, 1, blockb) if(blockb > 0) call davidson_add_task(handler, 1, blockb)
call davidson_run(handler, v_0, s_0) call davidson_run(handler, v_0, s_0, size(v_0,1))
!$OMP PARALLEL DEFAULT(NONE) & !$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,hij,s2,j,k,jj,vt,st,ii,sh,sh2,ni,exa,ext,org_i,org_j,endi,sorted_i,istate)& !$OMP PRIVATE(i,hij,s2,j,k,jj,vt,st,ii,sh,sh2,ni,exa,ext,org_i,org_j,endi,sorted_i,istate)&
@ -284,8 +277,8 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
do istate=1,n_st do istate=1,n_st
vt (istate,org_i) = vt (istate,org_i) + hij*ut(istate,org_j) vt (istate,org_i) = vt (istate,org_i) + hij*ut(istate,org_j)
vt (istate,org_j) = vt (istate,org_j) + hij*ut(istate,org_i) vt (istate,org_j) = vt (istate,org_j) + hij*ut(istate,org_i)
st (istate,org_i) = st (istate,org_i) + s2*ut(istate,org_j) st (istate,org_i) = st (istate,org_i) + s2 *ut(istate,org_j)
st (istate,org_j) = st (istate,org_j) + s2*ut(istate,org_i) st (istate,org_j) = st (istate,org_j) + s2 *ut(istate,org_i)
enddo enddo
end if end if
end do end do
@ -295,7 +288,7 @@ subroutine H_S2_u_0_nstates(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze_8)
!$OMP CRITICAL !$OMP CRITICAL
do istate=1,N_st do istate=1,N_st
do i=n,1,-1 do i=1,n
v_0(i,istate) = v_0(i,istate) + vt(istate,i) v_0(i,istate) = v_0(i,istate) + vt(istate,i)
s_0(i,istate) = s_0(i,istate) + st(istate,i) s_0(i,istate) = s_0(i,istate) + st(istate,i)
enddo enddo