LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-12-22 20:35:19 +01:00
Code Issues Releases Wiki Activity

Merge pull request #222 from scemama/master

Browse Source 
Bug fixes
This commit is contained in:
Thomas Applencourt 2018-01-10 14:20:24 -06:00 committed by GitHub
commit df3963d1f4
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
41 changed files with 260 additions and 373 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

16
configure vendored
View File

@ -480,24 +480,24 @@ def create_ninja_and_rc(l_installed):
'export QP_PYTHON={0}'.format(":".join(l_python)), "",
'export IRPF90={0}'.format(path_irpf90.replace(QP_ROOT,"${QP_ROOT}")),
'export NINJA={0}'.format(path_ninja.replace(QP_ROOT,"${QP_ROOT}")),
'qp_append_export () {',
' #Append path $2:${!1}. Add the semicolon only if ${!1} is defined',
' echo ${2}${!1:+:${!1}}',
'function qp_append_export () {',
' #Append path $2:${!1}. Add the semicolon only if ${!1} is defined',
' eval "value_1=\"\${$1}\""',
' echo ${2}${value_1:+:${value_1}}',
'}',
'export PYTHONPATH=$(qp_append_export "PYTHONPATH" "${QP_EZFIO}/Python":"${QP_PYTHON}")',
'export PATH=$(qp_append_export "PATH" "${QP_PYTHON}":"${QP_ROOT}"/bin:"${QP_ROOT}"/ocaml)',
'export LD_LIBRARY_PATH=$(qp_append_export "LD_LIBRARY_PATH" "${QP_ROOT}"/lib:"${QP_ROOT}"/lib64)',
'export LIBRARY_PATH=$(qp_append_export "LIBRARY_PATH" "${QP_ROOT}"/lib:"${QP_ROOT}"/lib64)',
'export C_INCLUDE_PATH=$(qp_append_export "C_INCLUDE_PATH" "${QP_ROOT}"/include)'
'export C_INCLUDE_PATH=$(qp_append_export "C_INCLUDE_PATH" "${QP_ROOT}"/include)',
'',
'source ${QP_ROOT}/install/EZFIO/Bash/ezfio.sh',
'if [[ $SHELL == "bash" ]] ; then',
' source ${QP_ROOT}/install/EZFIO/Bash/ezfio.sh',
'fi',
'',
'# Choose the correct network interface',
'# export QP_NIC=ib0',
'# export QP_NIC=eth0',
'',
'# Choose how to start MPI processes',
'# export QP_MPIRUN="mpirun"',
''
]

5
ocaml/Gto.ml
View File

@ -36,7 +36,10 @@ let of_prim_coef_list pc =
let read_one in_channel =
(* Fetch number of lines to read on first line *)
let buffer = input_line in_channel in
let buffer =
try input_line in_channel with
| End_of_file -> raise End_Of_Basis
in
if ( (String_ext.strip buffer) = "" ) then
raise End_Of_Basis;
let sym_str = String.sub buffer 0 2 in

43
ocaml/Input_determinants_by_hand.ml
View File

