PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-06-26 15:12:19 +02:00
Code Releases Activity

check TC energy after rotations

Browse Source
This commit is contained in:
Abdallah Ammar 2023-05-13 21:43:01 +02:00
parent 8cefe6eb44
commit 19f2ede59c
2 changed files with 31 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

31
src/tc_scf/routines_rotates.irp.f
View File

@ -93,14 +93,22 @@ subroutine routine_save_rotated_mos(thr_deg, good_angles)
integer :: i, j, k, n_degen_list, m, n, n_degen, ilast, ifirst
double precision :: max_angle, norm
double precision :: E_old, E_new, E_thr
integer, allocatable :: list_degen(:,:)
double precision, allocatable :: new_angles(:)
double precision, allocatable :: mo_r_coef_old(:,:), mo_l_coef_old(:,:)
double precision, allocatable :: mo_r_coef_good(:,:), mo_l_coef_good(:,:)
double precision, allocatable :: mo_r_coef_new(:,:)
double precision, allocatable :: fock_diag(:),s_mat(:,:)
double precision, allocatable :: fock_diag(:), s_mat(:,:)
double precision, allocatable :: stmp(:,:), T(:,:), Snew(:,:), smat2(:,:)
double precision, allocatable :: mo_l_coef_tmp(:,:), mo_r_coef_tmp(:,:), mo_l_coef_new(:,:)
E_thr = 1d-8
E_old = TC_HF_energy
allocate(mo_l_coef_old(ao_num,mo_num), mo_r_coef_old(ao_num,mo_num))
mo_r_coef_old = mo_r_coef
mo_l_coef_old = mo_l_coef
good_angles = .False.
allocate(mo_l_coef_good(ao_num, mo_num), mo_r_coef_good(ao_num,mo_num))
@ -253,11 +261,32 @@ subroutine routine_save_rotated_mos(thr_deg, good_angles)
call ezfio_set_bi_ortho_mos_mo_r_coef(mo_r_coef)
TOUCH mo_l_coef mo_r_coef
! check if TC energy has changed
E_new = TC_HF_energy
if(dabs(E_new - E_old) .gt. E_thr) then
mo_r_coef = mo_r_coef_old
mo_l_coef = mo_l_coef_old
deallocate(mo_l_coef_old, mo_r_coef_old)
call ezfio_set_bi_ortho_mos_mo_l_coef(mo_l_coef)
call ezfio_set_bi_ortho_mos_mo_r_coef(mo_r_coef)
TOUCH mo_l_coef mo_r_coef
print*, ' TC energy bef rotation = ', E_old
print*, ' TC energy aft rotation = ', E_new
print*, ' the rotation is refused'
stop
endif
allocate(new_angles(mo_num))
new_angles(1:mo_num) = dabs(angle_left_right(1:mo_num))
max_angle = maxval(new_angles)
good_angles = max_angle.lt.45.d0
print *, ' max_angle = ', max_angle
deallocate(new_angles)
deallocate(mo_l_coef_old, mo_r_coef_old)
deallocate(mo_l_coef_good, mo_r_coef_good)
deallocate(mo_r_coef_new)
end

2
src/tc_scf/tc_scf.irp.f
View File

@ -53,7 +53,7 @@ program tc_scf
stop
endif
!call minimize_tc_orb_angles()
call minimize_tc_orb_angles()
endif