10
0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-03 20:53:54 +01:00

Sigma-vector seems to be working.

This commit is contained in:
v1j4y 2022-06-06 17:29:18 +02:00
parent 48e48fbf16
commit 1c41bfdd17
2 changed files with 124 additions and 122 deletions

View File

@ -226,7 +226,7 @@ subroutine generate_all_singles_cfg(cfg,singles,n_singles,Nint)
enddo enddo
end end
subroutine generate_all_singles_cfg_with_type(cfgInp,singles,idxs_singles,pq_singles,ex_type_singles,n_singles,Nint) subroutine generate_all_singles_cfg_with_type(bit_tmp,cfgInp,singles,idxs_singles,pq_singles,ex_type_singles,n_singles,Nint)
implicit none implicit none
use bitmasks use bitmasks
BEGIN_DOC BEGIN_DOC
@ -238,6 +238,7 @@ subroutine generate_all_singles_cfg_with_type(cfgInp,singles,idxs_singles,pq_sin
! ex_type_singles : on output contains type of excitations : ! ex_type_singles : on output contains type of excitations :
! !
END_DOC END_DOC
integer*8, intent(in) :: bit_tmp(0:N_configuration+1)
integer, intent(in) :: Nint integer, intent(in) :: Nint
integer, intent(inout) :: n_singles integer, intent(inout) :: n_singles
integer, intent(out) :: idxs_singles(*) integer, intent(out) :: idxs_singles(*)
@ -251,6 +252,7 @@ subroutine generate_all_singles_cfg_with_type(cfgInp,singles,idxs_singles,pq_sin
integer(bit_kind) :: single(Nint,2) integer(bit_kind) :: single(Nint,2)
logical :: i_ok logical :: i_ok
n_singles = 0 n_singles = 0
!TODO !TODO
!Make list of Somo and Domo for holes !Make list of Somo and Domo for holes
@ -261,7 +263,7 @@ subroutine generate_all_singles_cfg_with_type(cfgInp,singles,idxs_singles,pq_sin
addcfg = -1 addcfg = -1
call do_single_excitation_cfg_with_type(cfgInp,single,i_hole,i_particle,ex_type,i_ok) call do_single_excitation_cfg_with_type(cfgInp,single,i_hole,i_particle,ex_type,i_ok)
if (i_ok) then if (i_ok) then
call binary_search_cfg(single,addcfg) call binary_search_cfg(single,addcfg,bit_tmp)
if(addcfg .EQ. -1) cycle if(addcfg .EQ. -1) cycle
n_singles = n_singles + 1 n_singles = n_singles + 1
do k=1,Nint do k=1,Nint

View File

@ -21,7 +21,7 @@
! required for the calculation of prototype arrays. ! required for the calculation of prototype arrays.
END_DOC END_DOC
NSOMOMax = min(elec_num, cfg_nsomo_max + 2) NSOMOMax = min(elec_num, cfg_nsomo_max + 2)
NSOMOMin = 0 NSOMOMin = max(0,cfg_nsomo_min-2)
! Note that here we need NSOMOMax + 2 sizes ! Note that here we need NSOMOMax + 2 sizes
NCSFMax = max(1,nint((binom(NSOMOMax,(NSOMOMax+1)/2)-binom(NSOMOMax,((NSOMOMax+1)/2)+1)))) ! TODO: NCSFs for MS=0 NCSFMax = max(1,nint((binom(NSOMOMax,(NSOMOMax+1)/2)-binom(NSOMOMax,((NSOMOMax+1)/2)+1)))) ! TODO: NCSFs for MS=0
NBFMax = NCSFMax NBFMax = NCSFMax
@ -143,6 +143,117 @@ end subroutine get_phase_qp_to_cfg
BEGIN_PROVIDER [ real*8, DetToCSFTransformationMatrix, (0:NSOMOMax,NBFMax,maxDetDimPerBF)]
&BEGIN_PROVIDER [ real*8, psi_coef_config, (n_CSF,1)]
&BEGIN_PROVIDER [ integer, psi_config_data, (N_configuration,2)]
&BEGIN_PROVIDER [ integer, psi_csf_to_config_data, (n_CSF)]
use cfunctions
implicit none
BEGIN_DOC
! Documentation for DetToCSFTransformationMatrix
! Provides the matrix of transformatons for the
! conversion between determinant to CSF basis (in BFs)
END_DOC
integer*8 :: Isomo, Idomo
integer(bit_kind) :: Ialpha(N_int),Ibeta(N_int)
integer :: rows, cols, i, j, k
integer :: startdet, enddet, idx
integer*8 MS
integer ndetI
integer :: getNSOMO
real*8,dimension(:,:),allocatable :: tempBuffer
real*8,dimension(:),allocatable :: tempCoeff
real*8 :: norm_det1, phasedet
integer :: nt
norm_det1 = 0.d0
MS = elec_alpha_num - elec_beta_num
print *,"Maxbfdim=",NBFMax
print *,"Maxdetdim=",maxDetDimPerBF
print *,"n_CSF=",n_CSF
print *,"N_configurations=",N_configuration
print *,"n_core_orb=",n_core_orb
! initialization
psi_coef_config = 0.d0
DetToCSFTransformationMatrix(0,:,:) = 1.d0
do i = 2-iand(elec_alpha_num-elec_beta_num,1), NSOMOMax,2
Isomo = IBSET(0_8, i) - 1_8
! rows = Ncsfs
! cols = Ndets
bfIcfg = max(1,nint((binom(i,(i+1)/2)-binom(i,((i+1)/2)+1))))
ndetI = max(1,nint((binom(i,(i+1)/2))))
allocate(tempBuffer(bfIcfg,ndetI))
call getCSFtoDETTransformationMatrix(Isomo, MS, NBFMax, maxDetDimPerBF, tempBuffer)
DetToCSFTransformationMatrix(i,1:bfIcfg,1:ndetI) = tempBuffer(1:bfIcfg,1:ndetI)
deallocate(tempBuffer)
enddo
integer s, bfIcfg
integer countcsf
countcsf = 0
integer countdet
countdet = 0
integer istate
istate = 1
psi_csf_to_config_data(1) = 1
phasedet = 1.0d0
call omp_set_max_active_levels(1)
!$OMP PARALLEL
!$OMP MASTER
do i = 1,N_configuration
startdet = psi_configuration_to_psi_det(1,i)
enddet = psi_configuration_to_psi_det(2,i)
ndetI = enddet-startdet+1
allocate(tempCoeff(ndetI))
countdet = 1
do j = startdet, enddet
idx = psi_configuration_to_psi_det_data(j)
Ialpha(:) = psi_det(:,1,idx)
Ibeta(:) = psi_det(:,2,idx)
call get_phase_qp_to_cfg(Ialpha, Ibeta, phasedet)
tempCoeff(countdet) = psi_coef(idx, istate)*phasedet
norm_det1 += tempCoeff(countdet)*tempCoeff(countdet)
countdet += 1
enddo
!print *,"dimcoef=",bfIcfg,norm_det1
!call printMatrix(tempCoeff,ndetI,1)
s = 0
do k=1,N_int
if (psi_configuration(k,1,i) == 0_bit_kind) cycle
s = s + popcnt(psi_configuration(k,1,i))
enddo
bfIcfg = max(1,nint((binom(s,(s+1)/2)-binom(s,((s+1)/2)+1))))
! perhaps blocking with CFGs of same seniority
! can be more efficient
allocate(tempBuffer(bfIcfg,ndetI))
tempBuffer = DetToCSFTransformationMatrix(s,:bfIcfg,:ndetI)
call dgemm('N','N', bfIcfg, 1, ndetI, 1.d0, tempBuffer, size(tempBuffer,1), tempCoeff, size(tempCoeff,1), 0.d0, psi_coef_config(countcsf+1,1), size(psi_coef_config,1))
!call dgemv('N', NBFMax, maxDetDimPerBF, 1.d0, tempBuffer, size(tempBuffer,1), tempCoeff, 1, 0.d0, psi_coef_config(countcsf), 1)
deallocate(tempCoeff)
deallocate(tempBuffer)
psi_config_data(i,1) = countcsf + 1
do k=1,bfIcfg
psi_csf_to_config_data(countcsf+k) = i
enddo
countcsf += bfIcfg
psi_config_data(i,2) = countcsf
enddo
print *,"Norm det=",norm_det1, size(psi_coef_config,1), " Dim csf=", countcsf
!$OMP END MASTER
!$OMP END PARALLEL
call omp_set_max_active_levels(4)
END_PROVIDER
BEGIN_PROVIDER [ integer, AIJpqMatrixDimsList, (NSOMOMin:NSOMOMax,4,NSOMOMax+1,NSOMOMax+1,2)] BEGIN_PROVIDER [ integer, AIJpqMatrixDimsList, (NSOMOMin:NSOMOMax,4,NSOMOMax+1,NSOMOMax+1,2)]
&BEGIN_PROVIDER [ integer, rowsmax] &BEGIN_PROVIDER [ integer, rowsmax]
&BEGIN_PROVIDER [ integer, colsmax] &BEGIN_PROVIDER [ integer, colsmax]
@ -781,117 +892,6 @@ subroutine calculate_preconditioner_cfg(diag_energies)
end subroutine calculate_preconditioner_cfg end subroutine calculate_preconditioner_cfg
BEGIN_PROVIDER [ real*8, DetToCSFTransformationMatrix, (0:NSOMOMax,NBFMax,maxDetDimPerBF)]
&BEGIN_PROVIDER [ real*8, psi_coef_config, (n_CSF,1)]
&BEGIN_PROVIDER [ integer, psi_config_data, (N_configuration,2)]
&BEGIN_PROVIDER [ integer, psi_csf_to_config_data, (n_CSF)]
use cfunctions
implicit none
BEGIN_DOC
! Documentation for DetToCSFTransformationMatrix
! Provides the matrix of transformatons for the
! conversion between determinant to CSF basis (in BFs)
END_DOC
integer*8 :: Isomo, Idomo
integer(bit_kind) :: Ialpha(N_int),Ibeta(N_int)
integer :: rows, cols, i, j, k
integer :: startdet, enddet, idx
integer*8 MS
integer ndetI
integer :: getNSOMO
real*8,dimension(:,:),allocatable :: tempBuffer
real*8,dimension(:),allocatable :: tempCoeff
real*8 :: norm_det1, phasedet
integer :: nt
norm_det1 = 0.d0
MS = elec_alpha_num - elec_beta_num
print *,"Maxbfdim=",NBFMax
print *,"Maxdetdim=",maxDetDimPerBF
print *,"n_CSF=",n_CSF
print *,"N_configurations=",N_configuration
print *,"n_core_orb=",n_core_orb
! initialization
psi_coef_config = 0.d0
DetToCSFTransformationMatrix(0,:,:) = 1.d0
do i = 2-iand(elec_alpha_num-elec_beta_num,1), NSOMOMax,2
Isomo = IBSET(0_8, i) - 1_8
! rows = Ncsfs
! cols = Ndets
bfIcfg = max(1,nint((binom(i,(i+1)/2)-binom(i,((i+1)/2)+1))))
ndetI = max(1,nint((binom(i,(i+1)/2))))
allocate(tempBuffer(bfIcfg,ndetI))
call getCSFtoDETTransformationMatrix(Isomo, MS, NBFMax, maxDetDimPerBF, tempBuffer)
DetToCSFTransformationMatrix(i,1:bfIcfg,1:ndetI) = tempBuffer(1:bfIcfg,1:ndetI)
deallocate(tempBuffer)
enddo
integer s, bfIcfg
integer countcsf
countcsf = 0
integer countdet
countdet = 0
integer istate
istate = 1
psi_csf_to_config_data(1) = 1
phasedet = 1.0d0
call omp_set_max_active_levels(1)
!$OMP PARALLEL
!$OMP MASTER
do i = 1,N_configuration
startdet = psi_configuration_to_psi_det(1,i)
enddet = psi_configuration_to_psi_det(2,i)
ndetI = enddet-startdet+1
allocate(tempCoeff(ndetI))
countdet = 1
do j = startdet, enddet
idx = psi_configuration_to_psi_det_data(j)
Ialpha(:) = psi_det(:,1,idx)
Ibeta(:) = psi_det(:,2,idx)
call get_phase_qp_to_cfg(Ialpha, Ibeta, phasedet)
tempCoeff(countdet) = psi_coef(idx, istate)*phasedet
norm_det1 += tempCoeff(countdet)*tempCoeff(countdet)
countdet += 1
enddo
!print *,"dimcoef=",bfIcfg,norm_det1
!call printMatrix(tempCoeff,ndetI,1)
s = 0
do k=1,N_int
if (psi_configuration(k,1,i) == 0_bit_kind) cycle
s = s + popcnt(psi_configuration(k,1,i))
enddo
bfIcfg = max(1,nint((binom(s,(s+1)/2)-binom(s,((s+1)/2)+1))))
! perhaps blocking with CFGs of same seniority
! can be more efficient
allocate(tempBuffer(bfIcfg,ndetI))
tempBuffer = DetToCSFTransformationMatrix(s,:bfIcfg,:ndetI)
call dgemm('N','N', bfIcfg, 1, ndetI, 1.d0, tempBuffer, size(tempBuffer,1), tempCoeff, size(tempCoeff,1), 0.d0, psi_coef_config(countcsf+1,1), size(psi_coef_config,1))
!call dgemv('N', NBFMax, maxDetDimPerBF, 1.d0, tempBuffer, size(tempBuffer,1), tempCoeff, 1, 0.d0, psi_coef_config(countcsf), 1)
deallocate(tempCoeff)
deallocate(tempBuffer)
psi_config_data(i,1) = countcsf + 1
do k=1,bfIcfg
psi_csf_to_config_data(countcsf+k) = i
enddo
countcsf += bfIcfg
psi_config_data(i,2) = countcsf
enddo
print *,"Norm det=",norm_det1, size(psi_coef_config,1), " Dim csf=", countcsf
!$OMP END MASTER
!$OMP END PARALLEL
call omp_set_max_active_levels(4)
END_PROVIDER
subroutine obtain_connected_I_foralpha_fromfilterdlist(idxI, nconnectedJ, idslistconnectedJ, listconnectedJ, Ialpha, connectedI, idxs_connectedI, nconnectedI, excitationIds, excitationTypes, diagfactors) subroutine obtain_connected_I_foralpha_fromfilterdlist(idxI, nconnectedJ, idslistconnectedJ, listconnectedJ, Ialpha, connectedI, idxs_connectedI, nconnectedI, excitationIds, excitationTypes, diagfactors)
implicit none implicit none
use bitmasks use bitmasks
@ -1348,10 +1348,10 @@ subroutine calculate_sigma_vector_cfg_nst_naive_store(psi_out, psi_in, n_st, sze
!$OMP AIJpqMatrixDimsList, diag_energies, n_CSF, lock, NBFmax,nconnectedtotalmax, nconnectedmaxJ,maxnalphas,& !$OMP AIJpqMatrixDimsList, diag_energies, n_CSF, lock, NBFmax,nconnectedtotalmax, nconnectedmaxJ,maxnalphas,&
!$OMP num_threads_max) !$OMP num_threads_max)
allocate(singlesI(N_INT,2,max(sze,100))) allocate(singlesI(N_INT,2,max(sze,1000)))
allocate(idxs_singlesI(max(sze,100))) allocate(idxs_singlesI(max(sze,1000)))
allocate(excitationIds_single(2,max(sze,100))) allocate(excitationIds_single(2,max(sze,1000)))
allocate(excitationTypes_single(max(sze,100))) allocate(excitationTypes_single(max(sze,1000)))
! !
!!! Single Excitations !!! !!! Single Excitations !!!
@ -1488,11 +1488,8 @@ subroutine calculate_sigma_vector_cfg_nst_naive_store(psi_out, psi_in, n_st, sze
allocate(CCmattmp(n_st,NBFmax)) allocate(CCmattmp(n_st,NBFmax))
! Loop over all selected configurations ! Loop over all selected configurations
!$OMP DO SCHEDULE(dynamic,16) !$OMP DO SCHEDULE(static)
do i = istart_cfg,iend_cfg do i = istart_cfg,iend_cfg
! TKI = 0.d0
! GIJpqrs = 0.d0
! TKIGIJ = 0.d0
! if Seniority_range > 8 then ! if Seniority_range > 8 then
! continue ! continue
@ -1552,6 +1549,9 @@ subroutine calculate_sigma_vector_cfg_nst_naive_store(psi_out, psi_in, n_st, sze
allocate(GIJpqrs(totcolsTKI,nconnectedI)) ! gpqrs allocate(GIJpqrs(totcolsTKI,nconnectedI)) ! gpqrs
allocate(TKIGIJ(n_st,rowsTKI,nconnectedI)) ! TKI * gpqrs allocate(TKIGIJ(n_st,rowsTKI,nconnectedI)) ! TKI * gpqrs
!print *,"\t---rowsTKI=",rowsTKI," totCols=",totcolsTKI !print *,"\t---rowsTKI=",rowsTKI," totCols=",totcolsTKI
TKI = 0.d0
GIJpqrs = 0.d0
TKIGIJ = 0.d0
totcolsTKI = 0 totcolsTKI = 0
do j = 1,nconnectedI do j = 1,nconnectedI