@ -11,6 +11,7 @@ module Determinants_by_hand : sig
expected_s2 : Positive_float.t;
psi_coef : Det_coef.t array;
psi_det : Determinant.t array;
state_average_weight : Positive_float.t array;
} [@@deriving sexp]
val read : unit -> t
val read_maybe : unit -> t option
@ -30,6 +31,7 @@ end = struct
expected_s2 : Positive_float.t;
psi_coef : Det_coef.t array;
psi_det : Determinant.t array;
state_average_weight : Positive_float.t array;
} [@@deriving sexp]
;;
@ -82,7 +84,6 @@ end = struct
|> Ezfio.set_determinants_n_det
;;
let read_n_states () =
if not (Ezfio.has_determinants_n_states ()) then
Ezfio.set_determinants_n_states 1
@ -96,6 +97,36 @@ end = struct
|> Ezfio.set_determinants_n_states
;;
let write_state_average_weight data =
let n_states =
read_n_states ()
|> States_number.to_int
in
let data =
Array.map ~f:Positive_float.to_float data
|> Array.to_list
in
Ezfio.ezfio_array_of_list ~rank:1 ~dim:[| n_states |] ~data
|> Ezfio.set_determinants_state_average_weight
;;
let read_state_average_weight () =
if not (Ezfio.has_determinants_state_average_weight ()) then
begin
let n_states =
read_n_states ()
|> States_number.to_int
in
let data =
Array.init n_states (fun _ -> 1./.(float_of_int n_states))
|> Array.map ~f:Positive_float.of_float
in
write_state_average_weight data;
end;
Ezfio.get_determinants_state_average_weight ()
|> Ezfio.flattened_ezfio
|> Array.map ~f:Positive_float.of_float
;;
let read_expected_s2 () =
if not (Ezfio.has_determinants_expected_s2 ()) then
@ -205,6 +236,7 @@ end = struct
psi_coef = read_psi_coef () ;
psi_det = read_psi_det () ;
n_states = read_n_states () ;
state_average_weight = read_state_average_weight () ;
}
else
failwith "No molecular orbitals, so no determinants"
@ -228,6 +260,7 @@ end = struct
psi_coef ;
psi_det ;
n_states ;
state_average_weight ;
} =
write_n_int n_int ;
write_bit_kind bit_kind;
@ -236,6 +269,7 @@ end = struct
write_expected_s2 expected_s2;
write_psi_coef ~n_det:n_det ~n_states:n_states psi_coef ;
write_psi_det ~n_int:n_int ~n_det:n_det psi_det;
write_state_average_weight state_average_weight;
;;
@ -288,12 +322,17 @@ Number of determinants ::
n_det = %s
State average weights ::
state_average_weight = (%s)
Determinants ::
%s
"
(b.expected_s2 |> Positive_float.to_string)
(b.n_det |> Det_number.to_string)
(b.state_average_weight |> Array.to_list |> List.map ~f:Positive_float.to_string |> String.concat ~sep:"\t")
det_text
|> Rst_string.of_string
;;
@ -307,6 +346,7 @@ bit_kind = %s
n_det = %s
n_states = %s
expected_s2 = %s
state_average_weight = %s
psi_coef = %s
psi_det = %s
"
@ -315,6 +355,7 @@ psi_det = %s
(b.n_det |> Det_number.to_string)
(b.n_states |> States_number.to_string)
(b.expected_s2 |> Positive_float.to_string)
(b.state_average_weight |> Array.to_list |> List.map ~f:Positive_float.to_string |> String.concat ~sep:",")
(b.psi_coef |> Array.to_list |> List.map ~f:Det_coef.to_string
|> String.concat ~sep:", ")
(b.psi_det |> Array.to_list |> List.map ~f:(Determinant.to_string

11
ocaml/Zmatrix.ml
View File

@ -196,17 +196,6 @@ let rotation_matrix axis angle =
(2. *. (b *. d -. a *. c),
2. *. (c *. d +. a *. b),
a *. a +. d *. d -. b *. b -. c *. c)]
(*
[(a *. a +. b *. b -. c *. c -. d *. d,
2. *. (b *. c +. a *. d),
2. *. (b *. d -. a *. c));
(2. *. (b *. c -. a *. d),
a *. a +. c *. c -.b *. b -. d *. d,
2. *. (c *. d +. a *. b));
(2. *. (b *. d +. a *. c),
2. *. (c *. d -. a *. b),
a *. a +. d *. d -. b *. b -. c *. c)]
*)

14
plugins/CIS/super_ci.irp.f
View File

@ -12,12 +12,12 @@ subroutine super_CI
integer :: k
character :: save_char
call write_time(output_hartree_fock)
write(output_hartree_fock,'(A4,X,A16, X, A16, X, A16 )') &
call write_time(6)
write(6,'(A4,X,A16, X, A16, X, A16 )') &
'====','================','================','================'
write(output_hartree_fock,'(A4,X,A16, X, A16, X, A16 )') &
write(6,'(A4,X,A16, X, A16, X, A16 )') &
' N ', 'Energy ', 'Energy diff ', 'Save '
write(output_hartree_fock,'(A4,X,A16, X, A16, X, A16 )') &
write(6,'(A4,X,A16, X, A16, X, A16 )') &
'====','================','================','================'
E = HF_energy + 1.d0
@ -39,7 +39,7 @@ subroutine super_CI
save_char = ' '
endif
E_min = min(E,E_min)
write(output_hartree_fock,'(I4,X,F16.10, X, F16.10, X, A8 )') &
write(6,'(I4,X,F16.10, X, F16.10, X, A8 )') &
k, E, delta_E, save_char
if ( (delta_E < 0.d0).and.(dabs(delta_E) < thresh_scf) ) then
exit
@ -55,8 +55,8 @@ subroutine super_CI
TOUCH mo_coef
enddo
write(output_hartree_fock,'(A4,X,A16, X, A16, X, A16 )') &
write(6,'(A4,X,A16, X, A16, X, A16 )') &
'====','================','================','================'
call write_time(output_hartree_fock)
call write_time(6)
end

13
plugins/Full_CI_ZMQ/pt2_stoch_routines.irp.f
View File

@ -25,8 +25,8 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
double precision :: sumabove(comb_teeth), sum2above(comb_teeth), Nabove(comb_teeth)
double precision, external :: omp_get_wtime
double precision :: state_average_weight_save(N_states), w(N_states)
double precision :: time
double precision :: w(N_states)
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
if (N_det < max(10,N_states)) then
@ -35,11 +35,12 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
error(:) = 0.d0
else
state_average_weight_save(:) = state_average_weight(:)
do pt2_stoch_istate=1,N_states
SOFT_TOUCH pt2_stoch_istate
w(:) = 0.d0
w(pt2_stoch_istate) = 1.d0
call update_psi_average_norm_contrib(w)
state_average_weight(:) = 0.d0
state_average_weight(pt2_stoch_istate) = 1.d0
TOUCH state_average_weight
allocate(pt2_detail(N_states,N_det_generators+1), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
sumabove = 0d0
@ -141,7 +142,9 @@ subroutine ZMQ_pt2(E, pt2,relative_error, absolute_error, error)
deallocate(pt2_detail, comb, computed, tbc)
enddo
FREE psi_average_norm_contrib pt2_stoch_istate
FREE pt2_stoch_istate
state_average_weight(:) = state_average_weight_save(:)
TOUCH state_average_weight
endif
do k=N_det+1,N_states
pt2(k) = 0.d0

4
plugins/Generators_CAS/generators.irp.f
View File

@ -7,7 +7,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
END_DOC
integer :: i,k,l
logical :: good
call write_time(output_determinants)
call write_time(6)
N_det_generators = 0
do i=1,N_det
do l=1,n_cas_bitmask
@ -28,7 +28,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
endif
enddo
N_det_generators = max(N_det_generators,1)
call write_int(output_determinants,N_det_generators,'Number of generators')
call write_int(6,N_det_generators,'Number of generators')
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ]

15
plugins/Generators_full/generators.irp.f
View File

@ -8,7 +8,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
END_DOC
integer :: i
double precision :: norm
call write_time(output_determinants)
call write_time(6)
norm = 0.d0
N_det_generators = N_det
do i=1,N_det
@ -19,7 +19,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
endif
enddo
N_det_generators = max(N_det_generators,1)
call write_int(output_determinants,N_det_generators,'Number of generators')
call write_int(6,N_det_generators,'Number of generators')
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ]
@ -30,15 +30,8 @@ END_PROVIDER
! Hartree-Fock determinant
END_DOC
integer :: i, k
psi_coef_generators = 0.d0
psi_det_generators = 0_bit_kind
do i=1,N_det_generators
do k=1,N_int
psi_det_generators(k,1,i) = psi_det_sorted(k,1,i)
psi_det_generators(k,2,i) = psi_det_sorted(k,2,i)
enddo
psi_coef_generators(i,:) = psi_coef_sorted(i,:)
enddo
psi_det_generators(1:N_int,1:2,1:N_det) = psi_det_sorted(1:N_int,1:2,1:N_det)
psi_coef_generators(1:N_det,1:N_states) = psi_coef_sorted(1:N_det,1:N_states)
END_PROVIDER

2
plugins/Generators_restart/generators.irp.f
View File

@ -14,7 +14,7 @@ BEGIN_PROVIDER [ integer, N_det_generators ]
else
print*,'PB in generators restart !!!'
endif
call write_int(output_determinants,N_det_generators,'Number of generators')
call write_int(6,N_det_generators,'Number of generators')
END_PROVIDER

20
plugins/Hartree_Fock/Roothaan_Hall_SCF.irp.f
View File

@ -15,18 +15,20 @@ END_DOC
integer :: i,j
double precision, allocatable :: mo_coef_save(:,:)
PROVIDE ao_md5 mo_occ level_shift
allocate(mo_coef_save(ao_num,mo_tot_num), &
Fock_matrix_DIIS (ao_num,ao_num,max_dim_DIIS), &
error_matrix_DIIS(ao_num,ao_num,max_dim_DIIS) &
)
call write_time(output_hartree_fock)
call write_time(6)
write(output_hartree_fock,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
'====','================','================','================'
write(output_hartree_fock,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
' N ', 'Energy ', 'Energy diff ', 'DIIS error '
write(output_hartree_fock,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
'====','================','================','================'
! Initialize energies and density matrices
@ -115,7 +117,7 @@ END_DOC
! Print results at the end of each iteration
write(output_hartree_fock,'(I4, 1X, F16.10, 1X, F16.10, 1X, F16.10, 1X, I3)') &
write(6,'(I4, 1X, F16.10, 1X, F16.10, 1X, F16.10, 1X, I3)') &
iteration_SCF, energy_SCF, Delta_energy_SCF, max_error_DIIS, dim_DIIS
if (Delta_energy_SCF < 0.d0) then
@ -128,18 +130,18 @@ END_DOC
! End of Main SCF loop
!
write(output_hartree_fock,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
'====','================','================','================'
write(output_hartree_fock,*)
write(6,*)
if(.not.no_oa_or_av_opt)then
call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1),size(Fock_matrix_mo,2),mo_label,1,.true.)
endif
call write_double(output_hartree_fock, Energy_SCF, 'Hartree-Fock energy')
call write_double(6, Energy_SCF, 'Hartree-Fock energy')
call ezfio_set_hartree_fock_energy(Energy_SCF)
call write_time(output_hartree_fock)
call write_time(6)
end

2
plugins/Hartree_Fock/SCF.irp.f
View File

@ -23,7 +23,7 @@ subroutine create_guess
mo_coef = ao_ortho_lowdin_coef
TOUCH mo_coef
mo_label = 'Guess'
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral,size(mo_mono_elec_integral,1),size(mo_mono_elec_integral,2),mo_label,.false.)
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral,size(mo_mono_elec_integral,1),size(mo_mono_elec_integral,2),mo_label,1,.false.)
SOFT_TOUCH mo_coef mo_label
else if (mo_guess_type == "Huckel") then
call huckel_guess

