From 05c88a79ba1d8ade596f7a4acd916db5b6cd2868 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Thu, 23 Mar 2017 15:13:11 +0100 Subject: [PATCH] MRCC optimizations --- src/Determinants/filter_connected.irp.f | 16 +++++++++---- src/Determinants/slater_rules.irp.f | 31 +++++++++++++++---------- 2 files changed, 31 insertions(+), 16 deletions(-) diff --git a/src/Determinants/filter_connected.irp.f b/src/Determinants/filter_connected.irp.f index b76540f7..ed6aa6d2 100644 --- a/src/Determinants/filter_connected.irp.f +++ b/src/Determinants/filter_connected.irp.f @@ -285,7 +285,9 @@ subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_micro do i=1, N_minilist - do j=1,Nint + mobileMask(1,1) = iand(key_mask_neg(1,1), minilist(1,1,i)) + mobileMask(1,2) = iand(key_mask_neg(1,2), minilist(1,2,i)) + do j=2,Nint mobileMask(j,1) = iand(key_mask_neg(j,1), minilist(j,1,i)) mobileMask(j,2) = iand(key_mask_neg(j,2), minilist(j,2,i)) end do @@ -296,7 +298,9 @@ subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_micro if(n_element(1) + n_element(2) /= 4) then idx_microlist(cur_microlist(0)) = i - do k=1,Nint + microlist(1,1,cur_microlist(0)) = minilist(1,1,i) + microlist(1,2,cur_microlist(0)) = minilist(1,2,i) + do k=2,Nint microlist(k,1,cur_microlist(0)) = minilist(k,1,i) microlist(k,2,cur_microlist(0)) = minilist(k,2,i) enddo @@ -305,8 +309,10 @@ subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_micro do j=1,n_element(1) nt = list(j,1) idx_microlist(cur_microlist(nt)) = i + microlist(1,1,cur_microlist(nt)) = minilist(1,1,i) + microlist(1,2,cur_microlist(nt)) = minilist(1,2,i) ! TODO : Page faults - do k=1,Nint + do k=2,Nint microlist(k,1,cur_microlist(nt)) = minilist(k,1,i) microlist(k,2,cur_microlist(nt)) = minilist(k,2,i) enddo @@ -316,7 +322,9 @@ subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_micro do j=1,n_element(2) nt = list(j,2) + mo_tot_num idx_microlist(cur_microlist(nt)) = i - do k=1,Nint + microlist(1,1,cur_microlist(nt)) = minilist(1,1,i) + microlist(1,2,cur_microlist(nt)) = minilist(1,2,i) + do k=2,Nint microlist(k,1,cur_microlist(nt)) = minilist(k,1,i) microlist(k,2,cur_microlist(nt)) = minilist(k,2,i) enddo diff --git a/src/Determinants/slater_rules.irp.f b/src/Determinants/slater_rules.irp.f index 7556e2b9..eebf9611 100644 --- a/src/Determinants/slater_rules.irp.f +++ b/src/Determinants/slater_rules.irp.f @@ -925,22 +925,29 @@ subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullLis N_miniList = 0 + integer :: e_ab + e_ab = n_a+n_b do i=1,N_fullList - e_a = n_a - popcnt(iand(fullList(1, 1, i), key_mask(1, 1))) - e_b = n_b - popcnt(iand(fullList(1, 2, i), key_mask(1, 2))) + e_a = e_ab - popcnt(iand(fullList(1, 1, i), key_mask(1, 1))) & + - popcnt(iand(fullList(1, 2, i), key_mask(1, 2))) do ni=2,nint - e_a -= popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1))) - e_b -= popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2))) + e_a = e_a - popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1))) & + - popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2))) end do - if(e_a + e_b <= 2) then - N_miniList = N_miniList + 1 - do ni=1,Nint - miniList(ni,1,N_miniList) = fullList(ni,1,i) - miniList(ni,2,N_miniList) = fullList(ni,2,i) - enddo - idx_miniList(N_miniList) = i - end if + if(e_a > 2) then + cycle + endif + + N_miniList = N_miniList + 1 + miniList(1,1,N_miniList) = fullList(1,1,i) + miniList(1,2,N_miniList) = fullList(1,2,i) + do ni=2,Nint + miniList(ni,1,N_miniList) = fullList(ni,1,i) + miniList(ni,2,N_miniList) = fullList(ni,2,i) + enddo + idx_miniList(N_miniList) = i + end do end subroutine