From 14a0ff375ae4426be2127e8341ad7f483160ca47 Mon Sep 17 00:00:00 2001 From: v1j4y Date: Tue, 14 Jun 2022 12:51:35 +0200 Subject: [PATCH] Added missing file. --- src/csf/obtain_I_foralpha.irp.f | 383 ++++++++++++++++++++++++++++++++ 1 file changed, 383 insertions(+) create mode 100644 src/csf/obtain_I_foralpha.irp.f diff --git a/src/csf/obtain_I_foralpha.irp.f b/src/csf/obtain_I_foralpha.irp.f new file mode 100644 index 00000000..d6fbe386 --- /dev/null +++ b/src/csf/obtain_I_foralpha.irp.f @@ -0,0 +1,383 @@ +subroutine obtain_connected_J_givenI(idxI, givenI, connectedI, idxs_connectedI, nconnectedI,ntotalconnectedI) + implicit none + use bitmasks + BEGIN_DOC + ! Documentation for obtain_connected_I_foralpha + ! This function returns all those selected configurations + ! which are connected to the input configuration + ! givenI by a single excitation. + ! + ! The type of excitations are ordered as follows: + ! Type 1 - SOMO -> SOMO + ! Type 2 - DOMO -> VMO + ! Type 3 - SOMO -> VMO + ! Type 4 - DOMO -> SOMO + ! + ! Order of operators + ! \alpha> = a^\dag_p a_q |I> = E_pq |I> + END_DOC + integer ,intent(in) :: idxI + integer(bit_kind),intent(in) :: givenI(N_int,2) + integer(bit_kind),intent(out) :: connectedI(N_int,2,*) + integer ,intent(out) :: idxs_connectedI(*) + integer,intent(out) :: nconnectedI + integer,intent(out) :: ntotalconnectedI + integer*8 :: Idomo + integer*8 :: Isomo + integer*8 :: Jdomo + integer*8 :: Jsomo + integer*8 :: IJsomo + integer*8 :: diffSOMO + integer*8 :: diffDOMO + integer*8 :: xordiffSOMODOMO + integer :: ndiffSOMO + integer :: ndiffDOMO + integer :: nxordiffSOMODOMO + integer :: ii,i,j,k,l,p,q,nsomoJ,nsomoalpha,starti,endi,extyp,nholes + integer :: listholes(mo_num) + integer :: holetype(mo_num) + integer :: end_index + integer :: Nsomo_I + + ! + ! 2 2 1 1 0 0 : 1 1 0 0 0 0 + ! 0 0 1 1 0 0 + ! + ! 2 1 1 1 1 0 : 1 0 0 0 0 0 + ! 0 1 1 1 1 0 + !xorS 0 1 0 0 1 0 : 2 + !xorD 0 1 0 0 0 0 : 1 + !xorSD 0 0 0 0 1 0 : 1 + ! ----- + ! 4 + ! 1 1 1 1 1 1 : 0 0 0 0 0 0 + ! 1 1 1 1 1 1 + ! 1 1 0 0 1 1 : 4 + ! 1 1 0 0 0 0 : 2 + ! 0 0 0 0 1 1 : 2 + ! ----- + ! 8 + ! + + nconnectedI = 0 + ntotalconnectedI = 0 + end_index = N_configuration + + ! Since CFGs are sorted wrt to seniority + ! we don't have to search the full CFG list + Isomo = givenI(1,1) + Idomo = givenI(1,2) + Nsomo_I = POPCNT(Isomo) + end_index = min(N_configuration,cfg_seniority_index(min(Nsomo_I+6,elec_num))-1) + if(end_index .LT. 0) end_index= N_configuration + !end_index = N_configuration + !print *,"Start and End = ",idxI, end_index + + + p = 0 + q = 0 + do i=idxI,end_index + !if(.True.) then + ! nconnectedI += 1 + ! connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + ! idxs_connectedI(nconnectedI)=i + ! cycle + !endif + Isomo = givenI(1,1) + Idomo = givenI(1,2) + Jsomo = psi_configuration(1,1,i) + Jdomo = psi_configuration(1,2,i) + diffSOMO = IEOR(Isomo,Jsomo) + ndiffSOMO = POPCNT(diffSOMO) + diffDOMO = IEOR(Idomo,Jdomo) + xordiffSOMODOMO = IEOR(diffSOMO,diffDOMO) + ndiffDOMO = POPCNT(diffDOMO) + nxordiffSOMODOMO = POPCNT(xordiffSOMODOMO) + nxordiffSOMODOMO += ndiffSOMO + ndiffDOMO + if((nxordiffSOMODOMO .EQ. 4) .AND. ndiffSOMO .EQ. 2) then + !------- + ! MONO | + !------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=i + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)) + else if((nxordiffSOMODOMO .EQ. 8) .AND. ndiffSOMO .EQ. 4) then + !---------------------------- + ! DOMO -> VMO + DOMO -> VMO | + !---------------------------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=i + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)) + else if((nxordiffSOMODOMO .EQ. 6) .AND. ndiffSOMO .EQ. 2) then + !---------------------------- + ! DOUBLE + !---------------------------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=i + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)) + else if((nxordiffSOMODOMO .EQ. 2) .AND. ndiffSOMO .EQ. 3) then + !----------------- + ! DOUBLE + !----------------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=i + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)) + else if((nxordiffSOMODOMO .EQ. 4) .AND. ndiffSOMO .EQ. 0) then + !----------------- + ! DOUBLE + !----------------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=i + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)) + else if((ndiffSOMO + ndiffDOMO) .EQ. 0) then + !-------- + ! I = I | + !-------- + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)= i + ! find out all pq holes possible + nholes = 0 + ! holes in SOMO + Isomo = psi_configuration(1,1,i) + Idomo = psi_configuration(1,2,i) + do ii = 1,mo_num + if(POPCNT(IAND(Isomo,IBSET(0_8,ii-1))) .EQ. 1) then + nholes += 1 + listholes(nholes) = ii + holetype(nholes) = 1 + endif + end do + ! holes in DOMO + do ii = 1,mo_num + if(POPCNT(IAND(Idomo,IBSET(0_8,ii-1))) .EQ. 1) then + nholes += 1 + listholes(nholes) = ii + holetype(nholes) = 2 + endif + end do + ntotalconnectedI += max(1,(psi_config_data(i,2)-psi_config_data(i,1)+1)*nholes) + endif + end do + +end subroutine obtain_connected_J_givenI + +subroutine obtain_connected_I_foralpha(idxI, Ialpha, connectedI, idxs_connectedI, nconnectedI, excitationIds, excitationTypes, diagfactors) + implicit none + use bitmasks + BEGIN_DOC + ! Documentation for obtain_connected_I_foralpha + ! This function returns all those selected configurations + ! which are connected to the input configuration + ! Ialpha by a single excitation. + ! + ! The type of excitations are ordered as follows: + ! Type 1 - SOMO -> SOMO + ! Type 2 - DOMO -> VMO + ! Type 3 - SOMO -> VMO + ! Type 4 - DOMO -> SOMO + ! + ! Order of operators + ! \alpha> = a^\dag_p a_q |I> = E_pq |I> + END_DOC + integer ,intent(in) :: idxI + integer(bit_kind),intent(in) :: Ialpha(N_int,2) + integer(bit_kind),intent(out) :: connectedI(N_int,2,*) + integer ,intent(out) :: idxs_connectedI(*) + integer,intent(out) :: nconnectedI + integer,intent(out) :: excitationIds(2,*) + integer,intent(out) :: excitationTypes(*) + real*8 ,intent(out) :: diagfactors(*) + integer*8 :: Idomo + integer*8 :: Isomo + integer*8 :: Jdomo + integer*8 :: Jsomo + integer*8 :: IJsomo + integer*8 :: diffSOMO + integer*8 :: diffDOMO + integer*8 :: xordiffSOMODOMO + integer :: ndiffSOMO + integer :: ndiffDOMO + integer :: nxordiffSOMODOMO + integer :: ii,i,j,k,l,p,q,nsomoJ,nsomoalpha,starti,endi,extyp,nholes + integer :: listholes(mo_num) + integer :: holetype(mo_num) + integer :: end_index + integer :: Nsomo_alpha + + nconnectedI = 0 + end_index = N_configuration + + ! Since CFGs are sorted wrt to seniority + ! we don't have to search the full CFG list + Isomo = Ialpha(1,1) + Idomo = Ialpha(1,2) + Nsomo_alpha = POPCNT(Isomo) + end_index = min(N_configuration,cfg_seniority_index(min(Nsomo_alpha+4,elec_num))-1) + if(end_index .LT. 0) end_index= N_configuration + end_index = N_configuration + + + p = 0 + q = 0 + if (N_int > 1) stop 'obtain_connected_i_foralpha : N_int > 1' + do i=idxI,end_index + Isomo = Ialpha(1,1) + Idomo = Ialpha(1,2) + Jsomo = psi_configuration(1,1,i) + Jdomo = psi_configuration(1,2,i) + diffSOMO = IEOR(Isomo,Jsomo) + ndiffSOMO = POPCNT(diffSOMO) + !if(idxI.eq.1)then + ! print *," \t idxI=",i," diffS=",ndiffSOMO," popJs=", POPCNT(Jsomo)," popIs=",POPCNT(Isomo) + !endif + diffDOMO = IEOR(Idomo,Jdomo) + xordiffSOMODOMO = IEOR(diffSOMO,diffDOMO) + ndiffDOMO = POPCNT(diffDOMO) + nxordiffSOMODOMO = POPCNT(xordiffSOMODOMO) + nxordiffSOMODOMO += ndiffSOMO + ndiffDOMO + if((nxordiffSOMODOMO .EQ. 4) .AND. ndiffSOMO .EQ. 2) then + select case(ndiffDOMO) + case (0) + ! SOMO -> VMO + !print *,"obt SOMO -> VMO" + extyp = 3 + IJsomo = IEOR(Isomo, Jsomo) +!IRP_IF WITHOUT_TRAILZ +! p = (popcnt(ieor( IAND(Isomo,IJsomo) , IAND(Isomo,IJsomo) -1))-1) + 1 +!IRP_ELSE + p = TRAILZ(IAND(Isomo,IJsomo)) + 1 +!IRP_ENDIF + IJsomo = IBCLR(IJsomo,p-1) +!IRP_IF WITHOUT_TRAILZ +! q = (popcnt(ieor(IJsomo,IJsomo-1))-1) + 1 +!IRP_ELSE + q = TRAILZ(IJsomo) + 1 +!IRP_ENDIF + case (1) + ! DOMO -> VMO + ! or + ! SOMO -> SOMO + nsomoJ = POPCNT(Jsomo) + nsomoalpha = POPCNT(Isomo) + if(nsomoJ .GT. nsomoalpha) then + ! DOMO -> VMO + !print *,"obt DOMO -> VMO" + extyp = 2 +!IRP_IF WITHOUT_TRAILZ +! p = (popcnt(ieor( IEOR(Idomo,Jdomo),IEOR(Idomo,Jdomo) -1))-1) + 1 +!IRP_ELSE + p = TRAILZ(IEOR(Idomo,Jdomo)) + 1 +!IRP_ENDIF + Isomo = IEOR(Isomo, Jsomo) + Isomo = IBCLR(Isomo,p-1) +!IRP_IF WITHOUT_TRAILZ +! q = (popcnt(ieor(Isomo,Isomo-1))-1) + 1 +!IRP_ELSE + q = TRAILZ(Isomo) + 1 +!IRP_ENDIF + else + ! SOMO -> SOMO + !print *,"obt SOMO -> SOMO" + extyp = 1 +!IRP_IF WITHOUT_TRAILZ +! q = (popcnt(ieor( IEOR(Idomo,Jdomo), IEOR(Idomo,Jdomo)-1))-1) + 1 +!IRP_ELSE + q = TRAILZ(IEOR(Idomo,Jdomo)) + 1 +!IRP_ENDIF + Isomo = IEOR(Isomo, Jsomo) + Isomo = IBCLR(Isomo,q-1) +!IRP_IF WITHOUT_TRAILZ +! p = (popcnt(ieor(Isomo,Isomo-1))-1) + 1 +!IRP_ELSE + p = TRAILZ(Isomo) + 1 +!IRP_ENDIF + end if + case (2) + ! DOMO -> SOMO + !print *,"obt DOMO -> SOMO" + extyp = 4 + IJsomo = IEOR(Isomo, Jsomo) +!IRP_IF WITHOUT_TRAILZ +! p = (popcnt(ieor( IAND(Jsomo,IJsomo), IAND(Jsomo,IJsomo)-1))-1) + 1 +!IRP_ELSE + p = TRAILZ(IAND(Jsomo,IJsomo)) + 1 +!IRP_ENDIF + IJsomo = IBCLR(IJsomo,p-1) +!IRP_IF WITHOUT_TRAILZ +! q = (popcnt(ieor( IJsomo , IJsomo -1))-1) + 1 +!IRP_ELSE + q = TRAILZ(IJsomo) + 1 +!IRP_ENDIF + case default + print *,"something went wront in get connectedI" + end select + starti = psi_config_data(i,1) + endi = psi_config_data(i,2) + nconnectedI += 1 + do k=1,N_int + connectedI(k,1,nconnectedI) = psi_configuration(k,1,i) + connectedI(k,2,nconnectedI) = psi_configuration(k,2,i) + enddo + idxs_connectedI(nconnectedI)=starti + excitationIds(1,nconnectedI)=p + excitationIds(2,nconnectedI)=q + excitationTypes(nconnectedI) = extyp + diagfactors(nconnectedI) = 1.0d0 + else if((ndiffSOMO + ndiffDOMO) .EQ. 0) then + ! find out all pq holes possible + nholes = 0 + ! holes in SOMO + Isomo = psi_configuration(1,1,i) + Idomo = psi_configuration(1,2,i) + do ii = 1,mo_num + if(POPCNT(IAND(Isomo,IBSET(0_8,ii-1))) .EQ. 1) then + nholes += 1 + listholes(nholes) = ii + holetype(nholes) = 1 + endif + end do + ! holes in DOMO + do ii = 1,mo_num + if(POPCNT(IAND(Idomo,IBSET(0_8,ii-1))) .EQ. 1) then + nholes += 1 + listholes(nholes) = ii + holetype(nholes) = 2 + endif + end do + + do k=1,nholes + p = listholes(k) + q = p + extyp = 1 + if(holetype(k) .EQ. 1) then + starti = psi_config_data(i,1) + endi = psi_config_data(i,2) + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=starti + excitationIds(1,nconnectedI)=p + excitationIds(2,nconnectedI)=q + excitationTypes(nconnectedI) = extyp + diagfactors(nconnectedI) = 1.0d0 + else + starti = psi_config_data(i,1) + endi = psi_config_data(i,2) + nconnectedI += 1 + connectedI(:,:,nconnectedI) = psi_configuration(:,:,i) + idxs_connectedI(nconnectedI)=starti + excitationIds(1,nconnectedI)=p + excitationIds(2,nconnectedI)=q + excitationTypes(nconnectedI) = extyp + diagfactors(nconnectedI) = 2.0d0 + endif + enddo + endif + end do + +end subroutine obtain_connected_I_foralpha