18
plugins/Hartree_Fock/damping_SCF.irp.f
View File

@ -28,13 +28,13 @@ subroutine damping_SCF
enddo
call write_time(output_hartree_fock)
call write_time(6)
write(output_hartree_fock,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
write(6,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
'====','================','================','================', '===='
write(output_hartree_fock,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
write(6,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
' N ', 'Energy ', 'Energy diff ', 'Density diff ', 'Save'
write(output_hartree_fock,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
write(6,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') &
'====','================','================','================', '===='
E = HF_energy + 1.d0
@ -58,7 +58,7 @@ subroutine damping_SCF
save_char = ' '
endif
write(output_hartree_fock,'(I4,1X,F16.10, 1X, F16.10, 1X, F16.10, 3X, A )') &
write(6,'(I4,1X,F16.10, 1X, F16.10, 1X, F16.10, 3X, A )') &
k, E, delta_E, delta_D, save_char
D_alpha = HF_density_matrix_ao_alpha
@ -115,17 +115,17 @@ subroutine damping_SCF
TOUCH mo_coef
enddo
write(output_hartree_fock,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') '====','================','================','================', '===='
write(output_hartree_fock,*)
write(6,'(A4,1X,A16, 1X, A16, 1X, A16, 1X, A4 )') '====','================','================','================', '===='
write(6,*)
if(.not.no_oa_or_av_opt)then
call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1),size(Fock_matrix_mo,2),mo_label,1,.true.)
endif
call write_double(output_hartree_fock, E_min, 'Hartree-Fock energy')
call write_double(6, E_min, 'Hartree-Fock energy')
call ezfio_set_hartree_fock_energy(E_min)
call write_time(output_hartree_fock)
call write_time(6)
deallocate(D_alpha,D_beta,F_new,D_new_alpha,D_new_beta,delta_alpha,delta_beta)
end

8
plugins/MRCC_Utils/mrcc_utils.irp.f
View File

@ -211,7 +211,7 @@ END_PROVIDER
call davidson_diag_mrcc_HS2(psi_det,eigenvectors, &
size(eigenvectors,1), &
eigenvalues,N_det,N_states,N_states_diag,N_int, &
output_determinants,mrcc_state)
6,mrcc_state)
CI_eigenvectors_dressed(1:N_det,mrcc_state) = eigenvectors(1:N_det,mrcc_state)
CI_electronic_energy_dressed(mrcc_state) = eigenvalues(mrcc_state)
enddo
@ -316,12 +316,12 @@ BEGIN_PROVIDER [ double precision, CI_energy_dressed, (N_states_diag) ]
integer :: j
character*(8) :: st
call write_time(output_determinants)
call write_time(6)
do j=1,min(N_det,N_states)
write(st,'(I4)') j
CI_energy_dressed(j) = CI_electronic_energy_dressed(j) + nuclear_repulsion
call write_double(output_determinants,CI_energy_dressed(j),'Energy of state '//trim(st))
call write_double(output_determinants,CI_eigenvectors_s2_dressed(j),'S^2 of state '//trim(st))
call write_double(6,CI_energy_dressed(j),'Energy of state '//trim(st))
call write_double(6,CI_eigenvectors_s2_dressed(j),'S^2 of state '//trim(st))
enddo
END_PROVIDER

6
plugins/MRPT_Utils/mrpt_utils.irp.f
View File

@ -354,12 +354,12 @@ BEGIN_PROVIDER [ double precision, CI_dressed_pt2_new_energy, (N_states) ]
integer :: j
character*(8) :: st
call write_time(output_determinants)
call write_time(6)
do j=1,N_states
CI_dressed_pt2_new_energy(j) = CI_electronic_dressed_pt2_new_energy(j) + nuclear_repulsion
write(st,'(I4)') j
call write_double(output_determinants,CI_dressed_pt2_new_energy(j),'Energy of state '//trim(st))
call write_double(output_determinants,CI_dressed_pt2_new_eigenvectors_s2(j),'S^2 of state '//trim(st))
call write_double(6,CI_dressed_pt2_new_energy(j),'Energy of state '//trim(st))
call write_double(6,CI_dressed_pt2_new_eigenvectors_s2(j),'S^2 of state '//trim(st))
enddo
END_PROVIDER

2
plugins/Perturbation/selection.irp.f
View File

@ -125,7 +125,7 @@ subroutine remove_small_contributions
if (N_removed > 0) then
N_det = N_det - N_removed
SOFT_TOUCH N_det psi_det psi_coef
call write_int(output_determinants,N_removed, 'Removed determinants')
call write_int(6,N_removed, 'Removed determinants')
endif
end

1
plugins/Selectors_Utils/zmq.irp.f
View File

@ -45,6 +45,7 @@ integer function zmq_get_$X(zmq_to_qp_run_socket, worker_id)
integer :: rc
character*(256) :: msg
PROVIDE zmq_state
zmq_get_$X = 0
if (mpi_master) then

4
plugins/Selectors_full/selectors.irp.f
View File

@ -8,7 +8,7 @@ BEGIN_PROVIDER [ integer, N_det_selectors]
END_DOC
integer :: i
double precision :: norm, norm_max
call write_time(output_determinants)
call write_time(6)
N_det_selectors = N_det
if (threshold_generators < 1.d0) then
norm = 0.d0
@ -21,7 +21,7 @@ BEGIN_PROVIDER [ integer, N_det_selectors]
enddo
N_det_selectors = max(N_det_selectors,N_det_generators)
endif
call write_int(output_determinants,N_det_selectors,'Number of selectors')
call write_int(6,N_det_selectors,'Number of selectors')
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_selectors, (N_int,2,psi_selectors_size) ]

4
plugins/Selectors_no_sorted/selectors.irp.f
View File

@ -8,11 +8,11 @@ BEGIN_PROVIDER [ integer, N_det_selectors]
END_DOC
integer :: i
double precision :: norm
call write_time(output_determinants)
call write_time(6)
norm = 0.d0
N_det_selectors = N_det
N_det_selectors = max(N_det_selectors,1)
call write_int(output_determinants,N_det_selectors,'Number of selectors')
call write_int(6,N_det_selectors,'Number of selectors')
END_PROVIDER

11
plugins/Symmetry/Symmetry.main.irp.f
View File

@ -3,9 +3,17 @@ program Symmetry
BEGIN_DOC
! TODO
END_DOC
integer :: i, j
integer :: i, j, k
character*8 :: sym
do k=1,n_irrep
print *, sym_operation(k)
do i=1,mo_tot_num
print '(1000(F8.4,X))', mo_symm(i,:,k), sum(mo_symm(i,:,k))
enddo
print *, ''
enddo
print *, 'Molecule is linear: ', molecule_is_linear
print *, 'Has center of inversion: ', molecule_has_center_of_inversion
print *, 'Has S2n improper rotation: ', molecule_has_improper_rotation
@ -17,5 +25,4 @@ program Symmetry
do i=1,n_irrep
print *, i, real(character_table(i,:))
enddo
PROVIDE mo_sym
end

20
plugins/Symmetry/aos.irp.f
View File

@ -22,15 +22,19 @@ subroutine generate_sym_coord(n_sym_points,result)
BEGIN_DOC
! xyz coordinates of points to check the symmetry, drawn uniformly in the molecular box.
END_DOC
integer :: i, xyz
integer :: i, iop
do i=1,n_sym_points
call random_number(result(1,i))
call random_number(result(2,i))
call random_number(result(3,i))
enddo
do xyz=1,3
result(xyz,1:n_sym_points) = sym_box(xyz,1) + result(xyz,:) * (sym_box(xyz,2)-sym_box(xyz,1))
double precision, external :: halton_ranf
do i=1,n_sym_points,n_irrep
result(1,i) = sym_box(1,1) + halton_ranf(1) * (sym_box(1,2)-sym_box(1,1))
result(2,i) = sym_box(1,1) + halton_ranf(2) * (sym_box(2,2)-sym_box(2,1))
result(3,i) = sym_box(1,1) + halton_ranf(3) * (sym_box(3,2)-sym_box(3,1))
do iop=2,n_irrep
if (iop-1+i > n_sym_points) exit
call dgemm('N','N',3,1,3,1.d0,sym_transformation_matrices(1,1,iop), &
size(sym_transformation_matrices,1),&
result(1,i),size(result,1),0.d0,result(1,i+iop-1),size(result,1))
enddo
enddo
end

108
plugins/Symmetry/find_sym.irp.f
View File

@ -47,6 +47,7 @@ BEGIN_PROVIDER [ integer, sym_rotation_axis, (3) ]
logical :: found
double precision, external :: u_dot_u
integer :: iorder, iaxis
do iaxis=1,3
do iorder=12,2,-1
sym_rotation_axis(iaxis) = iorder
@ -300,14 +301,6 @@ BEGIN_PROVIDER [ character*16, point_group ]
END_PROVIDER
BEGIN_PROVIDER [ character*8, mo_symmetry ]
implicit none
BEGIN_DOC
! Symmetry of the MOs
END_DOC
integer :: i,j
END_PROVIDER
BEGIN_PROVIDER [ integer, n_irrep ]
implicit none
@ -369,105 +362,6 @@ BEGIN_PROVIDER [ integer, mo_sym, (mo_tot_num) ]
double precision :: sym_operations_on_mos(mo_tot_num)
logical :: possible_irrep(n_irrep,mo_tot_num)
n_sym_points = 10
allocate(val(n_sym_points,mo_tot_num,2), sym_points(3,n_sym_points), ref_points(3,n_sym_points))
call generate_sym_coord(n_sym_points,ref_points)
call compute_sym_mo_values(ref_points,n_sym_points,val(1,1,2))
possible_irrep = .True.
do iop=1,n_irrep
if (sym_operation(iop) == 'E') then
cycle
endif
if (sym_operation(iop) == 'i') then
do ipoint=1,n_sym_points
call sym_apply_inversion(ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 'sh') then
do ipoint=1,n_sym_points
call sym_apply_reflexion(molecule_principal_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 's') then
do ipoint=1,n_sym_points
call sym_apply_reflexion(molecule_principal_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 'sv') then
do ipoint=1,n_sym_points
call sym_apply_reflexion(molecule_ternary_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 'sd') then
angle = dble(maxval(sym_rotation_axis))
do ipoint=1,n_sym_points
call sym_apply_diagonal_reflexion(angle,molecule_principal_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 'C2''') then
angle = 2.d0
do ipoint=1,n_sym_points
call sym_apply_rotation(angle,molecule_secondary_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop) == 'C2"') then
angle = 2.d0
do ipoint=1,n_sym_points
call sym_apply_rotation(angle,molecule_ternary_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else
do l=2,len(sym_operation(iop))
if (sym_operation(iop)(l:l) == '^') exit
enddo
read(sym_operation(iop)(2:l-1), *) iangle
if (l == len(sym_operation(iop))+1) then
l=1
else
read(sym_operation(iop)(l+1:), *, err=10, end=10) l
10 continue
endif
angle = dble(iangle)/(dble(l))
if (sym_operation(iop)(1:1) == 'C') then
do ipoint=1,n_sym_points
call sym_apply_rotation(angle,molecule_principal_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
else if (sym_operation(iop)(1:1) == 'S') then
do ipoint=1,n_sym_points
call sym_apply_improper_rotation(angle,molecule_principal_axis,ref_points(1,ipoint),sym_points(1,ipoint))
enddo
endif
endif
call compute_sym_mo_values(sym_points,n_sym_points,val(1,1,1))
print *, sym_operation(iop)
double precision :: icount
do imo=1,mo_tot_num
sym_operations_on_mos(imo) = 0.d0
icount = 0
do ipoint=1,n_sym_points
double precision :: x
if (dabs(val(ipoint,imo,1)) < 1.d-5) cycle
icount += 1.d0
x = val(ipoint,imo,1)/val(ipoint,imo,2)
if (dabs(x) > 1.d0) then
x = 1.d0/x
endif
sym_operations_on_mos(imo) += x
enddo
sym_operations_on_mos(imo) *= 1.d0/icount
if (dabs(sym_operations_on_mos(imo) - 1.d0) < 1.d-2) then
sym_operations_on_mos(imo) = 1.d0
else if (dabs(sym_operations_on_mos(imo) + 1.d0) < 1.d-2) then
sym_operations_on_mos(imo) = -1.d0
else if (dabs(sym_operations_on_mos(imo)) < 1.d-2) then
sym_operations_on_mos(imo) = 0.d0
endif
print *, imo, sym_operations_on_mos(imo)
do i=1,n_irrep
if (dabs(character_table(i,iop) - sym_operations_on_mos(imo)) > 1.d-2) then
possible_irrep(i,imo) = .False.
endif
enddo
enddo
enddo
do imo=1,mo_tot_num
print *, 'MO ', imo
do i=1,n_irrep

22
plugins/Symmetry/nuclei.irp.f
View File

@ -69,26 +69,26 @@ BEGIN_PROVIDER [ double precision, nucl_coord_sym, (nucl_num,3) ]
character*(64), parameter :: ft= '(A16, 4(1X,A12 ))'
double precision, parameter :: a0= 0.529177249d0
call write_time(output_Nuclei)
write(output_Nuclei,'(A)') ''
write(output_Nuclei,'(A)') 'Nuclear Coordinates in standard orientation (Angstroms)'
write(output_Nuclei,'(A)') '======================================================='
write(output_Nuclei,'(A)') ''
write(output_Nuclei,ft) &
call write_time(6)
write(6,'(A)') ''
write(6,'(A)') 'Nuclear Coordinates in standard orientation (Angstroms)'
write(6,'(A)') '======================================================='
write(6,'(A)') ''
write(6,ft) &
'================','============','============','============','============'
write(output_Nuclei,*) &
write(6,*) &
' Atom Charge X Y Z '
write(output_Nuclei,ft) &
write(6,ft) &
'================','============','============','============','============'
do i=1,nucl_num
write(output_Nuclei,f) nucl_label(i), nucl_charge(i), &
write(6,f) nucl_label(i), nucl_charge(i), &
nucl_coord_sym(i,1)*a0, &
nucl_coord_sym(i,2)*a0, &
nucl_coord_sym(i,3)*a0
enddo
write(output_Nuclei,ft) &
write(6,ft) &
'================','============','============','============','============'
write(output_Nuclei,'(A)') ''
write(6,'(A)') ''
endif

2
plugins/Symmetry/sym_operation.irp.f
View File

@ -23,7 +23,7 @@ subroutine sym_apply_diagonal_reflexion(angle,iaxis,point_in,point_out)
double precision :: point_tmp1(3), point_tmp2(3)
integer :: iaxis2
iaxis2 = mod(iaxis,3)+1
iaxis2 = mod(iaxis2,3)+1
! iaxis2 = mod(iaxis2,3)+1
call sym_apply_rotation(-angle,iaxis,point_in,point_tmp1)
call sym_apply_reflexion(iaxis2,point_tmp1,point_tmp2)
call sym_apply_rotation(angle,iaxis,point_tmp2,point_out)

13
plugins/mrcepa0/dressing.irp.f
View File

@ -963,6 +963,19 @@ END_PROVIDER
enddo
end do
end do
! else if(mrmode == 10) then
! do i = 1, N_det_ref
! do i_state = 1, N_states
! delta_ii(i_state,i)= delta_ii_mrsc2(i_state,i)
! delta_ii_s2(i_state,i)= delta_ii_s2_mrsc2(i_state,i)
! enddo
! do j = 1, N_det_non_ref
! do i_state = 1, N_states
! delta_ij(i_state,j,i) = delta_ij_mrsc2(i_state,j,i)
! delta_ij_s2(i_state,j,i) = delta_ij_s2_mrsc2(i_state,j,i)
! enddo
! end do
! end do
else if(mrmode == 5) then
do i = 1, N_det_ref
do i_state = 1, N_states

10
plugins/mrcepa0/mrcc_stoch_routines.irp.f
View File

@ -32,17 +32,15 @@ subroutine ZMQ_mrcc(E, mrcc, delta, delta_s2, relative_error)
double precision, external :: omp_get_wtime
double precision :: time
double precision :: w(N_states)
state_average_weight(:) = 0.d0
state_average_weight(mrcc_stoch_istate) = 1.d0
TOUCH state_average_weight
provide nproc fragment_first fragment_count mo_bielec_integrals_in_map mo_mono_elec_integral mrcc_weight psi_selectors
w(:) = 0.d0
w(mrcc_stoch_istate) = 1.d0
call update_psi_average_norm_contrib(w)
print *, '========== ================= ================= ================='

3
scripts/ezfio_interface/ei_handler.py
View File

@ -329,7 +329,8 @@ def create_ezfio_provider(dict_ezfio_cfg):
ez_p.set_doc(dict_info['doc'])
ez_p.set_ezfio_dir(dict_info['ezfio_dir'])
ez_p.set_ezfio_name(dict_info['ezfio_name'])
ez_p.set_output("output_%s" % dict_info['module'].lower)
ez_p.set_output("6")
# ez_p.set_output("output_%s" % dict_info['module'].lower)
# (nuclei.nucl_num,pseudo.klocmax) => (nucl_num,klocmax)
ez_p.set_size(re.sub(r'\w+\.', "", dict_info['size']))

10
scripts/generate_h_apply.py
View File

@ -142,7 +142,7 @@ class H_apply(object):
endif
SOFT_TOUCH psi_det psi_coef N_det
"""
s["printout_now"] = """write(output_determinants,*) &
s["printout_now"] = """write(6,*) &
100.*float(i_generator)/float(N_det_generators), '% in ', wall_1-wall_0, 's'"""
self.data = s
@ -370,9 +370,9 @@ class H_apply(object):
delta_pt2(k) = 0.d0
pt2_old(k) = 0.d0
enddo
write(output_determinants,'(A12, 1X, A8, 3(2X, A9), 2X, A8, 2X, A8, 2X, A8)') &
write(6,'(A12, 1X, A8, 3(2X, A9), 2X, A8, 2X, A8, 2X, A8)') &
'N_generators', 'Norm', 'Delta PT2', 'PT2', 'Est. PT2', 'secs'
write(output_determinants,'(A12, 1X, A8, 3(2X, A9), 2X, A8, 2X, A8, 2X, A8)') &
write(6,'(A12, 1X, A8, 3(2X, A9), 2X, A8, 2X, A8, 2X, A8)') &
'============', '========', '=========', '=========', '=========', &
'========='
"""
@ -385,7 +385,7 @@ class H_apply(object):
"""
self.data["printout_now"] = """
do k=1,N_st
write(output_determinants,'(I10, 4(2X, F9.6), 2X, F8.1)') &
write(6,'(I10, 4(2X, F9.6), 2X, F8.1)') &
i_generator, norm_psi(k), delta_pt2(k), pt2(k), &
pt2(k)/(norm_psi(k)*norm_psi(k)), &
wall_1-wall_0
@ -416,7 +416,7 @@ class H_apply(object):
SOFT_TOUCH psi_det psi_coef N_det
selection_criterion_min = min(selection_criterion_min, maxval(select_max))*0.1d0
selection_criterion = selection_criterion_min
call write_double(output_determinants,selection_criterion,'Selection criterion')
call write_double(6,selection_criterion,'Selection criterion')
"""
self.data["keys_work"] = """
e_2_pert_buffer = 0.d0

2
src/Bitmask/bitmasks.irp.f
View File

@ -100,7 +100,7 @@ BEGIN_PROVIDER [ integer, N_generators_bitmask ]
! Number of bitmasks for generators
END_DOC
logical :: exists
PROVIDE ezfio_filename
PROVIDE ezfio_filename N_int
if (mpi_master) then
call ezfio_has_bitmasks_N_mask_gen(exists)

10
src/Davidson/diagonalize_CI.irp.f
View File

@ -7,14 +7,14 @@ BEGIN_PROVIDER [ double precision, CI_energy, (N_states_diag) ]
integer :: j
character*(8) :: st
call write_time(output_determinants)
call write_time(6)
do j=1,min(N_det,N_states_diag)
CI_energy(j) = CI_electronic_energy(j) + nuclear_repulsion
enddo
do j=1,min(N_det,N_states)
write(st,'(I4)') j
call write_double(output_determinants,CI_energy(j),'Energy of state '//trim(st))
call write_double(output_determinants,CI_eigenvectors_s2(j),'S^2 of state '//trim(st))
call write_double(6,CI_energy(j),'Energy of state '//trim(st))
call write_double(6,CI_eigenvectors_s2(j),'S^2 of state '//trim(st))
enddo
END_PROVIDER
@ -58,14 +58,14 @@ END_PROVIDER
! call davidson_diag(psi_det,CI_eigenvectors,CI_electronic_energy, &
! size(CI_eigenvectors,1), &
! N_det,min(N_det,N_states),min(N_det,N_states_diag),N_int,output_determinants)
! N_det,min(N_det,N_states),min(N_det,N_states_diag),N_int,6)
!
! call u_0_S2_u_0(CI_eigenvectors_s2,CI_eigenvectors,N_det,psi_det,N_int,&
! min(N_det,N_states_diag),size(CI_eigenvectors,1))
call davidson_diag_HS2(psi_det,CI_eigenvectors, CI_eigenvectors_s2, &
size(CI_eigenvectors,1),CI_electronic_energy, &
N_det,min(N_det,N_states),min(N_det,N_states_diag),N_int,output_determinants)
N_det,min(N_det,N_states),min(N_det,N_states_diag),N_int,6)
else if (diag_algorithm == "Lapack") then

2
src/Determinants/H_apply.irp.f
View File

@ -265,7 +265,7 @@ subroutine remove_duplicates_in_psi_det(found_duplicates)
endif
enddo
N_det = k
call write_bool(output_determinants,found_duplicates,'Found duplicate determinants')
call write_bool(6,found_duplicates,'Found duplicate determinants')
SOFT_TOUCH N_det psi_det psi_coef
endif
deallocate (duplicate,bit_tmp)

6
src/Determinants/density_matrix.irp.f
View File

@ -368,13 +368,13 @@ BEGIN_PROVIDER [ double precision, state_average_weight, (N_states) ]
END_DOC
logical :: exists
state_average_weight = 1.d0
state_average_weight(:) = 1.d0
call ezfio_has_determinants_state_average_weight(exists)
if (exists) then
call ezfio_get_determinants_state_average_weight(state_average_weight)
endif
state_average_weight = state_average_weight+1.d-31
state_average_weight = state_average_weight/(sum(state_average_weight))
state_average_weight(:) = state_average_weight(:)+1.d-31
state_average_weight(:) = state_average_weight(:)/(sum(state_average_weight(:)))
END_PROVIDER

48
src/Determinants/determinants.irp.f
View File

@ -43,7 +43,7 @@ BEGIN_PROVIDER [ integer, N_det ]
else
N_det = 1
endif
call write_int(output_determinants,N_det,'Number of determinants')
call write_int(6,N_det,'Number of determinants')
endif
IRP_IF MPI
include 'mpif.h'
@ -77,7 +77,7 @@ BEGIN_PROVIDER [ integer, psi_det_size ]
BEGIN_DOC
! Size of the psi_det/psi_coef arrays
END_DOC
PROVIDE ezfio_filename output_determinants
PROVIDE ezfio_filename
logical :: exists
if (mpi_master) then
call ezfio_has_determinants_n_det(exists)
@ -87,7 +87,7 @@ BEGIN_PROVIDER [ integer, psi_det_size ]
psi_det_size = 1
endif
psi_det_size = max(psi_det_size,100000)
call write_int(output_determinants,psi_det_size,'Dimension of the psi arrays')
call write_int(6,psi_det_size,'Dimension of the psi arrays')
endif
IRP_IF MPI
include 'mpif.h'
@ -225,34 +225,6 @@ BEGIN_PROVIDER [ double precision, psi_coef, (psi_det_size,N_states) ]
END_PROVIDER
subroutine update_psi_average_norm_contrib(w)
implicit none
BEGIN_DOC
! Compute psi_average_norm_contrib for different state average weights w(:)
END_DOC
double precision, intent(in) :: w(N_states)
double precision :: w0(N_states), f
w0(:) = w(:)/sum(w(:))
integer :: i,j,k
do i=1,N_det
psi_average_norm_contrib(i) = psi_coef(i,1)*psi_coef(i,1)*w(1)
enddo
do k=2,N_states
do i=1,N_det
psi_average_norm_contrib(i) = psi_average_norm_contrib(i) + &
psi_coef(i,k)*psi_coef(i,k)*w(k)
enddo
enddo
f = 1.d0/sum(psi_average_norm_contrib(1:N_det))
do i=1,N_det
psi_average_norm_contrib(i) = psi_average_norm_contrib(i)*f
enddo
SOFT_TOUCH psi_average_norm_contrib
end subroutine
BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (psi_det_size) ]
implicit none
BEGIN_DOC
@ -260,14 +232,12 @@ BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (psi_det_size) ]
END_DOC
integer :: i,j,k
double precision :: f
f = 1.d0/dble(N_states)
do i=1,N_det
psi_average_norm_contrib(i) = psi_coef(i,1)*psi_coef(i,1)*f
enddo
do k=2,N_states
psi_average_norm_contrib(:) = 0.d0
do k=1,N_states
do i=1,N_det
psi_average_norm_contrib(i) = psi_average_norm_contrib(i) + &
psi_coef(i,k)*psi_coef(i,k)*f
psi_coef(i,k)*psi_coef(i,k)*state_average_weight(k)
enddo
enddo
f = 1.d0/sum(psi_average_norm_contrib(1:N_det))
@ -562,7 +532,7 @@ subroutine save_wavefunction_general(ndet,nstates,psidet,dim_psicoef,psicoef)
call ezfio_set_determinants_psi_coef(psi_coef_save)
deallocate (psi_coef_save)
call write_int(output_determinants,ndet,'Saved determinants')
call write_int(6,ndet,'Saved determinants')
endif
end
@ -634,7 +604,7 @@ subroutine save_wavefunction_specified(ndet,nstates,psidet,psicoef,ndetsave,inde
enddo
call ezfio_set_determinants_psi_coef(psi_coef_save)
call write_int(output_determinants,ndet,'Saved determinants')
call write_int(6,ndet,'Saved determinants')
deallocate (psi_coef_save)
end

2
src/Determinants/psi_cas.irp.f
View File

@ -41,7 +41,7 @@ use bitmasks
enddo
endif
enddo
call write_int(output_determinants,N_det_cas, 'Number of determinants in the CAS')
call write_int(6,N_det_cas, 'Number of determinants in the CAS')
END_PROVIDER

3
src/Determinants/zmq.irp.f
View File

@ -87,6 +87,7 @@ integer function zmq_get_$X(zmq_to_qp_run_socket, worker_id)
integer :: rc
character*(256) :: msg
PROVIDE zmq_state
zmq_get_$X = 0
if (mpi_master) then
write(msg,'(A,1X,I8,1X,A200)') 'get_data '//trim(zmq_state), worker_id, '$X'
@ -272,6 +273,7 @@ integer function zmq_get_psi_det(zmq_to_qp_run_socket, worker_id)
integer*8 :: rc8
character*(256) :: msg
PROVIDE zmq_state
zmq_get_psi_det = 0
if (mpi_master) then
write(msg,'(A,1X,I8,1X,A200)') 'get_data '//trim(zmq_state), worker_id, 'psi_det'
@ -319,6 +321,7 @@ integer function zmq_get_psi_coef(zmq_to_qp_run_socket, worker_id)
integer*8 :: rc8
character*(256) :: msg
PROVIDE zmq_state psi_det_size
zmq_get_psi_coef = 0
if (mpi_master) then
write(msg,'(A,1X,I8,1X,A200)') 'get_data '//trim(zmq_state), worker_id, 'psi_coef'

30
src/Ezfio_files/output.irp.f
View File

@ -8,36 +8,6 @@
call wall_time(output_wall_time_0)
END_PROVIDER
BEGIN_SHELL [ /bin/bash ]
for NAME in $(\ls -d ${QP_ROOT}/src/*/)
do
NAME=$(basename ${NAME})
cat << EOF
BEGIN_PROVIDER [ integer, output_$NAME ]
implicit none
BEGIN_DOC
! Output file for $NAME
END_DOC
PROVIDE output_wall_time_0 output_cpu_time_0 ezfio_filename
! integer :: getUnitAndOpen
! call ezfio_set_output_empty(.False.)
IRP_IF COARRAY
if (this_image() == 1) then
output_$NAME = 6 !getUnitAndOpen(trim(ezfio_filename)//'/output/'//'$NAME.rst','a')
else
output_$NAME = getUnitAndOpen('/dev/null','w')
endif
IRP_ELSE
output_$NAME = 6 !getUnitAndOpen(trim(ezfio_filename)//'/output/'//'$NAME.rst','a')
IRP_ENDIF
write(output_$NAME,'(A)') &
'--------------------------------------------------------------------------------'
END_PROVIDER
EOF
done
END_SHELL
subroutine write_time(iunit)
implicit none

10
src/Integrals_Bielec/mo_bi_integrals.irp.f
View File

@ -122,17 +122,17 @@ BEGIN_PROVIDER [ logical, mo_bielec_integrals_in_map ]
endif
else
! call add_integrals_to_map(full_ijkl_bitmask_4)
call add_integrals_to_map(full_ijkl_bitmask_4)
! call four_index_transform_zmq(ao_integrals_map,mo_integrals_map, &
! mo_coef, size(mo_coef,1), &
! 1, 1, 1, 1, ao_num, ao_num, ao_num, ao_num, &
! 1, 1, 1, 1, mo_num, mo_num, mo_num, mo_num)
!
call four_index_transform_block(ao_integrals_map,mo_integrals_map, &
mo_coef, size(mo_coef,1), &
1, 1, 1, 1, ao_num, ao_num, ao_num, ao_num, &
1, 1, 1, 1, mo_num, mo_num, mo_num, mo_num)
! call four_index_transform_block(ao_integrals_map,mo_integrals_map, &
! mo_coef, size(mo_coef,1), &
! 1, 1, 1, 1, ao_num, ao_num, ao_num, ao_num, &
! 1, 1, 1, 1, mo_num, mo_num, mo_num, mo_num)
!
! call four_index_transform(ao_integrals_map,mo_integrals_map, &
! mo_coef, size(mo_coef,1), &

2
src/MO_Basis/ao_ortho_canonical.irp.f
View File

@ -129,7 +129,7 @@ END_PROVIDER
enddo
ao_ortho_canonical_num = ao_cart_to_sphe_num
call ortho_canonical (ao_cart_to_sphe_overlap, size(ao_cart_to_sphe_overlap,1), &
call ortho_canonical(ao_cart_to_sphe_overlap, size(ao_cart_to_sphe_overlap,1), &
ao_cart_to_sphe_num, S, size(S,1), ao_ortho_canonical_num)
call dgemm('N','N', ao_num, ao_ortho_canonical_num, ao_cart_to_sphe_num, 1.d0, &

72
src/MO_Basis/utils.irp.f
View File

@ -55,7 +55,7 @@ subroutine mo_as_eigvectors_of_mo_matrix(matrix,n,m,label,sign,output)
double precision, allocatable :: mo_coef_new(:,:), R(:,:),eigvalues(:), A(:,:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, R
call write_time(output_mo_basis)
call write_time(6)
if (m /= mo_tot_num) then
print *, irp_here, ': Error : m/= mo_tot_num'
stop 1
@ -78,12 +78,12 @@ subroutine mo_as_eigvectors_of_mo_matrix(matrix,n,m,label,sign,output)
call lapack_diag(eigvalues,R,A,n,m)
if (output) then
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)') '======== ================'
write (6,'(A)') 'MOs are now **'//trim(label)//'**'
write (6,'(A)') ''
write (6,'(A)') 'Eigenvalues'
write (6,'(A)') '-----------'
write (6,'(A)') ''
write (6,'(A)') '======== ================'
endif
if (sign == -1) then
do i=1,m
@ -92,15 +92,15 @@ subroutine mo_as_eigvectors_of_mo_matrix(matrix,n,m,label,sign,output)
endif
if (output) then
do i=1,m
write (output_mo_basis,'(I8,1X,F16.10)') i,eigvalues(i)
write (6,'(I8,1X,F16.10)') i,eigvalues(i)
enddo
write (output_mo_basis,'(A)') '======== ================'
write (output_mo_basis,'(A)') ''
write (6,'(A)') '======== ================'
write (6,'(A)') ''
endif
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)
call write_time(6)
mo_label = label
end
@ -115,7 +115,7 @@ subroutine mo_as_svd_vectors_of_mo_matrix(matrix,lda,m,n,label)
double precision, allocatable :: mo_coef_new(:,:), U(:,:),D(:), A(:,:), Vt(:,:), work(:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, U, Vt, A
call write_time(output_mo_basis)
call write_time(6)
if (m /= mo_tot_num) then
print *, irp_here, ': Error : m/= mo_tot_num'
stop 1
@ -132,22 +132,22 @@ subroutine mo_as_svd_vectors_of_mo_matrix(matrix,lda,m,n,label)
call svd(A,lda,U,lda,D,Vt,lda,m,n)
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)') '======== ================'
write (6,'(A)') 'MOs are now **'//trim(label)//'**'
write (6,'(A)') ''
write (6,'(A)') 'Eigenvalues'
write (6,'(A)') '-----------'
write (6,'(A)') ''
write (6,'(A)') '======== ================'
do i=1,m
write (output_mo_basis,'(I8,1X,F16.10)') i,D(i)
write (6,'(I8,1X,F16.10)') i,D(i)
enddo
write (output_mo_basis,'(A)') '======== ================'
write (output_mo_basis,'(A)') ''
write (6,'(A)') '======== ================'
write (6,'(A)') ''
call dgemm('N','N',ao_num,m,m,1.d0,mo_coef_new,size(mo_coef_new,1),U,size(U,1),0.d0,mo_coef,size(mo_coef,1))
deallocate(A,mo_coef_new,U,Vt,D)
call write_time(output_mo_basis)
call write_time(6)
mo_label = label
end
@ -162,7 +162,7 @@ subroutine mo_as_eigvectors_of_mo_matrix_sort_by_observable(matrix,observable,n,
integer,allocatable :: iorder(:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, R
call write_time(output_mo_basis)
call write_time(6)
if (m /= mo_tot_num) then
print *, irp_here, ': Error : m/= mo_tot_num'
stop 1
@ -213,21 +213,21 @@ subroutine mo_as_eigvectors_of_mo_matrix_sort_by_observable(matrix,observable,n,
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)') '======== ================'
write (6,'(A)') 'MOs are now **'//trim(label)//'**'
write (6,'(A)') ''
write (6,'(A)') 'Eigenvalues'
write (6,'(A)') '-----------'
write (6,'(A)') ''
write (6,'(A)') '======== ================'
do i = 1, m
write (output_mo_basis,'(I8,1X,F16.10)') i,eigvalues(i)
write (6,'(I8,1X,F16.10)') i,eigvalues(i)
enddo
write (output_mo_basis,'(A)') '======== ================'
write (output_mo_basis,'(A)') ''
write (6,'(A)') '======== ================'
write (6,'(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)
call write_time(6)
mo_label = label
SOFT_TOUCH mo_coef mo_label
@ -261,12 +261,12 @@ subroutine mo_sort_by_observable(observable,label)
enddo
enddo
write (output_mo_basis,'(A)') 'MOs are now **'//trim(label)//'**'
write (output_mo_basis,'(A)') ''
write (6,'(A)') 'MOs are now **'//trim(label)//'**'
write (6,'(A)') ''
deallocate(mo_coef_new,value)
! call write_time(output_mo_basis)
! call write_time(6)
mo_label = label
SOFT_TOUCH mo_coef mo_label

37
src/Nuclei/nuclei.irp.f
View File

@ -31,26 +31,26 @@ BEGIN_PROVIDER [ double precision, nucl_coord, (nucl_num,3) ]
character*(64), parameter :: ft= '(A16, 4(1X,A12 ))'
double precision, parameter :: a0= 0.529177249d0
call write_time(output_Nuclei)
write(output_Nuclei,'(A)') ''
write(output_Nuclei,'(A)') 'Nuclear Coordinates (Angstroms)'
write(output_Nuclei,'(A)') '==============================='
write(output_Nuclei,'(A)') ''
write(output_Nuclei,ft) &
call write_time(6)
write(6,'(A)') ''
write(6,'(A)') 'Nuclear Coordinates (Angstroms)'
write(6,'(A)') '==============================='
write(6,'(A)') ''
write(6,ft) &
'================','============','============','============','============'
write(output_Nuclei,*) &
write(6,*) &
' Atom Charge X Y Z '
write(output_Nuclei,ft) &
write(6,ft) &
'================','============','============','============','============'
do i=1,nucl_num
write(output_Nuclei,f) nucl_label(i), nucl_charge(i), &
write(6,f) nucl_label(i), nucl_charge(i), &
nucl_coord(i,1)*a0, &
nucl_coord(i,2)*a0, &
nucl_coord(i,3)*a0
enddo
write(output_Nuclei,ft) &
write(6,ft) &
'================','============','============','============','============'
write(output_Nuclei,'(A)') ''
write(6,'(A)') ''
endif
@ -190,9 +190,8 @@ BEGIN_PROVIDER [ double precision, nuclear_repulsion ]
nuclear_repulsion *= 0.5d0
end if
call write_time(output_Nuclei)
call write_double(output_Nuclei,nuclear_repulsion, &
'Nuclear repulsion energy')
call write_time(6)
call write_double(6,nuclear_repulsion,'Nuclear repulsion energy')
if (disk_access_nuclear_repulsion.EQ.'Write') then
if (mpi_master) then
@ -203,8 +202,8 @@ BEGIN_PROVIDER [ double precision, nuclear_repulsion ]
END_PROVIDER
BEGIN_PROVIDER [ character*(4), element_name, (0:128)]
&BEGIN_PROVIDER [ double precision, element_mass, (0:128) ]
BEGIN_PROVIDER [ character*(4), element_name, (0:127)]
&BEGIN_PROVIDER [ double precision, element_mass, (0:127) ]
BEGIN_DOC
! Array of the name of element, sorted by nuclear charge (integer)
END_DOC
@ -216,7 +215,7 @@ END_PROVIDER
filename = trim(filename)//'/data/list_element.txt'
iunit = getUnitAndOpen(filename,'r')
element_mass(:) = 0.d0
do i=0,128
do i=0,127
write(element_name(i),'(I4)') i
enddo
character*(80) :: buffer, dummy
@ -232,11 +231,11 @@ END_PROVIDER
IRP_IF MPI
include 'mpif.h'
integer :: ierr
call MPI_BCAST( element_name, size(element_name)*4, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr)
call MPI_BCAST( element_name, 128*4, MPI_CHARACTER, 0, MPI_COMM_WORLD, ierr)
if (ierr /= MPI_SUCCESS) then
stop 'Unable to read element_name with MPI'
endif
call MPI_BCAST( element_mass, size(element_mass), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr)
call MPI_BCAST( element_mass, 128, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr)
if (ierr /= MPI_SUCCESS) then
stop 'Unable to read element_name with MPI'
endif

17
src/Utils/LinearAlgebra.irp.f
View File

@ -81,10 +81,11 @@ subroutine ortho_canonical(overlap,LDA,N,C,LDC,m)
call svd(overlap,lda,U,ldc,D,Vt,lda,n,n)
D(:) = dsqrt(D(:))
m=n
do i=1,n
if ( D(i) >= 1.d-6 ) then
D(i) = 1.d0/dsqrt(D(i))
D(i) = 1.d0/D(i)
else
m = i-1
print *, 'Removed Linear dependencies below:', 1.d0/D(m)
@ -101,25 +102,19 @@ subroutine ortho_canonical(overlap,LDA,N,C,LDC,m)
endif
enddo
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP SHARED(S,U,D,Vt,n,C,m) &
!$OMP PRIVATE(i,j)
!$OMP DO
do j=1,n
do i=1,n
S(i,j) = U(i,j)*D(j)
enddo
enddo
do j=1,n
do i=1,n
U(i,j) = C(i,j)
enddo
enddo
!$OMP END DO
!$OMP END PARALLEL
call dgemm('N','N',n,m,n,1.d0,U,size(U,1),S,size(S,1),0.d0,C,size(C,1))
call dgemm('N','N',n,n,n,1.d0,U,size(U,1),S,size(S,1),0.d0,C,size(C,1))
deallocate (U, Vt, D, S)
end

1
src/ZMQ/put_get.irp.f
View File

@ -51,6 +51,7 @@ integer function zmq_get_dvector(zmq_to_qp_run_socket, worker_id, name, x, size_
integer*8 :: rc8
character*(256) :: msg
PROVIDE zmq_state
! Success
zmq_get_dvector = 0