subroutine save_mos implicit none double precision, allocatable :: buffer(:,:) integer :: i,j call system('$QP_ROOT/scripts/save_current_mos.sh '//trim(ezfio_filename)) call ezfio_set_mo_basis_mo_tot_num(mo_tot_num) call ezfio_set_mo_basis_mo_label(mo_label) call ezfio_set_mo_basis_ao_md5(ao_md5) allocate ( buffer(ao_num,mo_tot_num) ) buffer = 0.d0 do j = 1, mo_tot_num do i = 1, ao_num buffer(i,j) = mo_coef(i,j) enddo enddo call ezfio_set_mo_basis_mo_coef(buffer) call ezfio_set_mo_basis_mo_occ(mo_occ) deallocate (buffer) end subroutine save_mos_truncated(n) implicit none double precision, allocatable :: buffer(:,:) integer :: i,j,n call system('$QP_ROOT/scripts/save_current_mos.sh '//trim(ezfio_filename)) call ezfio_set_mo_basis_mo_tot_num(n) call ezfio_set_mo_basis_mo_label(mo_label) call ezfio_set_mo_basis_ao_md5(ao_md5) allocate ( buffer(ao_num,n) ) buffer = 0.d0 do j = 1, n do i = 1, ao_num buffer(i,j) = mo_coef(i,j) enddo enddo call ezfio_set_mo_basis_mo_coef(buffer) call ezfio_set_mo_basis_mo_occ(mo_occ) deallocate (buffer) end subroutine mo_as_eigvectors_of_mo_matrix(matrix,n,m,label,sign) implicit none integer,intent(in) :: n,m, sign character*(64), intent(in) :: label double precision, intent(in) :: matrix(n,m) integer :: i,j double precision, allocatable :: mo_coef_new(:,:), R(:,:),eigvalues(:), A(:,:) !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, R call write_time(output_mo_basis) if (m /= mo_tot_num) then print *, irp_here, ': Error : m/= mo_tot_num' stop 1 endif allocate(A(n,m),R(n,m),mo_coef_new(ao_num_align,m),eigvalues(m)) if (sign == -1) then do j=1,m do i=1,n A(i,j) = -matrix(i,j) enddo enddo else do j=1,m do i=1,n A(i,j) = matrix(i,j) enddo enddo endif mo_coef_new = mo_coef call lapack_diag(eigvalues,R,A,n,m) write (output_mo_basis,'(A)'), 'MOs are now **'//trim(label)//'**' write (output_mo_basis,'(A)'), '' write (output_mo_basis,'(A)'), 'Eigenvalues' write (output_mo_basis,'(A)'), '-----------' write (output_mo_basis,'(A)'), '' write (output_mo_basis,'(A)'), '======== ================' if (sign == -1) then do i=1,m eigvalues(i) = -eigvalues(i) enddo endif do i=1,m write (output_mo_basis,'(I8,X,F16.10)'), i,eigvalues(i) enddo write (output_mo_basis,'(A)'), '======== ================' write (output_mo_basis,'(A)'), '' call dgemm('N','N',ao_num,m,m,1.d0,mo_coef_new,size(mo_coef_new,1),R,size(R,1),0.d0,mo_coef,size(mo_coef,1)) deallocate(A,mo_coef_new,R,eigvalues) call write_time(output_mo_basis) mo_label = label end subroutine mo_as_eigvectors_of_mo_matrix_sort_by_observable(matrix,observable,n,m,label) implicit none integer,intent(in) :: n,m character*(64), intent(in) :: label double precision, intent(in) :: matrix(n,m),observable(n,n) double precision, allocatable :: mo_coef_new(:,:), R(:,:),eigvalues(:),value(:) integer,allocatable :: iorder(:) !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, R call write_time(output_mo_basis) if (m /= mo_tot_num) then print *, irp_here, ': Error : m/= mo_tot_num' stop 1 endif allocate(R(n,m),mo_coef_new(ao_num_align,m),eigvalues(m),value(m),iorder(m)) mo_coef_new = mo_coef call lapack_diag(eigvalues,R,matrix,size(matrix,1),size(matrix,2)) do i = 1, mo_tot_num value(i) = 0.d0 iorder(i) = i do j = 1, mo_tot_num do k = 1, mo_tot_num value(i) += R(k,i) * R(j,i) * observable(k,j) enddo enddo ! print*,'value(i) = ',i,value(i) enddo integer :: i,j,k,index double precision :: R_tmp(m,m),obs(m) print*,'sort ....' call dsort(value,iorder,n) do i = 1, mo_tot_num index = iorder(i) ! print*,'iorder(i) = ',iorder(i) obs(i) = eigvalues(iorder(i)) do j = 1, mo_tot_num R_tmp(j,i) = R(j,index) enddo enddo do i = 1, mo_tot_num do j = 1, mo_tot_num R(j,i) = R_tmp(j,i) enddo enddo do i = 1, mo_tot_num value(i) = 0.d0 do j = 1, mo_tot_num do k = 1, mo_tot_num value(i) += R(k,i) * R(j,i) * observable(k,j) enddo enddo print*,'i = ',i print*,'value(i) = ',value(i) print*,'obs(i) = ',obs(i) print*,'' enddo write (output_mo_basis,'(A)'), 'MOs are now **'//trim(label)//'**' write (output_mo_basis,'(A)'), '' write (output_mo_basis,'(A)'), 'Eigenvalues' write (output_mo_basis,'(A)'), '-----------' write (output_mo_basis,'(A)'), '' write (output_mo_basis,'(A)'), '======== ================' do i = 1, m write (output_mo_basis,'(I8,X,F16.10)'), i,eigvalues(i) enddo write (output_mo_basis,'(A)'), '======== ================' write (output_mo_basis,'(A)'), '' call dgemm('N','N',ao_num,m,m,1.d0,mo_coef_new,size(mo_coef_new,1),R,size(R,1),0.d0,mo_coef,size(mo_coef,1)) deallocate(mo_coef_new,R,eigvalues) call write_time(output_mo_basis) mo_label = label SOFT_TOUCH mo_coef mo_label end subroutine mo_sort_by_observable(observable,label) implicit none character*(64), intent(in) :: label double precision, intent(in) :: observable(mo_tot_num) double precision, allocatable :: mo_coef_new(:,:),value(:) integer,allocatable :: iorder(:) allocate(mo_coef_new(ao_num_align,mo_tot_num),value(mo_tot_num),iorder(mo_tot_num)) print*,'allocate !' mo_coef_new = mo_coef do i = 1, mo_tot_num iorder(i) = i value(i) = observable(i) enddo integer :: i,j,k,index print*,'sort ....' call dsort(value,iorder,mo_tot_num) do i = 1, mo_tot_num index = iorder(i) do j = 1, mo_tot_num mo_coef(j,i) = mo_coef_new(j,index) enddo enddo write (output_mo_basis,'(A)'), 'MOs are now **'//trim(label)//'**' write (output_mo_basis,'(A)'), '' deallocate(mo_coef_new,value) ! call write_time(output_mo_basis) mo_label = label SOFT_TOUCH mo_coef mo_label end