From 90f78a9b782842429c970213018381c9516ef241 Mon Sep 17 00:00:00 2001 From: v1j4y Date: Thu, 3 Nov 2022 15:36:07 +0100 Subject: [PATCH 1/3] Remove QR for CSFs because it does not work. --- src/utils/linear_algebra.irp.f | 191 +++++++++++++++++---------------- 1 file changed, 96 insertions(+), 95 deletions(-) diff --git a/src/utils/linear_algebra.irp.f b/src/utils/linear_algebra.irp.f index 3463370f..16f67182 100644 --- a/src/utils/linear_algebra.irp.f +++ b/src/utils/linear_algebra.irp.f @@ -1133,102 +1133,103 @@ subroutine ortho_svd(A,LDA,m,n) end -subroutine ortho_qr_withB(A,LDA,B,m,n) - implicit none - BEGIN_DOC - ! Orthogonalization using Q.R factorization - ! - ! A : Overlap Matrix - ! - ! LDA : leftmost dimension of A - ! - ! m : Number of rows of A - ! - ! n : Number of columns of A - ! - ! B : Output orthogonal basis - ! - END_DOC - integer, intent(in) :: m,n, LDA - double precision, intent(inout) :: A(LDA,n) - double precision, intent(inout) :: B(LDA,n) +! QR to orthonormalize CSFs does not work :-( +!subroutine ortho_qr_withB(A,LDA,B,m,n) +! implicit none +! BEGIN_DOC +! ! Orthogonalization using Q.R factorization +! ! +! ! A : Overlap Matrix +! ! +! ! LDA : leftmost dimension of A +! ! +! ! m : Number of rows of A +! ! +! ! n : Number of columns of A +! ! +! ! B : Output orthogonal basis +! ! +! END_DOC +! integer, intent(in) :: m,n, LDA +! double precision, intent(inout) :: A(LDA,n) +! double precision, intent(inout) :: B(LDA,n) +! +! integer :: LWORK, INFO +! integer, allocatable :: jpvt(:) +! double precision, allocatable :: TAU(:), WORK(:) +! double precision, allocatable :: C(:,:) +! double precision :: norm +! integer :: i,j +! +! allocate (TAU(min(m,n)), WORK(1)) +! allocate (jpvt(n)) +! !print *," In function ortho" +! B = A +! +! jpvt(1:n)=1 +! +! LWORK=-1 +! call dgeqp3( m, n, A, LDA, jpvt, TAU, WORK, LWORK, INFO ) +! +! ! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 +! LWORK=max(n,int(WORK(1))) +! +! deallocate(WORK) +! allocate(WORK(LWORK)) +! call dgeqp3(m, n, A, LDA, jpvt, TAU, WORK, LWORK, INFO ) +! print *,A +! print *,jpvt +! deallocate(WORK,TAU) +! !stop +! +! !LWORK=-1 +! !call dgeqrf( m, n, A, LDA, TAU, WORK, LWORK, INFO ) +! !! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 +! !LWORK=max(n,int(WORK(1))) +! +! !deallocate(WORK) +! !allocate(WORK(LWORK)) +! !call dgeqrf(m, n, A, LDA, TAU, WORK, LWORK, INFO ) +! +! !LWORK=-1 +! !call dorgqr(m, n, n, A, LDA, TAU, WORK, LWORK, INFO) +! !! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 +! !LWORK=max(n,int(WORK(1))) +! +! !deallocate(WORK) +! !allocate(WORK(LWORK)) +! !call dorgqr(m, n, n, A, LDA, TAU, WORK, LWORK, INFO) +! ! +! !allocate(C(LDA,n)) +! !call dgemm('N','N',m,n,n,1.0d0,B,LDA,A,LDA,0.0d0,C,LDA) +! !norm = 0.0d0 +! !B = 0.0d0 +! !!print *,C +! !do i=1,m +! ! norm = 0.0d0 +! ! do j=1,n +! ! norm = norm + C(j,i)*C(j,i) +! ! end do +! ! norm = 1.0d0/dsqrt(norm) +! ! do j=1,n +! ! B(j,i) = C(j,i) +! ! end do +! !end do +! !print *,B +! +! +! !deallocate(WORK,TAU) +!end - integer :: LWORK, INFO - integer, allocatable :: jpvt(:) - double precision, allocatable :: TAU(:), WORK(:) - double precision, allocatable :: C(:,:) - double precision :: norm - integer :: i,j - - allocate (TAU(min(m,n)), WORK(1)) - allocate (jpvt(n)) - !print *," In function ortho" - B = A - - jpvt(1:n)=1 - - LWORK=-1 - call dgeqp3( m, n, A, LDA, jpvt, TAU, WORK, LWORK, INFO ) - - ! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 - LWORK=max(n,int(WORK(1))) - - deallocate(WORK) - allocate(WORK(LWORK)) - call dgeqp3(m, n, A, LDA, jpvt, TAU, WORK, LWORK, INFO ) - print *,A - print *,jpvt - deallocate(WORK,TAU) - !stop - - !LWORK=-1 - !call dgeqrf( m, n, A, LDA, TAU, WORK, LWORK, INFO ) - !! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 - !LWORK=max(n,int(WORK(1))) - - !deallocate(WORK) - !allocate(WORK(LWORK)) - !call dgeqrf(m, n, A, LDA, TAU, WORK, LWORK, INFO ) - - !LWORK=-1 - !call dorgqr(m, n, n, A, LDA, TAU, WORK, LWORK, INFO) - !! /!\ int(WORK(1)) becomes negative when WORK(1) > 2147483648 - !LWORK=max(n,int(WORK(1))) - - !deallocate(WORK) - !allocate(WORK(LWORK)) - !call dorgqr(m, n, n, A, LDA, TAU, WORK, LWORK, INFO) - ! - !allocate(C(LDA,n)) - !call dgemm('N','N',m,n,n,1.0d0,B,LDA,A,LDA,0.0d0,C,LDA) - !norm = 0.0d0 - !B = 0.0d0 - !!print *,C - !do i=1,m - ! norm = 0.0d0 - ! do j=1,n - ! norm = norm + C(j,i)*C(j,i) - ! end do - ! norm = 1.0d0/dsqrt(norm) - ! do j=1,n - ! B(j,i) = C(j,i) - ! end do - !end do - !print *,B - - - !deallocate(WORK,TAU) -end - -subroutine ortho_qr_csf(A, LDA, B, m, n) bind(C, name="ortho_qr_csf") - use iso_c_binding - integer(c_int32_t), value :: LDA - integer(c_int32_t), value :: m - integer(c_int32_t), value :: n - integer(c_int16_t) :: A(LDA,n) - integer(c_int16_t) :: B(LDA,n) - call ortho_qr_withB(A,LDA,B,m,n) -end subroutine ortho_qr_csf +!subroutine ortho_qr_csf(A, LDA, B, m, n) bind(C, name="ortho_qr_csf") +! use iso_c_binding +! integer(c_int32_t), value :: LDA +! integer(c_int32_t), value :: m +! integer(c_int32_t), value :: n +! integer(c_int16_t) :: A(LDA,n) +! integer(c_int16_t) :: B(LDA,n) +! call ortho_qr_withB(A,LDA,B,m,n) +!end subroutine ortho_qr_csf subroutine ortho_qr(A,LDA,m,n) implicit none From 3e44e79be5000dbc6ae791f78ffec721a953c503 Mon Sep 17 00:00:00 2001 From: v1j4y Date: Thu, 3 Nov 2022 15:38:05 +0100 Subject: [PATCH 2/3] Deactive qr also in CFG utils function. --- src/csf/cfgCI_utils.c | 35 ++++++++++++++++++----------------- 1 file changed, 18 insertions(+), 17 deletions(-) diff --git a/src/csf/cfgCI_utils.c b/src/csf/cfgCI_utils.c index 76b64dd0..3510db37 100644 --- a/src/csf/cfgCI_utils.c +++ b/src/csf/cfgCI_utils.c @@ -255,23 +255,24 @@ void generateAllBFs(int64_t Isomo, int64_t MS, Tree *bftree, int *NBF, int *NSOM void ortho_qr_csf(double *overlapMatrix, int lda, double *orthoMatrix, int rows, int cols); -void gramSchmidt_qp(double *overlapMatrix, int rows, int cols, double *orthoMatrix){ - int i,j; - //for(j=0;j Date: Thu, 3 Nov 2022 15:40:05 +0100 Subject: [PATCH 3/3] Remove debug printing and fix compilation issue. --- src/csf/sigma_vector.irp.f | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/src/csf/sigma_vector.irp.f b/src/csf/sigma_vector.irp.f index 21c19aaa..3fa437a5 100644 --- a/src/csf/sigma_vector.irp.f +++ b/src/csf/sigma_vector.irp.f @@ -112,7 +112,7 @@ ! endif !enddo n_CSF = 0 - print *," -9(((((((((((((( NSOMOMin=",NSOMOMin + !print *," -9(((((((((((((( NSOMOMin=",NSOMOMin ncfgprev = cfg_seniority_index(NSOMOMin) ! can be -1 if(ncfgprev.eq.-1)then ncfgprev=1 @@ -145,10 +145,10 @@ endif endif n_CSF += ncfg*dimcsfpercfg - print *," i=",i," dimcsf=",dimcsfpercfg," ncfg=",ncfg, " ncfgprev=",ncfgprev, " senor=",cfg_seniority_index(i) + !print *," i=",i," dimcsf=",dimcsfpercfg," ncfg=",ncfg, " ncfgprev=",ncfgprev, " senor=",cfg_seniority_index(i) ncfgprev = cfg_seniority_index(i+2) end do - print *," ^^^^^ N_CSF = ",n_CSF," N_CFG=",N_configuration + !print *," ^^^^^ N_CSF = ",n_CSF," N_CFG=",N_configuration END_PROVIDER @@ -240,7 +240,7 @@ end subroutine get_phase_qp_to_cfg ! initialization psi_coef_config = 0.d0 DetToCSFTransformationMatrix(0,:,:) = 1.d0 - do i = 2-iand(MS,1), NSOMOMax,2 + do i = 2-iand(MS,1_8), NSOMOMax,2 Isomo = IBSET(0_8, i) - 1_8 ! rows = Ncsfs ! cols = Ndets @@ -342,8 +342,8 @@ end subroutine get_phase_qp_to_cfg nsomomin = elec_alpha_num-elec_beta_num rowsmax = 0 colsmax = 0 - print *,"NSOMOMax = ",NSOMOMax - print *,"NSOMOMin = ",NSOMOMin + !print *,"NSOMOMax = ",NSOMOMax + !print *,"NSOMOMin = ",NSOMOMin !allocate(AIJpqMatrixDimsList(NSOMOMax,NSOMOMax,4,NSOMOMax,NSOMOMax,2)) ! Type ! 1. SOMO -> SOMO @@ -525,7 +525,7 @@ end subroutine get_phase_qp_to_cfg end do end do end do - print *,"Rowsmax=",rowsmax," Colsmax=",colsmax + !print *,"Rowsmax=",rowsmax," Colsmax=",colsmax END_PROVIDER BEGIN_PROVIDER [ real*8, AIJpqContainer, (NBFMax,NBFmax,NSOMOMax+1,NSOMOMax+1,4,NSOMOMin:NSOMOMax)]