mirror of
https://github.com/QuantumPackage/qp2.git
synced 2024-12-22 19:43:32 +01:00
add dgemm in max overlap
This commit is contained in:
parent
b68fe9feb2
commit
e9c900fe00
@ -8,26 +8,28 @@
|
|||||||
double precision, allocatable :: proj(:)
|
double precision, allocatable :: proj(:)
|
||||||
integer, allocatable :: iorder(:)
|
integer, allocatable :: iorder(:)
|
||||||
double precision, allocatable :: mo_coef_tmp(:,:)
|
double precision, allocatable :: mo_coef_tmp(:,:)
|
||||||
allocate(overlap(mo_num,mo_num),proj(mo_num),iorder(mo_num),mo_coef_tmp(ao_num,mo_num))
|
double precision, allocatable :: tmp(:,:)
|
||||||
|
allocate(overlap(mo_num,mo_num),proj(mo_num),iorder(mo_num),mo_coef_tmp(ao_num,mo_num),tmp(mo_num,ao_num))
|
||||||
|
|
||||||
overlap(:,:) = 0d0
|
overlap(:,:) = 0d0
|
||||||
mo_coef_tmp(:,:) = 0d0
|
mo_coef_tmp(:,:) = 0d0
|
||||||
proj(:) = 0d0
|
proj(:) = 0d0
|
||||||
iorder(:) = 0d0
|
iorder(:) = 0d0
|
||||||
|
tmp(:,:) = 0d0
|
||||||
|
|
||||||
|
! These matrix products compute the overlap bewteen the initial and the current MOs
|
||||||
|
call dgemm('T','N', mo_num, ao_num, ao_num, 1.d0, &
|
||||||
|
mo_coef_begin_iteration, size(mo_coef_begin_iteration,1), &
|
||||||
|
ao_overlap, size(ao_overlap,1), 0.d0, &
|
||||||
|
tmp, size(tmp,1))
|
||||||
|
|
||||||
|
call dgemm('N','N', mo_num, mo_num, ao_num, 1.d0, &
|
||||||
|
tmp, size(tmp,1), &
|
||||||
|
mo_coef, size(mo_coef, 1), 0.d0, &
|
||||||
|
overlap, size(overlap,1) )
|
||||||
|
|
||||||
! this loop compute the overlap between the initial and the current MOS
|
|
||||||
do i = 1, mo_num ! old mo
|
|
||||||
do j = 1, mo_num ! curent mo
|
|
||||||
do k = 1, ao_num
|
|
||||||
do l = 1, ao_num
|
|
||||||
overlap(i,j) += mo_coef_begin_iteration(k,i)* ao_overlap(k,l) * mo_coef(l,j)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
|
|
||||||
! for each orbital compute the best overlap
|
! for each orbital compute the best overlap
|
||||||
|
|
||||||
do i = 1, mo_num
|
do i = 1, mo_num
|
||||||
iorder(i) = i ! initialize the iorder list as we need it to sort later
|
iorder(i) = i ! initialize the iorder list as we need it to sort later
|
||||||
do j = 1, elec_alpha_num
|
do j = 1, elec_alpha_num
|
||||||
@ -52,17 +54,21 @@
|
|||||||
integer, allocatable :: iorder_alpha(:)
|
integer, allocatable :: iorder_alpha(:)
|
||||||
allocate(overlap_alpha(mo_num,elec_alpha_num),proj_alpha(elec_alpha_num),iorder_alpha(elec_alpha_num))
|
allocate(overlap_alpha(mo_num,elec_alpha_num),proj_alpha(elec_alpha_num),iorder_alpha(elec_alpha_num))
|
||||||
overlap_alpha(:,:) = 0d0
|
overlap_alpha(:,:) = 0d0
|
||||||
|
mo_coef_tmp(:,:) = 0d0
|
||||||
proj_alpha(:) = 0d0
|
proj_alpha(:) = 0d0
|
||||||
iorder_alpha(:) = 0d0
|
iorder_alpha(:) = 0d0
|
||||||
do i = 1, mo_num ! old mo
|
tmp(:,:) = 0d0
|
||||||
do j = 1, elec_alpha_num ! curent mo
|
! These matrix products compute the overlap bewteen the initial and the current MOs
|
||||||
do k = 1, ao_num
|
call dgemm('T','N', mo_num, ao_num, ao_num, 1.d0, &
|
||||||
do l = 1, ao_num
|
mo_coef_begin_iteration, size(mo_coef_begin_iteration,1), &
|
||||||
overlap_alpha(i,j) += mo_coef_begin_iteration(k,i) * ao_overlap(k,l) * mo_coef(l,j)
|
ao_overlap, size(ao_overlap,1), 0.d0, &
|
||||||
enddo
|
tmp, size(tmp,1))
|
||||||
enddo
|
|
||||||
enddo
|
call dgemm('N','N', mo_num, elec_alpha_num, ao_num, 1.d0, &
|
||||||
enddo
|
tmp, size(tmp,1), &
|
||||||
|
mo_coef, size(mo_coef, 1), 0.d0, &
|
||||||
|
overlap_alpha, size(overlap_alpha,1) )
|
||||||
|
|
||||||
do i = 1, elec_alpha_num
|
do i = 1, elec_alpha_num
|
||||||
iorder_alpha(i) = i ! initialize the iorder list as we need it to sort later
|
iorder_alpha(i) = i ! initialize the iorder list as we need it to sort later
|
||||||
do j = 1, elec_beta_num
|
do j = 1, elec_beta_num
|
||||||
@ -70,6 +76,7 @@
|
|||||||
enddo
|
enddo
|
||||||
proj_alpha(i) = dsqrt(proj_alpha(i))
|
proj_alpha(i) = dsqrt(proj_alpha(i))
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
! sort the list of projection to find the mos with the largest overlap
|
! sort the list of projection to find the mos with the largest overlap
|
||||||
call dsort(proj_alpha(:),iorder_alpha(:),elec_alpha_num)
|
call dsort(proj_alpha(:),iorder_alpha(:),elec_alpha_num)
|
||||||
! reorder orbitals according to projection
|
! reorder orbitals according to projection
|
||||||
@ -79,10 +86,11 @@
|
|||||||
do i=1,elec_alpha_num
|
do i=1,elec_alpha_num
|
||||||
mo_coef(:,i) = mo_coef_tmp(:,i)
|
mo_coef(:,i) = mo_coef_tmp(:,i)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
|
deallocate(overlap_alpha, proj_alpha, iorder_alpha)
|
||||||
endif
|
endif
|
||||||
|
|
||||||
deallocate(overlap, proj, iorder, mo_coef_tmp)
|
deallocate(overlap, proj, iorder, mo_coef_tmp, tmp)
|
||||||
|
|
||||||
|
|
||||||
end
|
end
|
||||||
|
|
||||||
|
@ -99,8 +99,8 @@ END_DOC
|
|||||||
call reorder_mo_max_overlap
|
call reorder_mo_max_overlap
|
||||||
endif
|
endif
|
||||||
if(frozen_orb_scf)then
|
if(frozen_orb_scf)then
|
||||||
call reorder_core_orb
|
call reorder_core_orb
|
||||||
call initialize_mo_coef_begin_iteration
|
call initialize_mo_coef_begin_iteration
|
||||||
endif
|
endif
|
||||||
|
|
||||||
TOUCH MO_coef
|
TOUCH MO_coef
|
||||||
@ -112,12 +112,12 @@ END_DOC
|
|||||||
|
|
||||||
energy_SCF = SCF_energy
|
energy_SCF = SCF_energy
|
||||||
Delta_Energy_SCF = energy_SCF - energy_SCF_previous
|
Delta_Energy_SCF = energy_SCF - energy_SCF_previous
|
||||||
! if ( (SCF_algorithm == 'DIIS').and.(Delta_Energy_SCF > 0.d0) ) then
|
if ( (SCF_algorithm == 'DIIS').and.(Delta_Energy_SCF > 0.d0).and.(.not.do_mom) ) then
|
||||||
! Fock_matrix_AO(1:ao_num,1:ao_num) = Fock_matrix_DIIS (1:ao_num,1:ao_num,index_dim_DIIS)
|
Fock_matrix_AO(1:ao_num,1:ao_num) = Fock_matrix_DIIS (1:ao_num,1:ao_num,index_dim_DIIS)
|
||||||
! Fock_matrix_AO_alpha = Fock_matrix_AO*0.5d0
|
Fock_matrix_AO_alpha = Fock_matrix_AO*0.5d0
|
||||||
! Fock_matrix_AO_beta = Fock_matrix_AO*0.5d0
|
Fock_matrix_AO_beta = Fock_matrix_AO*0.5d0
|
||||||
! TOUCH Fock_matrix_AO_alpha Fock_matrix_AO_beta
|
TOUCH Fock_matrix_AO_alpha Fock_matrix_AO_beta
|
||||||
! endif
|
endif
|
||||||
|
|
||||||
if (.not.do_mom) then
|
if (.not.do_mom) then
|
||||||
double precision :: level_shift_save
|
double precision :: level_shift_save
|
||||||
|
Loading…
Reference in New Issue
Block a user