From 3815c4ef35c655fa71ad8877823115abd180d7aa Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Mon, 7 Jun 2021 23:41:37 +0200 Subject: [PATCH] Swaps in determinants --- make.config.ifort | 2 +- src/det.irp.f | 248 ++++++++++++++++++++++++++++++++++++++-------- 2 files changed, 209 insertions(+), 41 deletions(-) diff --git a/make.config.ifort b/make.config.ifort index d241749..7321337 100644 --- a/make.config.ifort +++ b/make.config.ifort @@ -10,7 +10,7 @@ ALIGN="32" ## FORTRAN compiler FC="ifort" -FCFLAGS="-O2 -g -ip -ftz -finline ${CPU_TYPE}" #-traceback +FCFLAGS="-O2 -g -ip -ftz -finline ${CPU_TYPE} -qopenmp-simd" #-traceback LIB="-mkl=sequential" ## IRPF90 diff --git a/src/det.irp.f b/src/det.irp.f index 91e2034..6cfd595 100644 --- a/src/det.irp.f +++ b/src/det.irp.f @@ -1,18 +1,22 @@ -BEGIN_PROVIDER [ integer, det_i ] + BEGIN_PROVIDER [ integer, det_i ] +&BEGIN_PROVIDER [ integer, det_i_prev ] BEGIN_DOC ! Current running alpha determinant END_DOC det_i=det_alpha_order(1) + det_i_prev=det_alpha_order(1) END_PROVIDER -BEGIN_PROVIDER [ integer, det_j ] + BEGIN_PROVIDER [ integer, det_j ] +&BEGIN_PROVIDER [ integer, det_j_prev ] BEGIN_DOC ! Current running beta determinant END_DOC det_j=det_beta_order(1) + det_j_prev=det_beta_order(1) END_PROVIDER @@ -323,7 +327,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) double precision :: zj, zj1, zj2, zj3 !DIR$ NOPREFETCH - !OMP$ SIMD NOVECREMAINDER do i=1,$n u(i) = m(i) - S(i,l) enddo @@ -331,7 +334,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) NOVECREMAINDER do i=1,$n-1,4 zj = zj + S_inv(i,l)*u(i) + S_inv(i+1,l)*u(i+1) & + S_inv(i+2,l)*u(i+2) + S_inv(i+3,l)*u(i+3) @@ -353,7 +355,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj3 = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2,zj3) NOVECREMAINDER do i=1,$n zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -367,7 +368,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) enddo !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) NOVECREMAINDER do i=1,$n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) @@ -380,7 +380,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj3 = z(i+3) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2,zj3) NOVECREMAINDER do j=1,$n S_inv(j,i ) = S_inv(j,i )*lambda - w(j)*zj S_inv(j,i+1) = S_inv(j,i+1)*lambda - w(j)*zj1 @@ -459,7 +458,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) do i=1,$n-1,4 zj = zj + S_inv(i,l)*u(i) + S_inv(i+1,l)*u(i+1) & + S_inv(i+2,l)*u(i+2) + S_inv(i+3,l)*u(i+3) @@ -482,7 +480,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj3 = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2,zj3) NOVECREMAINDER do i=1,$n-1 zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -498,14 +495,12 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) NOVECREMAINDER do i=1,$n-1 zj = zj + S_inv(i,$n)*u(i) enddo z($n) = zj + S_inv($n,$n)*u($n) !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) NOVECREMAINDER do i=1,$n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) @@ -518,7 +513,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj3 = z(i+3) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2,zj3) NOVECREMAINDER do j=1,$n-1 S_inv(j,i ) = S_inv(j,i )*lambda - w(j)*zj S_inv(j,i+1) = S_inv(j,i+1)*lambda - w(j)*zj1 @@ -534,7 +528,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = z($n) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj) NOVECREMAINDER do i=1,$n S_inv(i,$n) = S_inv(i,$n)*lambda -w(i)*zj enddo @@ -607,7 +600,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) double precision :: zj, zj1, zj2, zj3 !DIR$ NOPREFETCH - !OMP$ SIMD NOVECREMAINDER do i=1,$n u(i) = m(i) - S(i,l) enddo @@ -615,7 +607,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) NOVECREMAINDER do i=1,$n-2,4 zj = zj + S_inv(i,l)*u(i) + S_inv(i+1,l)*u(i+1) & + S_inv(i+2,l)*u(i+2) + S_inv(i+3,l)*u(i+3) @@ -638,7 +629,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj2 = 0.d0 zj3 = 0.d0 !DIR$ VECTOR ALIGNED - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2,zj3) NOVECREMAINDER do i=1,$n-2 zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -660,7 +650,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj1 = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj,zj1) NOVECREMAINDER do i=1,$n-2 zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -671,7 +660,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) z(j+1) = z(j+1) + S_inv($n,j+1)*u($n) !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) NOVECREMAINDER do i=1,$n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) @@ -683,7 +671,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj2 = z(i+2) zj3 = z(i+3) !DIR$ VECTOR ALIGNED - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2,zj3) NOVECREMAINDER do j=1,$n-2 S_inv(j,i ) = S_inv(j,i )*lambda -zj *w(j) S_inv(j,i+1) = S_inv(j,i+1)*lambda -zj1*w(j) @@ -704,7 +691,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = z(i) zj1= z(i+1) !DIR$ VECTOR ALIGNED - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1) do j=1,$n-2 S_inv(j,i ) = S_inv(j,i )*lambda -zj*w(j) S_inv(j,i+1) = S_inv(j,i+1)*lambda -zj1*w(j) @@ -780,7 +766,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) double precision :: zj, zj1, zj2, zj3 - !OMP$ SIMD do i=1,$n u(i) = m(i) - S(i,l) enddo @@ -788,7 +773,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) NOVECREMAINDER do i=1,$n-3,4 zj = zj + S_inv(i,l)*u(i) + S_inv(i+1,l)*u(i+1) & + S_inv(i+2,l)*u(i+2) + S_inv(i+3,l)*u(i+3) @@ -812,7 +796,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj2 = 0.d0 zj3 = 0.d0 !DIR$ VECTOR ALIGNED - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2,zj3) do i=1,$n-3 zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -839,7 +822,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj2 = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2) do i=1,$n-3 zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -856,7 +838,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) z(j+2) = z(j+2) + S_inv($n,j+2)*u($n) !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) do i=1,$n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) @@ -869,7 +850,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj3 = z(i+3) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2,zj3) do j=1,$n-3 S_inv(j,i ) = S_inv(j,i )*lambda - w(j)*zj S_inv(j,i+1) = S_inv(j,i+1)*lambda - w(j)*zj1 @@ -896,7 +876,6 @@ subroutine det_update$n(n,LDS,m,l,S,S_inv,d) zj2 = z(i+2) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2) do j=1,$n S_inv(j,i ) = S_inv(j,i )*lambda - w(j)*zj S_inv(j,i+1) = S_inv(j,i+1)*lambda - w(j)*zj1 @@ -980,7 +959,6 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) zl = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zl) do i=1,n zl = zl + S_inv(i,l)*u(i) enddo @@ -988,7 +966,7 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) d_inv = 1.d0/d d = d+zl lambda = d*d_inv - +! if ( dabs(lambda) < 1.d-3 ) then d = 0.d0 return @@ -1006,7 +984,6 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) zj3 = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj,zj1,zj2,zj3) do i=1,n zj = zj + S_inv(i,j )*u(i) zj1 = zj1 + S_inv(i,j+1)*u(i) @@ -1023,7 +1000,6 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) zj = 0.d0 !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD REDUCTION(+:zj) do i=1,n zj = zj + S_inv(i,j)*u(i) enddo @@ -1031,14 +1007,12 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) enddo !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) do i=1,n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) enddo !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(d_inv) do i=1,n w(i) = S_inv(i,l)*d_inv S(i,l) = m(i) @@ -1051,7 +1025,6 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) zj3 = z(i+3) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj,zj1,zj2,zj3) do j=1,n S_inv(j,i ) = S_inv(j,i )*lambda -zj *w(j) S_inv(j,i+1) = S_inv(j,i+1)*lambda -zj1*w(j) @@ -1064,7 +1037,6 @@ subroutine det_update_general(n,LDS,m,l,S,S_inv,d) zj = z(i) !DIR$ VECTOR ALIGNED !DIR$ NOPREFETCH - !OMP$ SIMD FIRSTPRIVATE(lambda,zj) do j=1,n S_inv(j,i) = S_inv(j,i)*lambda -zj*w(j) enddo @@ -1100,8 +1072,61 @@ subroutine bitstring_to_list( string, list, n_elements, Nint) enddo end +subroutine get_excitation_degree_spin(key1,key2,degree,Nint) + implicit none + BEGIN_DOC + ! Returns the excitation degree between two determinants. + END_DOC + integer, intent(in) :: Nint + integer(8), intent(in) :: key1(Nint) + integer(8), intent(in) :: key2(Nint) + integer, intent(out) :: degree + + integer(8) :: xorvec(32) + integer :: l + + ASSERT (Nint > 0) + ASSERT (Nint <= 32) + + select case (Nint) + + case (1) + xorvec(1) = xor( key1(1), key2(1)) + degree = popcnt(xorvec(1)) + + case (2) + xorvec(1) = xor( key1(1), key2(1)) + xorvec(2) = xor( key1(2), key2(2)) + degree = popcnt(xorvec(1))+popcnt(xorvec(2)) + + case (3) + xorvec(1) = xor( key1(1), key2(1)) + xorvec(2) = xor( key1(2), key2(2)) + xorvec(3) = xor( key1(3), key2(3)) + degree = sum(popcnt(xorvec(1:3))) + + case (4) + xorvec(1) = xor( key1(1), key2(1)) + xorvec(2) = xor( key1(2), key2(2)) + xorvec(3) = xor( key1(3), key2(3)) + xorvec(4) = xor( key1(4), key2(4)) + degree = sum(popcnt(xorvec(1:4))) + + case default + do l=1,Nint + xorvec(l) = xor( key1(l), key2(l)) + enddo + degree = sum(popcnt(xorvec(1:Nint))) + + end select + + degree = shiftr(degree,1) + +end + BEGIN_PROVIDER [ integer, mo_list_alpha_curr, (elec_alpha_num) ] &BEGIN_PROVIDER [ integer, mo_list_alpha_prev, (elec_alpha_num) ] +&BEGIN_PROVIDER [ integer, mo_exc_alpha_curr ] implicit none BEGIN_DOC ! List of MOs in the current alpha determinant @@ -1117,11 +1142,16 @@ end if (l /= elec_alpha_num) then stop 'error in number of alpha electrons' endif + call get_excitation_degree_spin(psi_det_alpha(1,det_i), & + psi_det_alpha(1,det_i_prev), & + mo_exc_alpha_curr, N_int) END_PROVIDER + BEGIN_PROVIDER [ integer, mo_list_beta_curr, (elec_beta_num) ] &BEGIN_PROVIDER [ integer, mo_list_beta_prev, (elec_beta_num) ] +&BEGIN_PROVIDER [ integer, mo_exc_beta_curr ] implicit none BEGIN_DOC ! List of MOs in the current beta determinant @@ -1141,6 +1171,10 @@ END_PROVIDER if (l /= elec_beta_num) then stop 'error in number of beta electrons' endif + call get_excitation_degree_spin(psi_det_beta(1,det_j), & + psi_det_beta(1,det_j_prev), & + mo_exc_beta_curr, N_int) + END_PROVIDER BEGIN_PROVIDER [ double precision, det_alpha_value_curr ] @@ -1163,7 +1197,10 @@ END_PROVIDER double precision :: ddet integer :: i,j,k,imo,l integer :: to_do(elec_alpha_num), n_to_do_old, n_to_do + double precision :: tmp_inv(elec_alpha_num_8) + real :: tmp_det(elec_alpha_num_8) integer, save :: ifirst + !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: tmp_inv, tmp_det if (ifirst == 0) then ifirst = 1 @@ -1175,7 +1212,6 @@ END_PROVIDER PROVIDE mo_value if (det_i /= det_alpha_order(1) ) then -! if (det_i == -1 ) then n_to_do = 0 do k=1,elec_alpha_num @@ -1186,6 +1222,69 @@ END_PROVIDER endif enddo + ! make swaps and keep 1 update + if (n_to_do > 1 .and. mo_exc_alpha_curr == 1) then + + if (iand(n_to_do+1,1)==1) then + det_alpha_value_curr = -det_alpha_value_curr + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + slater_matrix_alpha_inv_det = - slater_matrix_alpha_inv_det + endif + + if (mo_list_alpha_curr(to_do(1)) == mo_list_alpha_prev(to_do(1)+1)) then + + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + tmp_det(l) = slater_matrix_alpha(l,to_do(1)) + tmp_inv(l) = slater_matrix_alpha_inv_det(l,to_do(1)) + enddo + + do k=to_do(1),to_do(n_to_do-1) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + slater_matrix_alpha(l,k) = slater_matrix_alpha(l,k+1) + slater_matrix_alpha_inv_det(l,k) = slater_matrix_alpha_inv_det(l,k+1) + enddo + enddo + k = to_do(n_to_do) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + slater_matrix_alpha(l,k) = tmp_det(l) + slater_matrix_alpha_inv_det(l,k) = tmp_inv(l) + enddo + to_do(1) = to_do(n_to_do) + + else if (mo_list_alpha_curr(to_do(n_to_do)) == mo_list_alpha_prev(to_do(n_to_do)-1)) then + k = to_do(n_to_do) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + tmp_det(l) = slater_matrix_alpha(l,k) + tmp_inv(l) = slater_matrix_alpha_inv_det(l,k) + enddo + do k=to_do(n_to_do),to_do(2),-1 + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + slater_matrix_alpha(l,k) = slater_matrix_alpha(l,k-1) + slater_matrix_alpha_inv_det(l,k) = slater_matrix_alpha_inv_det(l,k-1) + enddo + enddo + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_alpha_num_8 + slater_matrix_alpha(l,to_do(1)) = tmp_det(l) + slater_matrix_alpha_inv_det(l,to_do(1)) = tmp_inv(l) + enddo + + endif + n_to_do = 1 + endif + ddet = 0.d0 if (n_to_do < shiftl(elec_alpha_num,1)) then @@ -1270,8 +1369,10 @@ END_PROVIDER double precision :: ddet integer :: i,j,k,imo,l integer :: to_do(elec_alpha_num-mo_closed_num), n_to_do_old, n_to_do - + double precision :: tmp_inv(elec_alpha_num_8) + real :: tmp_det(elec_alpha_num_8) integer, save :: ifirst + if (elec_beta_num == 0) then det_beta_value_curr = 0.d0 return @@ -1285,7 +1386,6 @@ END_PROVIDER PROVIDE mo_value if (det_j /= det_beta_order(1)) then -! if (det_j == -1) then n_to_do = 0 do k=mo_closed_num+1,elec_beta_num @@ -1296,7 +1396,71 @@ END_PROVIDER endif enddo + ! make swaps and keep 1 update + if (n_to_do > 1 .and. mo_exc_beta_curr == 1) then + + if (iand(n_to_do+1,1)==1) then + det_beta_value_curr = -det_beta_value_curr + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + slater_matrix_beta_inv_det = - slater_matrix_beta_inv_det + endif + + if (mo_list_beta_curr(to_do(1)) == mo_list_beta_prev(to_do(1)+1)) then + + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + tmp_det(l) = slater_matrix_beta(l,to_do(1)) + tmp_inv(l) = slater_matrix_beta_inv_det(l,to_do(1)) + enddo + + do k=to_do(1),to_do(n_to_do-1) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + slater_matrix_beta(l,k) = slater_matrix_beta(l,k+1) + slater_matrix_beta_inv_det(l,k) = slater_matrix_beta_inv_det(l,k+1) + enddo + enddo + k = to_do(n_to_do) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + slater_matrix_beta(l,k) = tmp_det(l) + slater_matrix_beta_inv_det(l,k) = tmp_inv(l) + enddo + to_do(1) = to_do(n_to_do) + + else if (mo_list_beta_curr(to_do(n_to_do)) == mo_list_beta_prev(to_do(n_to_do)-1)) then + k = to_do(n_to_do) + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + tmp_det(l) = slater_matrix_beta(l,k) + tmp_inv(l) = slater_matrix_beta_inv_det(l,k) + enddo + do k=to_do(n_to_do),to_do(2),-1 + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + slater_matrix_beta(l,k) = slater_matrix_beta(l,k-1) + slater_matrix_beta_inv_det(l,k) = slater_matrix_beta_inv_det(l,k-1) + enddo + enddo + !DIR$ VECTOR ALWAYS + !DIR$ VECTOR ALIGNED + do l=1,elec_beta_num_8 + slater_matrix_beta(l,to_do(1)) = tmp_det(l) + slater_matrix_beta_inv_det(l,to_do(1)) = tmp_inv(l) + enddo + + endif + n_to_do = 1 + endif + ddet = 0.d0 + if (n_to_do < shiftl(elec_beta_num,1)) then do while ( n_to_do > 0 ) @@ -1402,6 +1566,7 @@ END_PROVIDER do j=1,det_alpha_num + det_i_prev = det_i det_i = det_alpha_order(j) if (j > 1) then TOUCH det_i @@ -1416,7 +1581,8 @@ END_PROVIDER enddo det_i = det_alpha_order(1) - SOFT_TOUCH det_i + det_i_prev = det_alpha_order(1) + SOFT_TOUCH det_i det_i_prev END_PROVIDER @@ -1449,6 +1615,7 @@ END_PROVIDER do j=1,det_beta_num + det_j_prev = det_j det_j = det_beta_order(j) if (j > 1) then TOUCH det_j @@ -1463,7 +1630,8 @@ END_PROVIDER enddo det_j = det_beta_order(1) - SOFT_TOUCH det_j + det_j_prev = det_beta_order(1) + SOFT_TOUCH det_j det_j_prev END_PROVIDER