added print_angle_tc_orb.irp.f

src/bi_ortho_mos/overlap.irp.f

@@ -118,3 +118,21 @@ END_PROVIDER
    enddo
  enddo
 END_PROVIDER 
+
+BEGIN_PROVIDER [ double precision, angle_left_right, (mo_num)]
+ implicit none
+ BEGIN_DOC
+  ! angle_left_right(i) = angle between the left-eigenvector chi_i and the right-eigenvector phi_i
+ END_DOC
+ integer :: i,j
+ double precision :: left,right,arg
+ do i = 1, mo_num
+  left = overlap_mo_l(i,i)
+  right = overlap_mo_r(i,i)
+  arg = min(overlap_bi_ortho(i,i)/(left*right),1.d0)
+  arg = max(arg,-1.d0)
+  angle_left_right(i) = dacos(arg) * 180.d0/dacos(-1.d0)
+ enddo
+END_PROVIDER 
+
+

src/non_hermit_dav/lapack_diag_non_hermit.irp.f

@@ -1393,7 +1393,7 @@ subroutine impose_orthog_degen_eigvec(n, e0, C0)
     deg_num(i) = 1
   enddo
 
-  de_thr = 1d-10
+  de_thr = thr_degen_tc
 
   do i = 1, n-1
     ei = e0(i)
@@ -1853,7 +1853,7 @@ subroutine impose_biorthog_degen_eigvec(n, e0, L0, R0)
     deg_num(i) = 1
   enddo
 
-  de_thr = 1d-10
+  de_thr = thr_degen_tc
 
   do i = 1, n-1
     ei = e0(i)
@@ -1955,7 +1955,7 @@ subroutine impose_orthog_biorthog_degen_eigvec(n, e0, L0, R0)
     deg_num(i) = 1
   enddo
 
-  de_thr = 1d-10
+  de_thr = thr_degen_tc
 
   do i = 1, n-1
     ei = e0(i)
@@ -2049,7 +2049,7 @@ subroutine impose_unique_biorthog_degen_eigvec(n, e0, C0, W0, L0, R0)
     deg_num(i) = 1
   enddo
 
-  de_thr = 1d-10
+  de_thr = thr_degen_tc
 
   do i = 1, n-1
     ei = e0(i)

src/tc_keywords/EZFIO.cfg

@@ -123,3 +123,10 @@ type: integer
 doc: type of 1-body Jastrow
 interface: ezfio, provider, ocaml
 default: 0
+
+
+[thr_degen_tc]
+type: Threshold
+doc: Threshold to determine if two orbitals are degenerate in TCSCF in order to avoid random quasi orthogonality between the right- and left-eigenvector for the same eigenvalue
+interface: ezfio,provider,ocaml
+default: 1.e-6

src/tc_scf/print_angle_tc_orb.irp.f

@@ -0,0 +1,20 @@
+program print_angles
+ implicit none
+  my_grid_becke  = .True.
+!  my_n_pt_r_grid = 30
+!  my_n_pt_a_grid = 50
+  my_n_pt_r_grid = 10 ! small grid for quick debug
+  my_n_pt_a_grid = 14 ! small grid for quick debug
+  call routine
+end
+subroutine routine
+ implicit none
+ integer :: i,j
+ double precision :: left,right
+ print*,'energy,product of norms, angle between vectors'                                                                  
+ do i = 1, mo_num
+  left = overlap_mo_l(i,i)
+  right = overlap_mo_r(i,i)
+  print*,Fock_matrix_tc_mo_tot(i,i),left*right,angle_left_right(i)
+ enddo
+end

src/tc_scf/rotate_tcscf_orbitals.irp.f

@@ -142,7 +142,7 @@ subroutine rotate_degen_eigvec(n, m, e0, C0, W0, L0, R0)
     deg_num(i) = 1
   enddo
 
-  de_thr = 1d-10
+  de_thr = thr_degen_tc
 
   do i = 1, m-1
     ei = e0(i)

src/tc_scf/tc_scf.irp.f

@@ -166,9 +166,10 @@ subroutine routine_scf()
   endif
 
   print*,'Energy converged !'
-  print*,'Diagonal Fock elements '
+  print*,'Final TC energy = ', TC_HF_energy
+  print*,'Diag Fock elem, product of left/right norm, angle left/right '
   do i = 1, mo_num
-   print*,i,Fock_matrix_tc_mo_tot(i,i)
+   write(*,'(I3,X,100(F16.10,X))')i,Fock_matrix_tc_mo_tot(i,i),overlap_mo_l(i,i)*overlap_mo_r(i,i),angle_left_right(i)
   enddo
 
   deallocate(rho_old, rho_new)