Merge pull request #111 from scemama/master
Accelerated pseudo + bug in Huckel guess
3
configure
vendored
@ -62,8 +62,9 @@ d_dependency = {
|
|||||||
"resultsFile": ["python"],
|
"resultsFile": ["python"],
|
||||||
"emsl": ["python"],
|
"emsl": ["python"],
|
||||||
"gcc": [],
|
"gcc": [],
|
||||||
|
"g++": [],
|
||||||
"python": [],
|
"python": [],
|
||||||
"ninja": ["gcc", "python"],
|
"ninja": ["g++", "python"],
|
||||||
"make": [],
|
"make": [],
|
||||||
"p_graphviz": ["python"]
|
"p_graphviz": ["python"]
|
||||||
}
|
}
|
||||||
|
@ -7,8 +7,6 @@ module Determinants_by_hand : sig
|
|||||||
{ n_int : N_int_number.t;
|
{ n_int : N_int_number.t;
|
||||||
bit_kind : Bit_kind.t;
|
bit_kind : Bit_kind.t;
|
||||||
n_det : Det_number.t;
|
n_det : Det_number.t;
|
||||||
n_states : States_number.t;
|
|
||||||
n_states_diag : States_number.t;
|
|
||||||
expected_s2 : Positive_float.t;
|
expected_s2 : Positive_float.t;
|
||||||
psi_coef : Det_coef.t array;
|
psi_coef : Det_coef.t array;
|
||||||
psi_det : Determinant.t array;
|
psi_det : Determinant.t array;
|
||||||
@ -23,8 +21,6 @@ end = struct
|
|||||||
{ n_int : N_int_number.t;
|
{ n_int : N_int_number.t;
|
||||||
bit_kind : Bit_kind.t;
|
bit_kind : Bit_kind.t;
|
||||||
n_det : Det_number.t;
|
n_det : Det_number.t;
|
||||||
n_states : States_number.t;
|
|
||||||
n_states_diag : States_number.t;
|
|
||||||
expected_s2 : Positive_float.t;
|
expected_s2 : Positive_float.t;
|
||||||
psi_coef : Det_coef.t array;
|
psi_coef : Det_coef.t array;
|
||||||
psi_det : Determinant.t array;
|
psi_det : Determinant.t array;
|
||||||
@ -146,11 +142,12 @@ end = struct
|
|||||||
|> Array.map ~f:Det_coef.of_float
|
|> Array.map ~f:Det_coef.of_float
|
||||||
;;
|
;;
|
||||||
|
|
||||||
let write_psi_coef ~n_det ~n_states c =
|
let write_psi_coef ~n_det c =
|
||||||
let n_det = Det_number.to_int n_det
|
let n_det = Det_number.to_int n_det
|
||||||
and c = Array.to_list c
|
and c = Array.to_list c
|
||||||
|> List.map ~f:Det_coef.to_float
|
|> List.map ~f:Det_coef.to_float
|
||||||
and n_states = States_number.to_int n_states
|
and n_states =
|
||||||
|
read_n_states () |> States_number.to_int
|
||||||
in
|
in
|
||||||
Ezfio.ezfio_array_of_list ~rank:2 ~dim:[| n_det ; n_states |] ~data:c
|
Ezfio.ezfio_array_of_list ~rank:2 ~dim:[| n_det ; n_states |] ~data:c
|
||||||
|> Ezfio.set_determinants_psi_coef
|
|> Ezfio.set_determinants_psi_coef
|
||||||
@ -214,8 +211,6 @@ end = struct
|
|||||||
{ n_int = read_n_int () ;
|
{ n_int = read_n_int () ;
|
||||||
bit_kind = read_bit_kind () ;
|
bit_kind = read_bit_kind () ;
|
||||||
n_det = read_n_det () ;
|
n_det = read_n_det () ;
|
||||||
n_states = read_n_states () ;
|
|
||||||
n_states_diag = read_n_states_diag () ;
|
|
||||||
expected_s2 = read_expected_s2 () ;
|
expected_s2 = read_expected_s2 () ;
|
||||||
psi_coef = read_psi_coef () ;
|
psi_coef = read_psi_coef () ;
|
||||||
psi_det = read_psi_det () ;
|
psi_det = read_psi_det () ;
|
||||||
@ -227,8 +222,6 @@ end = struct
|
|||||||
let write { n_int ;
|
let write { n_int ;
|
||||||
bit_kind ;
|
bit_kind ;
|
||||||
n_det ;
|
n_det ;
|
||||||
n_states ;
|
|
||||||
n_states_diag ;
|
|
||||||
expected_s2 ;
|
expected_s2 ;
|
||||||
psi_coef ;
|
psi_coef ;
|
||||||
psi_det ;
|
psi_det ;
|
||||||
@ -236,10 +229,8 @@ end = struct
|
|||||||
write_n_int n_int ;
|
write_n_int n_int ;
|
||||||
write_bit_kind bit_kind;
|
write_bit_kind bit_kind;
|
||||||
write_n_det n_det;
|
write_n_det n_det;
|
||||||
write_n_states n_states;
|
|
||||||
write_n_states_diag ~n_states:n_states n_states_diag;
|
|
||||||
write_expected_s2 expected_s2;
|
write_expected_s2 expected_s2;
|
||||||
write_psi_coef ~n_det:n_det psi_coef ~n_states:n_states;
|
write_psi_coef ~n_det:n_det psi_coef ;
|
||||||
write_psi_det ~n_int:n_int ~n_det:n_det psi_det;
|
write_psi_det ~n_int:n_int ~n_det:n_det psi_det;
|
||||||
;;
|
;;
|
||||||
|
|
||||||
@ -249,7 +240,7 @@ end = struct
|
|||||||
let mo_tot_num = MO_number.of_int mo_tot_num ~max:mo_tot_num in
|
let mo_tot_num = MO_number.of_int mo_tot_num ~max:mo_tot_num in
|
||||||
let det_text =
|
let det_text =
|
||||||
let nstates =
|
let nstates =
|
||||||
States_number.to_int b.n_states
|
read_n_states () |> States_number.to_int
|
||||||
and ndet =
|
and ndet =
|
||||||
Det_number.to_int b.n_det
|
Det_number.to_int b.n_det
|
||||||
in
|
in
|
||||||
@ -284,12 +275,6 @@ If true, input the expected value of S^2 ::
|
|||||||
|
|
||||||
expected_s2 = %s
|
expected_s2 = %s
|
||||||
|
|
||||||
Number of requested states, and number of states used for the
|
|
||||||
Davidson diagonalization ::
|
|
||||||
|
|
||||||
n_states = %s
|
|
||||||
n_states_diag = %s
|
|
||||||
|
|
||||||
Number of determinants ::
|
Number of determinants ::
|
||||||
|
|
||||||
n_det = %s
|
n_det = %s
|
||||||
@ -299,8 +284,6 @@ Determinants ::
|
|||||||
%s
|
%s
|
||||||
"
|
"
|
||||||
(b.expected_s2 |> Positive_float.to_string)
|
(b.expected_s2 |> Positive_float.to_string)
|
||||||
(b.n_states |> States_number.to_string)
|
|
||||||
(b.n_states_diag |> States_number.to_string)
|
|
||||||
(b.n_det |> Det_number.to_string)
|
(b.n_det |> Det_number.to_string)
|
||||||
det_text
|
det_text
|
||||||
|> Rst_string.of_string
|
|> Rst_string.of_string
|
||||||
@ -313,8 +296,6 @@ Determinants ::
|
|||||||
n_int = %s
|
n_int = %s
|
||||||
bit_kind = %s
|
bit_kind = %s
|
||||||
n_det = %s
|
n_det = %s
|
||||||
n_states = %s
|
|
||||||
n_states_diag = %s
|
|
||||||
expected_s2 = %s
|
expected_s2 = %s
|
||||||
psi_coef = %s
|
psi_coef = %s
|
||||||
psi_det = %s
|
psi_det = %s
|
||||||
@ -322,8 +303,6 @@ psi_det = %s
|
|||||||
(b.n_int |> N_int_number.to_string)
|
(b.n_int |> N_int_number.to_string)
|
||||||
(b.bit_kind |> Bit_kind.to_string)
|
(b.bit_kind |> Bit_kind.to_string)
|
||||||
(b.n_det |> Det_number.to_string)
|
(b.n_det |> Det_number.to_string)
|
||||||
(b.n_states |> States_number.to_string)
|
|
||||||
(b.n_states_diag |> States_number.to_string)
|
|
||||||
(b.expected_s2 |> Positive_float.to_string)
|
(b.expected_s2 |> Positive_float.to_string)
|
||||||
(b.psi_coef |> Array.to_list |> List.map ~f:Det_coef.to_string
|
(b.psi_coef |> Array.to_list |> List.map ~f:Det_coef.to_string
|
||||||
|> String.concat ~sep:", ")
|
|> String.concat ~sep:", ")
|
||||||
|
@ -12,7 +12,6 @@ let rec transpose = function
|
|||||||
;;
|
;;
|
||||||
*)
|
*)
|
||||||
|
|
||||||
|
|
||||||
let input_to_sexp s =
|
let input_to_sexp s =
|
||||||
let result =
|
let result =
|
||||||
String.split_lines s
|
String.split_lines s
|
||||||
|
@ -83,6 +83,12 @@ let input_data = "
|
|||||||
assert (x >= 0.) ;
|
assert (x >= 0.) ;
|
||||||
assert (x <= 1.) ;
|
assert (x <= 1.) ;
|
||||||
|
|
||||||
|
* Energy : float
|
||||||
|
assert (x <=0.) ;
|
||||||
|
|
||||||
|
* S2 : float
|
||||||
|
assert (x >=0.) ;
|
||||||
|
|
||||||
* PT2_energy : float
|
* PT2_energy : float
|
||||||
assert (x >=0.) ;
|
assert (x >=0.) ;
|
||||||
|
|
||||||
|
1
plugins/CISD/.gitignore
vendored
@ -20,6 +20,7 @@ Pseudo
|
|||||||
Selectors_full
|
Selectors_full
|
||||||
SingleRefMethod
|
SingleRefMethod
|
||||||
Utils
|
Utils
|
||||||
|
cisd
|
||||||
cisd_lapack
|
cisd_lapack
|
||||||
ezfio_interface.irp.f
|
ezfio_interface.irp.f
|
||||||
irpf90.make
|
irpf90.make
|
||||||
|
10
plugins/CISD/EZFIO.cfg
Normal file
@ -0,0 +1,10 @@
|
|||||||
|
[energy]
|
||||||
|
type: double precision
|
||||||
|
doc: Variational CISD energy
|
||||||
|
interface: ezfio
|
||||||
|
|
||||||
|
[energy_pt2]
|
||||||
|
type: double precision
|
||||||
|
doc: Estimated CISD energy (including PT2)
|
||||||
|
interface: ezfio
|
||||||
|
|
@ -59,10 +59,6 @@ Documentation
|
|||||||
.. by the `update_README.py` script.
|
.. by the `update_README.py` script.
|
||||||
|
|
||||||
|
|
||||||
`cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/cisd_lapack.irp.f#L1>`_
|
|
||||||
Undocumented
|
|
||||||
|
|
||||||
|
|
||||||
`h_apply_cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/H_apply.irp.f_shell_8#L414>`_
|
`h_apply_cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/H_apply.irp.f_shell_8#L414>`_
|
||||||
Calls H_apply on the HF determinant and selects all connected single and double
|
Calls H_apply on the HF determinant and selects all connected single and double
|
||||||
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
|
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
|
||||||
|
Before Width: | Height: | Size: 87 KiB After Width: | Height: | Size: 84 KiB |
Before Width: | Height: | Size: 64 KiB After Width: | Height: | Size: 59 KiB |
@ -23,6 +23,10 @@ Documentation
|
|||||||
.. by the `update_README.py` script.
|
.. by the `update_README.py` script.
|
||||||
|
|
||||||
|
|
||||||
|
`e_curve <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/e_curve.irp.f#L1>`_
|
||||||
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`full_ci <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/full_ci_no_skip.irp.f#L1>`_
|
`full_ci <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/full_ci_no_skip.irp.f#L1>`_
|
||||||
Undocumented
|
Undocumented
|
||||||
|
|
||||||
|
@ -27,6 +27,8 @@ program var_pt2_ratio_run
|
|||||||
call diagonalize_CI
|
call diagonalize_CI
|
||||||
ratio = (CI_energy(1) - HF_energy) / (CI_energy(1)+pt2(1) - HF_energy)
|
ratio = (CI_energy(1) - HF_energy) / (CI_energy(1)+pt2(1) - HF_energy)
|
||||||
if (N_det > 20000) then
|
if (N_det > 20000) then
|
||||||
|
N_det = 20000
|
||||||
|
TOUCH N_det
|
||||||
exit
|
exit
|
||||||
endif
|
endif
|
||||||
enddo
|
enddo
|
||||||
|
Before Width: | Height: | Size: 100 KiB After Width: | Height: | Size: 98 KiB |
Before Width: | Height: | Size: 74 KiB After Width: | Height: | Size: 73 KiB |
1
plugins/Hartree_Fock/.gitignore
vendored
@ -21,4 +21,5 @@ Utils
|
|||||||
ezfio_interface.irp.f
|
ezfio_interface.irp.f
|
||||||
irpf90.make
|
irpf90.make
|
||||||
irpf90_entities
|
irpf90_entities
|
||||||
|
simple_SCF
|
||||||
tags
|
tags
|
@ -161,6 +161,10 @@ Documentation
|
|||||||
optional: mo_basis.mo_coef
|
optional: mo_basis.mo_coef
|
||||||
|
|
||||||
|
|
||||||
|
`simple_scf <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/simple_SCF.irp.f#L1>`_
|
||||||
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`thresh_scf <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L46>`_
|
`thresh_scf <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L46>`_
|
||||||
Threshold on the convergence of the Hartree Fock energy
|
Threshold on the convergence of the Hartree Fock energy
|
||||||
|
|
||||||
|
@ -73,19 +73,19 @@ BEGIN_PROVIDER [double precision, diagonal_Fock_matrix_mo_sum, (mo_tot_num)]
|
|||||||
END_DOC
|
END_DOC
|
||||||
integer :: i,j
|
integer :: i,j
|
||||||
double precision :: accu
|
double precision :: accu
|
||||||
do i = 1,elec_alpha_num
|
do j = 1,elec_alpha_num
|
||||||
accu = 0.d0
|
accu = 0.d0
|
||||||
do j = 1, elec_alpha_num
|
do i = 1, elec_alpha_num
|
||||||
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
|
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
|
||||||
enddo
|
enddo
|
||||||
diagonal_Fock_matrix_mo_sum(i) = accu + mo_mono_elec_integral(i,i)
|
diagonal_Fock_matrix_mo_sum(j) = accu + mo_mono_elec_integral(j,j)
|
||||||
enddo
|
enddo
|
||||||
do i = elec_alpha_num+1,mo_tot_num
|
do j = elec_alpha_num+1,mo_tot_num
|
||||||
accu = 0.d0
|
accu = 0.d0
|
||||||
do j = 1, elec_alpha_num
|
do i = 1, elec_alpha_num
|
||||||
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
|
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
|
||||||
enddo
|
enddo
|
||||||
diagonal_Fock_matrix_mo_sum(i) = accu + mo_mono_elec_integral(i,i)
|
diagonal_Fock_matrix_mo_sum(j) = accu + mo_mono_elec_integral(j,j)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
@ -4,7 +4,7 @@ subroutine huckel_guess
|
|||||||
! Build the MOs using the extended Huckel model
|
! Build the MOs using the extended Huckel model
|
||||||
END_DOC
|
END_DOC
|
||||||
integer :: i,j
|
integer :: i,j
|
||||||
double precision :: tmp_matrix(ao_num_align,ao_num),accu
|
double precision :: accu
|
||||||
double precision :: c
|
double precision :: c
|
||||||
character*(64) :: label
|
character*(64) :: label
|
||||||
|
|
||||||
@ -19,14 +19,12 @@ subroutine huckel_guess
|
|||||||
c = 0.5d0 * 1.75d0
|
c = 0.5d0 * 1.75d0
|
||||||
|
|
||||||
do j=1,ao_num
|
do j=1,ao_num
|
||||||
|
!DIR$ VECTOR ALIGNED
|
||||||
do i=1,ao_num
|
do i=1,ao_num
|
||||||
if (i.ne.j) then
|
Fock_matrix_ao(i,j) = c*ao_overlap(i,j)*(ao_mono_elec_integral_diag(i) + &
|
||||||
Fock_matrix_ao(i,j) = c*ao_overlap(i,j)*(ao_mono_elec_integral(i,i) + &
|
ao_mono_elec_integral_diag(j))
|
||||||
ao_mono_elec_integral(j,j))
|
|
||||||
else
|
|
||||||
Fock_matrix_ao(i,j) = Fock_matrix_alpha_ao(i,j)
|
|
||||||
endif
|
|
||||||
enddo
|
enddo
|
||||||
|
Fock_matrix_ao(j,j) = Fock_matrix_alpha_ao(j,j)
|
||||||
enddo
|
enddo
|
||||||
TOUCH Fock_matrix_ao
|
TOUCH Fock_matrix_ao
|
||||||
mo_coef = eigenvectors_fock_matrix_mo
|
mo_coef = eigenvectors_fock_matrix_mo
|
||||||
|
Before Width: | Height: | Size: 59 KiB After Width: | Height: | Size: 57 KiB |
Before Width: | Height: | Size: 110 KiB After Width: | Height: | Size: 106 KiB |
@ -107,92 +107,156 @@ Documentation
|
|||||||
Undocumented
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`perturb_buffer_by_mono_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L161>`_
|
`perturb_buffer_by_mono_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L791>`_
|
||||||
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L686>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L266>`_
|
`perturb_buffer_by_mono_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L266>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L476>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L371>`_
|
`perturb_buffer_by_mono_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L371>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L581>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L476>`_
|
`perturb_buffer_by_mono_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L476>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L371>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L791>`_
|
`perturb_buffer_by_mono_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L791>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L266>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L686>`_
|
`perturb_buffer_by_mono_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L686>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L161>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L581>`_
|
`perturb_buffer_by_mono_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L581>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L56>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L896>`_
|
`perturb_buffer_by_mono_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L896>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_by_mono_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L896>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_by_mono_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L56>`_
|
`perturb_buffer_by_mono_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L56>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
`perturb_buffer_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L110>`_
|
`perturb_buffer_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L740>`_
|
||||||
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L635>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L215>`_
|
`perturb_buffer_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L215>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L425>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L320>`_
|
`perturb_buffer_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L320>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L530>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L425>`_
|
`perturb_buffer_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L425>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L320>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L740>`_
|
`perturb_buffer_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L740>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L215>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L635>`_
|
`perturb_buffer_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L635>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L110>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L530>`_
|
`perturb_buffer_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L530>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L5>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L845>`_
|
`perturb_buffer_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L845>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`perturb_buffer_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L845>`_
|
||||||
|
=======
|
||||||
`perturb_buffer_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L5>`_
|
`perturb_buffer_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L5>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
|
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
|
||||||
routine.
|
routine.
|
||||||
|
|
||||||
|
Before Width: | Height: | Size: 78 KiB After Width: | Height: | Size: 76 KiB |
Before Width: | Height: | Size: 62 KiB After Width: | Height: | Size: 60 KiB |
Before Width: | Height: | Size: 7.6 KiB After Width: | Height: | Size: 7.4 KiB |
Before Width: | Height: | Size: 3.5 KiB After Width: | Height: | Size: 3.4 KiB |
3
plugins/QmcChem/.gitignore
vendored
@ -16,8 +16,11 @@ Makefile.depend
|
|||||||
Nuclei
|
Nuclei
|
||||||
Pseudo
|
Pseudo
|
||||||
Utils
|
Utils
|
||||||
|
e_curve_qmc
|
||||||
ezfio_interface.irp.f
|
ezfio_interface.irp.f
|
||||||
irpf90.make
|
irpf90.make
|
||||||
irpf90_entities
|
irpf90_entities
|
||||||
save_for_qmcchem
|
save_for_qmcchem
|
||||||
tags
|
tags
|
||||||
|
target_pt2_qmc
|
||||||
|
test
|
@ -66,6 +66,14 @@ Documentation
|
|||||||
title="f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C" />
|
title="f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C" />
|
||||||
|
|
||||||
|
|
||||||
|
`compute_energy <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/target_pt2_qmc.irp.f#L80>`_
|
||||||
|
Compute an energy when a threshold is applied
|
||||||
|
|
||||||
|
|
||||||
|
`e_curve <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/test.irp.f#L1>`_
|
||||||
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`mo_pseudo_grid <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/pot_ao_pseudo_ints.irp.f#L56>`_
|
`mo_pseudo_grid <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/pot_ao_pseudo_ints.irp.f#L56>`_
|
||||||
Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \phi_i^{A} (r-r_A) d\Omega_C
|
Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \phi_i^{A} (r-r_A) d\Omega_C
|
||||||
.br
|
.br
|
||||||
|
121
plugins/QmcChem/target_pt2_qmc.irp.f
Normal file
@ -0,0 +1,121 @@
|
|||||||
|
program e_curve
|
||||||
|
use bitmasks
|
||||||
|
implicit none
|
||||||
|
integer :: i,j,k, nab, m, l, n_up, n_dn, n
|
||||||
|
double precision :: norm, E, hij, num, ci, cj
|
||||||
|
integer, allocatable :: iorder(:)
|
||||||
|
double precision , allocatable :: norm_sort(:), psi_bilinear_matrix_values_save(:)
|
||||||
|
nab = n_det_alpha_unique+n_det_beta_unique
|
||||||
|
|
||||||
|
allocate ( norm_sort(0:nab), iorder(0:nab), psi_bilinear_matrix_values_save(N_det) )
|
||||||
|
|
||||||
|
|
||||||
|
norm_sort(0) = 0.d0
|
||||||
|
iorder(0) = 0
|
||||||
|
do i=1,n_det_alpha_unique
|
||||||
|
norm_sort(i) = det_alpha_norm(i)
|
||||||
|
iorder(i) = i
|
||||||
|
enddo
|
||||||
|
|
||||||
|
do i=1,n_det_beta_unique
|
||||||
|
norm_sort(i+n_det_alpha_unique) = det_beta_norm(i)
|
||||||
|
iorder(i+n_det_alpha_unique) = -i
|
||||||
|
enddo
|
||||||
|
|
||||||
|
call dsort(norm_sort(1),iorder(1),nab)
|
||||||
|
|
||||||
|
if (.not.read_wf) then
|
||||||
|
stop 'Please set read_wf to true'
|
||||||
|
endif
|
||||||
|
|
||||||
|
psi_bilinear_matrix_values_save = psi_bilinear_matrix_values(:,1)
|
||||||
|
print *, '=========================================================='
|
||||||
|
print '(A8,2X,A8,2X,A12,2X,A10,2X,A12)', 'Thresh.', 'Ndet', 'Cost', 'Norm', 'E'
|
||||||
|
print *, '=========================================================='
|
||||||
|
integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
|
||||||
|
|
||||||
|
double precision :: thresh, E_min, E_max, E_prev
|
||||||
|
thresh = 0.d0
|
||||||
|
call compute_energy(psi_bilinear_matrix_values_save,E_max,m,norm)
|
||||||
|
call i_h_j(psi_det_sorted(1,1,1), psi_det_sorted(1,1,1), N_int, E_min)
|
||||||
|
print *, E_min, E_max
|
||||||
|
|
||||||
|
n_up = nab
|
||||||
|
n_dn = 0
|
||||||
|
do while (n_up > n_dn)
|
||||||
|
n = n_dn + (n_up-n_dn)/2
|
||||||
|
psi_bilinear_matrix_values_save = psi_bilinear_matrix_values(:,1)
|
||||||
|
do j=1,n
|
||||||
|
i = iorder(j)
|
||||||
|
if (i<0) then
|
||||||
|
do k=1,n_det
|
||||||
|
if (psi_bilinear_matrix_columns(k) == -i) then
|
||||||
|
psi_bilinear_matrix_values_save(k) = 0.d0
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
else
|
||||||
|
do k=1,n_det
|
||||||
|
if (psi_bilinear_matrix_rows(k) == i) then
|
||||||
|
psi_bilinear_matrix_values_save(k) = 0.d0
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
call compute_energy(psi_bilinear_matrix_values_save,E,m,norm)
|
||||||
|
print '(E9.1,2X,I8,2X,F10.2,2X,F10.8,2X,F12.6)', norm_sort(n), m, &
|
||||||
|
dble( elec_alpha_num**3 + elec_alpha_num**2 * m ) / &
|
||||||
|
dble( elec_alpha_num**3 + elec_alpha_num**2 * n ), norm, E
|
||||||
|
if (E < target_energy) then
|
||||||
|
n_dn = n+1
|
||||||
|
else
|
||||||
|
n_up = n
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
print *, '=========================================================='
|
||||||
|
print *, norm_sort(n), target_energy
|
||||||
|
|
||||||
|
deallocate (iorder, norm_sort, psi_bilinear_matrix_values_save)
|
||||||
|
end
|
||||||
|
|
||||||
|
subroutine compute_energy(psi_bilinear_matrix_values_save, E, m, norm)
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Compute an energy when a threshold is applied
|
||||||
|
END_DOC
|
||||||
|
double precision, intent(in) :: psi_bilinear_matrix_values_save(n_det)
|
||||||
|
integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
|
||||||
|
integer :: i,j, k, l, m
|
||||||
|
double precision :: num, norm, ci, cj, hij, E
|
||||||
|
|
||||||
|
|
||||||
|
num = 0.d0
|
||||||
|
norm = 0.d0
|
||||||
|
m = 0
|
||||||
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(k,l,det_i,det_j,ci,cj,hij) REDUCTION(+:norm,m,num)
|
||||||
|
allocate( det_i(N_int,2), det_j(N_int,2))
|
||||||
|
!$OMP DO
|
||||||
|
do k=1,n_det
|
||||||
|
if (psi_bilinear_matrix_values_save(k) == 0.d0) then
|
||||||
|
cycle
|
||||||
|
endif
|
||||||
|
ci = psi_bilinear_matrix_values_save(k)
|
||||||
|
det_i(:,1) = psi_det_alpha_unique(:,psi_bilinear_matrix_rows(k))
|
||||||
|
det_i(:,2) = psi_det_beta_unique(:,psi_bilinear_matrix_columns(k))
|
||||||
|
do l=1,n_det
|
||||||
|
if (psi_bilinear_matrix_values_save(l) == 0.d0) then
|
||||||
|
cycle
|
||||||
|
endif
|
||||||
|
cj = psi_bilinear_matrix_values_save(l)
|
||||||
|
det_j(:,1) = psi_det_alpha_unique(:,psi_bilinear_matrix_rows(l))
|
||||||
|
det_j(:,2) = psi_det_beta_unique(:,psi_bilinear_matrix_columns(l))
|
||||||
|
call i_h_j(det_i, det_j, N_int, hij)
|
||||||
|
num = num + ci*cj*hij
|
||||||
|
enddo
|
||||||
|
norm = norm + ci*ci
|
||||||
|
m = m+1
|
||||||
|
enddo
|
||||||
|
!$OMP END DO
|
||||||
|
deallocate (det_i,det_j)
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
E = num / norm + nuclear_repulsion
|
||||||
|
end
|
Before Width: | Height: | Size: 74 KiB After Width: | Height: | Size: 72 KiB |
Before Width: | Height: | Size: 33 KiB After Width: | Height: | Size: 33 KiB |
@ -707,7 +707,7 @@ def ninja_dot_tree_rule():
|
|||||||
l_string = [
|
l_string = [
|
||||||
"rule build_dot_tree", " command = {0}".format(" ; ".join(l_cmd)),
|
"rule build_dot_tree", " command = {0}".format(" ; ".join(l_cmd)),
|
||||||
" generator = 1",
|
" generator = 1",
|
||||||
" description = Generate Png representtion of the Tree Dependencies of $module_rel",
|
" description = Generating Png representation of the Tree Dependencies of $module_rel",
|
||||||
""
|
""
|
||||||
]
|
]
|
||||||
|
|
||||||
|
@ -13,7 +13,6 @@ then
|
|||||||
fi
|
fi
|
||||||
source ${QP_ROOT}/scripts/qp_include.sh
|
source ${QP_ROOT}/scripts/qp_include.sh
|
||||||
|
|
||||||
|
|
||||||
function do_gitingore()
|
function do_gitingore()
|
||||||
{
|
{
|
||||||
cat << EOF > .gitignore
|
cat << EOF > .gitignore
|
||||||
|
@ -180,6 +180,11 @@ class ModuleHandler():
|
|||||||
def create_png(self, l_module):
|
def create_png(self, l_module):
|
||||||
"""Create the png of the dependency tree for a l_module"""
|
"""Create the png of the dependency tree for a l_module"""
|
||||||
|
|
||||||
|
# Don't update if we are not in the main repository
|
||||||
|
from is_master_repository import is_master_repository
|
||||||
|
if not is_master_repository:
|
||||||
|
return
|
||||||
|
|
||||||
basename = "tree_dependency"
|
basename = "tree_dependency"
|
||||||
path = '{0}.png'.format(basename)
|
path = '{0}.png'.format(basename)
|
||||||
|
|
||||||
@ -289,6 +294,12 @@ if __name__ == '__main__':
|
|||||||
pass
|
pass
|
||||||
|
|
||||||
if arguments["create_git_ignore"]:
|
if arguments["create_git_ignore"]:
|
||||||
|
|
||||||
|
# Don't update if we are not in the main repository
|
||||||
|
from is_master_repository import is_master_repository
|
||||||
|
if not is_master_repository:
|
||||||
|
sys.exit()
|
||||||
|
|
||||||
path = os.path.join(module_abs, ".gitignore")
|
path = os.path.join(module_abs, ".gitignore")
|
||||||
|
|
||||||
with open(path, "w+") as f:
|
with open(path, "w+") as f:
|
||||||
|
@ -168,7 +168,15 @@ def update_documentation(d_readmen, root_module):
|
|||||||
|
|
||||||
d_readme[path]["documentation"] = "\n".join(l_doc_section)
|
d_readme[path]["documentation"] = "\n".join(l_doc_section)
|
||||||
|
|
||||||
|
|
||||||
if __name__ == '__main__':
|
if __name__ == '__main__':
|
||||||
|
|
||||||
|
# Update documentation only if the remote repository is
|
||||||
|
# the main repository
|
||||||
|
from is_master_repository import is_master_repository
|
||||||
|
if not is_master_repository:
|
||||||
|
sys.exit(0)
|
||||||
|
|
||||||
arguments = docopt(__doc__)
|
arguments = docopt(__doc__)
|
||||||
|
|
||||||
if arguments["--root_module"]:
|
if arguments["--root_module"]:
|
||||||
|
14
scripts/utility/is_master_repository.py
Executable file
@ -0,0 +1,14 @@
|
|||||||
|
#!/usr/bin/env python
|
||||||
|
|
||||||
|
import subprocess
|
||||||
|
pipe = subprocess.Popen("git config --local --get remote.origin.url", \
|
||||||
|
shell=True, stdout=subprocess.PIPE)
|
||||||
|
result = pipe.stdout.read()
|
||||||
|
is_master_repository = "LCPQ/quantum_package" in result
|
||||||
|
|
||||||
|
if __name__ == "__main__":
|
||||||
|
import sys
|
||||||
|
if is_master_repository:
|
||||||
|
sys.exit(0)
|
||||||
|
else:
|
||||||
|
sys.exit(-1)
|
Before Width: | Height: | Size: 14 KiB After Width: | Height: | Size: 14 KiB |
Before Width: | Height: | Size: 28 KiB After Width: | Height: | Size: 28 KiB |
@ -53,9 +53,10 @@ interface: ezfio,provider,ocaml
|
|||||||
default: 0.999
|
default: 0.999
|
||||||
|
|
||||||
[n_states_diag]
|
[n_states_diag]
|
||||||
type: integer
|
type: States_number
|
||||||
doc: n_states_diag
|
doc: n_states_diag
|
||||||
interface: ezfio,provider
|
default: 1
|
||||||
|
interface: ezfio,provider,ocaml
|
||||||
|
|
||||||
[n_int]
|
[n_int]
|
||||||
interface: ezfio
|
interface: ezfio
|
||||||
@ -89,24 +90,25 @@ doc: psi_det
|
|||||||
type: integer*8
|
type: integer*8
|
||||||
size: (determinants.n_int*determinants.bit_kind/8,2,determinants.n_det)
|
size: (determinants.n_int*determinants.bit_kind/8,2,determinants.n_det)
|
||||||
|
|
||||||
[det_num]
|
|
||||||
interface: ezfio,provider
|
|
||||||
doc: det_num
|
|
||||||
type: integer
|
|
||||||
|
|
||||||
[det_occ]
|
[det_occ]
|
||||||
interface: ezfio,provider
|
interface: ezfio,provider
|
||||||
doc: det_occ
|
doc: det_occ
|
||||||
type: integer
|
type: integer
|
||||||
size: (electrons.elec_alpha_num,determinants.det_num,2)
|
size: (electrons.elec_alpha_num,determinants.n_det,2)
|
||||||
|
|
||||||
[det_coef]
|
[det_coef]
|
||||||
interface: ezfio,provider
|
interface: ezfio,provider
|
||||||
doc: det_coef
|
doc: det_coef
|
||||||
type: double precision
|
type: double precision
|
||||||
size: (determinants.det_num)
|
size: (determinants.n_det)
|
||||||
|
|
||||||
[expected_s2]
|
[expected_s2]
|
||||||
interface: ezfio,provider
|
interface: ezfio,provider
|
||||||
doc: expcted_s2
|
doc: Expected value of S^2
|
||||||
type: double precision
|
type: double precision
|
||||||
|
|
||||||
|
[target_energy]
|
||||||
|
interface: ezfio,provider,ocaml
|
||||||
|
doc: Energy that should be obtained when truncating the wave function (optional)
|
||||||
|
type: Energy
|
||||||
|
default: 0.
|
||||||
|
@ -113,7 +113,7 @@ Documentation
|
|||||||
After calling this subroutine, N_det, psi_det and psi_coef need to be touched
|
After calling this subroutine, N_det, psi_det and psi_coef need to be touched
|
||||||
|
|
||||||
|
|
||||||
`create_wf_of_psi_bilinear_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L417>`_
|
`create_wf_of_psi_bilinear_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L445>`_
|
||||||
Generate a wave function containing all possible products
|
Generate a wave function containing all possible products
|
||||||
of alpha and beta determinants
|
of alpha and beta determinants
|
||||||
|
|
||||||
@ -186,6 +186,18 @@ Documentation
|
|||||||
degree : Degree of excitation
|
degree : Degree of excitation
|
||||||
|
|
||||||
|
|
||||||
|
`det_alpha_norm <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L353>`_
|
||||||
|
Norm of the alpha and beta spin determinants in the wave function:
|
||||||
|
.br
|
||||||
|
||Da||_i \sum_j C_{ij}**2
|
||||||
|
|
||||||
|
|
||||||
|
`det_beta_norm <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L354>`_
|
||||||
|
Norm of the alpha and beta spin determinants in the wave function:
|
||||||
|
.br
|
||||||
|
||Da||_i \sum_j C_{ij}**2
|
||||||
|
|
||||||
|
|
||||||
`det_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L138>`_
|
`det_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L138>`_
|
||||||
det_coef
|
det_coef
|
||||||
|
|
||||||
@ -194,10 +206,6 @@ Documentation
|
|||||||
Undocumented
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`det_num <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L248>`_
|
|
||||||
det_num
|
|
||||||
|
|
||||||
|
|
||||||
`det_occ <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L226>`_
|
`det_occ <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L226>`_
|
||||||
det_occ
|
det_occ
|
||||||
|
|
||||||
@ -306,7 +314,21 @@ Documentation
|
|||||||
to repeat the excitations
|
to repeat the excitations
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`filter_connected_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/filter_connected.irp.f#L101>`_
|
||||||
|
Filters out the determinants that are not connected by H
|
||||||
|
returns the array idx which contains the index of the
|
||||||
|
determinants in the array key1 that interact
|
||||||
|
via the H operator with key2.
|
||||||
|
idx(0) is the number of determinants that interact with key1
|
||||||
|
.br
|
||||||
|
Determinants are taken from the psi_det_sorted_ab array
|
||||||
|
|
||||||
|
|
||||||
|
`generate_all_alpha_beta_det_products <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L500>`_
|
||||||
|
=======
|
||||||
`generate_all_alpha_beta_det_products <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L472>`_
|
`generate_all_alpha_beta_det_products <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L472>`_
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
Create a wave function from all possible alpha x beta determinants
|
Create a wave function from all possible alpha x beta determinants
|
||||||
|
|
||||||
|
|
||||||
@ -588,22 +610,22 @@ Documentation
|
|||||||
Wave function sorted by determinants contribution to the norm (state-averaged)
|
Wave function sorted by determinants contribution to the norm (state-averaged)
|
||||||
|
|
||||||
|
|
||||||
`psi_bilinear_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L400>`_
|
`psi_bilinear_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L428>`_
|
||||||
Coefficient matrix if the wave function is expressed in a bilinear form :
|
Coefficient matrix if the wave function is expressed in a bilinear form :
|
||||||
D_a^t C D_b
|
D_a^t C D_b
|
||||||
|
|
||||||
|
|
||||||
`psi_bilinear_matrix_columns <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L362>`_
|
`psi_bilinear_matrix_columns <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L390>`_
|
||||||
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
||||||
D_a^t C D_b
|
D_a^t C D_b
|
||||||
|
|
||||||
|
|
||||||
`psi_bilinear_matrix_rows <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L361>`_
|
`psi_bilinear_matrix_rows <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L389>`_
|
||||||
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
||||||
D_a^t C D_b
|
D_a^t C D_b
|
||||||
|
|
||||||
|
|
||||||
`psi_bilinear_matrix_values <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L360>`_
|
`psi_bilinear_matrix_values <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L388>`_
|
||||||
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
Sparse coefficient matrix if the wave function is expressed in a bilinear form :
|
||||||
D_a^t C D_b
|
D_a^t C D_b
|
||||||
|
|
||||||
@ -757,7 +779,7 @@ Documentation
|
|||||||
Reads the determinants from the EZFIO file
|
Reads the determinants from the EZFIO file
|
||||||
|
|
||||||
|
|
||||||
`read_wf <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L160>`_
|
`read_wf <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L116>`_
|
||||||
If true, read the wave function from the EZFIO file
|
If true, read the wave function from the EZFIO file
|
||||||
|
|
||||||
|
|
||||||
@ -782,7 +804,7 @@ Documentation
|
|||||||
Undocumented
|
Undocumented
|
||||||
|
|
||||||
|
|
||||||
`s2_eig <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L116>`_
|
`s2_eig <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L248>`_
|
||||||
Force the wave function to be an eigenfunction of S^2
|
Force the wave function to be an eigenfunction of S^2
|
||||||
|
|
||||||
|
|
||||||
@ -869,8 +891,13 @@ Documentation
|
|||||||
Weights in the state-average calculation of the density matrix
|
Weights in the state-average calculation of the density matrix
|
||||||
|
|
||||||
|
|
||||||
|
<<<<<<< HEAD
|
||||||
|
`target_energy <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/ezfio_interface.irp.f#L160>`_
|
||||||
|
Energy that should be obtained when truncating the wave function (optional)
|
||||||
|
=======
|
||||||
`tamiser <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/davidson.irp.f#L91>`_
|
`tamiser <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/davidson.irp.f#L91>`_
|
||||||
Undocumented
|
Undocumented
|
||||||
|
>>>>>>> 2d3ba8003b05cb406674cc93a9bfc917d94db7fc
|
||||||
|
|
||||||
|
|
||||||
`threshold_convergence_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_SC2.irp.f#L18>`_
|
`threshold_convergence_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_SC2.irp.f#L18>`_
|
||||||
|
@ -350,6 +350,34 @@ subroutine write_spindeterminants
|
|||||||
|
|
||||||
end
|
end
|
||||||
|
|
||||||
|
BEGIN_PROVIDER [ double precision, det_alpha_norm, (N_det_alpha_unique) ]
|
||||||
|
&BEGIN_PROVIDER [ double precision, det_beta_norm, (N_det_beta_unique) ]
|
||||||
|
implicit none
|
||||||
|
BEGIN_DOC
|
||||||
|
! Norm of the alpha and beta spin determinants in the wave function:
|
||||||
|
!
|
||||||
|
! ||Da||_i \sum_j C_{ij}**2
|
||||||
|
END_DOC
|
||||||
|
|
||||||
|
integer :: i,j,k,l
|
||||||
|
double precision :: f
|
||||||
|
|
||||||
|
det_alpha_norm = 0.d0
|
||||||
|
det_beta_norm = 0.d0
|
||||||
|
do k=1,N_det
|
||||||
|
i = psi_bilinear_matrix_rows(k)
|
||||||
|
j = psi_bilinear_matrix_columns(k)
|
||||||
|
do l=1,N_states
|
||||||
|
f = psi_bilinear_matrix_values(k,l)*psi_bilinear_matrix_values(k,l)
|
||||||
|
enddo
|
||||||
|
det_alpha_norm(i) += f
|
||||||
|
det_beta_norm(j) += f
|
||||||
|
enddo
|
||||||
|
det_alpha_norm = det_alpha_norm / dble(N_states)
|
||||||
|
det_beta_norm = det_beta_norm / dble(N_states)
|
||||||
|
|
||||||
|
END_PROVIDER
|
||||||
|
|
||||||
|
|
||||||
!==============================================================================!
|
!==============================================================================!
|
||||||
! !
|
! !
|
||||||
|
Before Width: | Height: | Size: 57 KiB After Width: | Height: | Size: 55 KiB |
Before Width: | Height: | Size: 7.1 KiB After Width: | Height: | Size: 7.0 KiB |
Before Width: | Height: | Size: 3.4 KiB After Width: | Height: | Size: 3.3 KiB |
@ -481,7 +481,7 @@ IRP_ENDIF COARRAY
|
|||||||
ao_bielec_integrals_in_map = .True.
|
ao_bielec_integrals_in_map = .True.
|
||||||
if (write_ao_integrals) then
|
if (write_ao_integrals) then
|
||||||
call dump_ao_integrals(trim(ezfio_filename)//'/work/ao_integrals.bin')
|
call dump_ao_integrals(trim(ezfio_filename)//'/work/ao_integrals.bin')
|
||||||
call ezfio_set_integrals_bielec_disk_access_ao_integrals(.True.)
|
call ezfio_set_integrals_bielec_disk_access_ao_integrals("Read")
|
||||||
endif
|
endif
|
||||||
|
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
@ -312,7 +312,7 @@ IRP_ENDIF
|
|||||||
|
|
||||||
if (write_mo_integrals) then
|
if (write_mo_integrals) then
|
||||||
call dump_mo_integrals(trim(ezfio_filename)//'/work/mo_integrals.bin')
|
call dump_mo_integrals(trim(ezfio_filename)//'/work/mo_integrals.bin')
|
||||||
call ezfio_set_integrals_bielec_disk_access_mo_integrals(.True.)
|
call ezfio_set_integrals_bielec_disk_access_mo_integrals("Read")
|
||||||
endif
|
endif
|
||||||
|
|
||||||
|
|
||||||
|
Before Width: | Height: | Size: 41 KiB After Width: | Height: | Size: 40 KiB |
@ -93,6 +93,11 @@ Documentation
|
|||||||
: sum of the kinetic and nuclear electronic potential
|
: sum of the kinetic and nuclear electronic potential
|
||||||
|
|
||||||
|
|
||||||
|
`ao_mono_elec_integral_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Monoelec/ao_mono_ints.irp.f#L2>`_
|
||||||
|
array of the mono electronic hamiltonian on the AOs basis
|
||||||
|
: sum of the kinetic and nuclear electronic potential
|
||||||
|
|
||||||
|
|
||||||
`ao_nucl_elec_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Monoelec/pot_ao_ints.irp.f#L1>`_
|
`ao_nucl_elec_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Monoelec/pot_ao_ints.irp.f#L1>`_
|
||||||
interaction nuclear electron
|
interaction nuclear electron
|
||||||
|
|
||||||
|
@ -1,4 +1,5 @@
|
|||||||
BEGIN_PROVIDER [ double precision, ao_mono_elec_integral,(ao_num_align,ao_num)]
|
BEGIN_PROVIDER [ double precision, ao_mono_elec_integral,(ao_num_align,ao_num)]
|
||||||
|
&BEGIN_PROVIDER [ double precision, ao_mono_elec_integral_diag,(ao_num)]
|
||||||
implicit none
|
implicit none
|
||||||
integer :: i,j,n,l
|
integer :: i,j,n,l
|
||||||
BEGIN_DOC
|
BEGIN_DOC
|
||||||
@ -6,9 +7,11 @@ BEGIN_PROVIDER [ double precision, ao_mono_elec_integral,(ao_num_align,ao_num)]
|
|||||||
! : sum of the kinetic and nuclear electronic potential
|
! : sum of the kinetic and nuclear electronic potential
|
||||||
END_DOC
|
END_DOC
|
||||||
do j = 1, ao_num
|
do j = 1, ao_num
|
||||||
|
!DIR$ VECTOR ALIGNED
|
||||||
do i = 1, ao_num
|
do i = 1, ao_num
|
||||||
ao_mono_elec_integral(i,j) = ao_nucl_elec_integral(i,j) + ao_kinetic_integral(i,j) + ao_pseudo_integral(i,j)
|
ao_mono_elec_integral(i,j) = ao_nucl_elec_integral(i,j) + ao_kinetic_integral(i,j) + ao_pseudo_integral(i,j)
|
||||||
enddo
|
enddo
|
||||||
|
ao_mono_elec_integral_diag(j) = ao_mono_elec_integral(j,j)
|
||||||
enddo
|
enddo
|
||||||
END_PROVIDER
|
END_PROVIDER
|
||||||
|
|
||||||
|
@ -214,13 +214,14 @@ integer :: l,k, nkl_max
|
|||||||
! |_) | (_| (_| | | (_| \/
|
! |_) | (_| (_| | | (_| \/
|
||||||
! _| /
|
! _| /
|
||||||
|
|
||||||
double precision, allocatable :: array_coefs_A(:,:,:)
|
double precision, allocatable :: array_coefs_A(:,:)
|
||||||
double precision, allocatable :: array_coefs_B(:,:,:)
|
double precision, allocatable :: array_coefs_B(:,:)
|
||||||
|
|
||||||
double precision, allocatable :: array_R(:,:,:,:,:)
|
double precision, allocatable :: array_R(:,:,:,:,:)
|
||||||
double precision, allocatable :: array_I_A(:,:,:,:,:)
|
double precision, allocatable :: array_I_A(:,:,:,:,:)
|
||||||
double precision, allocatable :: array_I_B(:,:,:,:,:)
|
double precision, allocatable :: array_I_B(:,:,:,:,:)
|
||||||
|
|
||||||
|
double precision :: f1, f2, f3
|
||||||
|
|
||||||
! _
|
! _
|
||||||
! / _. | _ |
|
! / _. | _ |
|
||||||
@ -255,14 +256,14 @@ nkl_max=4
|
|||||||
! A l l o c a t e !
|
! A l l o c a t e !
|
||||||
!=!=!=!=!=!=!=!=!=!
|
!=!=!=!=!=!=!=!=!=!
|
||||||
|
|
||||||
allocate (array_coefs_A(0:ntot,0:ntot,0:ntot))
|
allocate (array_coefs_A(0:ntot,3))
|
||||||
allocate (array_coefs_B(0:ntot,0:ntot,0:ntot))
|
allocate (array_coefs_B(0:ntot,3))
|
||||||
|
|
||||||
allocate (array_R(0:ntot+nkl_max,kmax,0:lmax,0:lmax+ntot,0:lmax+ntot))
|
allocate (array_R(kmax,0:ntot+nkl_max,0:lmax,0:lmax+ntot,0:lmax+ntot))
|
||||||
|
|
||||||
allocate (array_I_A(0:lmax+ntot,-(lmax+ntot):lmax+ntot,0:ntot,0:ntot,0:ntot))
|
allocate (array_I_A(-(lmax+ntot):lmax+ntot,0:lmax+ntot,0:ntot,0:ntot,0:ntot))
|
||||||
|
|
||||||
allocate (array_I_B(0:lmax+ntot,-(lmax+ntot):lmax+ntot,0:ntot,0:ntot,0:ntot))
|
allocate (array_I_B(-(lmax+ntot):lmax+ntot,0:lmax+ntot,0:ntot,0:ntot,0:ntot))
|
||||||
|
|
||||||
if(ac.eq.0.d0.and.bc.eq.0.d0)then
|
if(ac.eq.0.d0.and.bc.eq.0.d0)then
|
||||||
|
|
||||||
@ -310,36 +311,36 @@ else if(ac.ne.0.d0.and.bc.ne.0.d0)then
|
|||||||
phi_BC0=datan2((b(2)-c(2))/bc,(b(1)-c(1))/bc)
|
phi_BC0=datan2((b(2)-c(2))/bc,(b(1)-c(1))/bc)
|
||||||
|
|
||||||
|
|
||||||
|
do lambdap=0,lmax+ntotB
|
||||||
|
|
||||||
do ktot=0,ntotA+ntotB+nkl_max
|
|
||||||
do lambda=0,lmax+ntotA
|
do lambda=0,lmax+ntotA
|
||||||
do lambdap=0,lmax+ntotB
|
do l=0,lmax
|
||||||
|
do ktot=0,ntotA+ntotB+nkl_max
|
||||||
do k=1,kmax
|
do k=1,kmax
|
||||||
do l=0,lmax
|
array_R(k,ktot,l,lambda,lambdap)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),lambda,lambdap,areal,breal,arg)
|
||||||
array_R(ktot,k,l,lambda,lambdap)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),lambda,lambdap,areal,breal,arg)
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
|
array_coefs_A(k1,1) = binom_func(n_a(1),k1)*(c(1)-a(1))**(n_a(1)-k1)
|
||||||
|
enddo
|
||||||
do k2=0,n_a(2)
|
do k2=0,n_a(2)
|
||||||
|
array_coefs_A(k2,2) = binom_func(n_a(2),k2)*(c(2)-a(2))**(n_a(2)-k2)
|
||||||
|
enddo
|
||||||
do k3=0,n_a(3)
|
do k3=0,n_a(3)
|
||||||
array_coefs_A(k1,k2,k3)=binom_func(n_a(1),k1)*binom_func(n_a(2),k2)*binom_func(n_a(3),k3) &
|
array_coefs_A(k3,3) = binom_func(n_a(3),k3)*(c(3)-a(3))**(n_a(3)-k3)
|
||||||
*(c(1)-a(1))**(n_a(1)-k1)*(c(2)-a(2))**(n_a(2)-k2)*(c(3)-a(3))**(n_a(3)-k3)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
do k1p=0,n_b(1)
|
do k1p=0,n_b(1)
|
||||||
|
array_coefs_B(k1p,1) = binom_func(n_b(1),k1p)*(c(1)-b(1))**(n_b(1)-k1p)
|
||||||
|
enddo
|
||||||
do k2p=0,n_b(2)
|
do k2p=0,n_b(2)
|
||||||
|
array_coefs_B(k2p,2) = binom_func(n_b(2),k2p)*(c(2)-b(2))**(n_b(2)-k2p)
|
||||||
|
enddo
|
||||||
do k3p=0,n_b(3)
|
do k3p=0,n_b(3)
|
||||||
array_coefs_B(k1p,k2p,k3p)=binom_func(n_b(1),k1p)*binom_func(n_b(2),k2p)*binom_func(n_b(3),k3p) &
|
array_coefs_B(k3p,3) = binom_func(n_b(3),k3p)*(c(3)-b(3))**(n_b(3)-k3p)
|
||||||
*(c(1)-b(1))**(n_b(1)-k1p)*(c(2)-b(2))**(n_b(2)-k2p)*(c(3)-b(3))**(n_b(3)-k3p)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!=!=!=!=!=!=!=!
|
!=!=!=!=!=!=!=!
|
||||||
@ -348,70 +349,72 @@ else if(ac.ne.0.d0.and.bc.ne.0.d0)then
|
|||||||
|
|
||||||
accu=0.d0
|
accu=0.d0
|
||||||
do l=0,lmax
|
do l=0,lmax
|
||||||
do m=-l,l
|
do m=-l,l
|
||||||
|
|
||||||
do lambda=0,l+ntotA
|
do k3=0,n_a(3)
|
||||||
do mu=-lambda,lambda
|
do k2=0,n_a(2)
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
do k2=0,n_a(2)
|
do lambda=0,l+ntotA
|
||||||
do k3=0,n_a(3)
|
do mu=-lambda,lambda
|
||||||
array_I_A(lambda,mu,k1,k2,k3)=bigI(lambda,mu,l,m,k1,k2,k3)
|
array_I_A(mu,lambda,k1,k2,k3)=bigI(lambda,mu,l,m,k1,k2,k3)
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
do lambdap=0,l+ntotB
|
do k3p=0,n_b(3)
|
||||||
do mup=-lambdap,lambdap
|
do k2p=0,n_b(2)
|
||||||
do k1p=0,n_b(1)
|
do k1p=0,n_b(1)
|
||||||
do k2p=0,n_b(2)
|
do lambdap=0,l+ntotB
|
||||||
do k3p=0,n_b(3)
|
do mup=-lambdap,lambdap
|
||||||
array_I_B(lambdap,mup,k1p,k2p,k3p)=bigI(lambdap,mup,l,m,k1p,k2p,k3p)
|
array_I_B(mup,lambdap,k1p,k2p,k3p)=bigI(lambdap,mup,l,m,k1p,k2p,k3p)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
do lambda=0,l+ntotA
|
do k3=0,n_a(3)
|
||||||
do mu=-lambda,lambda
|
do k2=0,n_a(2)
|
||||||
|
do k1=0,n_a(1)
|
||||||
|
|
||||||
do k1=0,n_a(1)
|
do lambda=0,l+ntotA
|
||||||
do k2=0,n_a(2)
|
do mu=-lambda,lambda
|
||||||
do k3=0,n_a(3)
|
|
||||||
|
|
||||||
prod=ylm(lambda,mu,theta_AC0,phi_AC0)*array_coefs_A(k1,k2,k3)*array_I_A(lambda,mu,k1,k2,k3)
|
prod=ylm(lambda,mu,theta_AC0,phi_AC0)*array_coefs_A(k1,1)*array_coefs_A(k2,2)*array_coefs_A(k3,3)*array_I_A(mu,lambda,k1,k2,k3)
|
||||||
|
|
||||||
do lambdap=0,l+ntotB
|
|
||||||
do mup=-lambdap,lambdap
|
|
||||||
|
|
||||||
do k1p=0,n_b(1)
|
do k3p=0,n_b(3)
|
||||||
do k2p=0,n_b(2)
|
do k2p=0,n_b(2)
|
||||||
do k3p=0,n_b(3)
|
do k1p=0,n_b(1)
|
||||||
|
do lambdap=0,l+ntotB
|
||||||
|
do mup=-lambdap,lambdap
|
||||||
|
|
||||||
prodp=ylm(lambdap,mup,theta_BC0,phi_BC0)*array_coefs_B(k1p,k2p,k3p)*array_I_B(lambdap,mup,k1p,k2p,k3p)
|
prodp=prod*ylm(lambdap,mup,theta_BC0,phi_BC0)* &
|
||||||
|
array_coefs_B(k1p,1)*array_coefs_B(k2p,2)*array_coefs_B(k3p,3)* &
|
||||||
|
array_I_B(mup,lambdap,k1p,k2p,k3p)
|
||||||
|
|
||||||
do k=1,kmax
|
do k=1,kmax
|
||||||
ktot=k1+k2+k3+k1p+k2p+k3p+n_kl(k,l)
|
ktot=k1+k2+k3+k1p+k2p+k3p+n_kl(k,l)
|
||||||
accu=accu+prod*prodp*v_kl(k,l)*array_R(ktot,k,l,lambda,lambdap)
|
accu=accu+prodp*v_kl(k,l)*array_R(k,ktot,l,lambda,lambdap)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
enddo
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!=!=!=!=!
|
!=!=!=!=!
|
||||||
@ -434,24 +437,24 @@ else if(ac.eq.0.d0.and.bc.ne.0.d0)then
|
|||||||
breal=2.d0*g_b*bc
|
breal=2.d0*g_b*bc
|
||||||
freal=dexp(-g_a*ac**2-g_b*bc**2)
|
freal=dexp(-g_a*ac**2-g_b*bc**2)
|
||||||
|
|
||||||
do ktot=0,ntotA+ntotB+nkl_max
|
do lambdap=0,lmax+ntotB
|
||||||
do lambdap=0,lmax+ntotB
|
do l=0,lmax
|
||||||
do k=1,kmax
|
do ktot=0,ntotA+ntotB+nkl_max
|
||||||
do l=0,lmax
|
do k=1,kmax
|
||||||
array_R(ktot,k,l,0,lambdap)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),0,lambdap,areal,breal,arg)
|
array_R(k,ktot,l,0,lambdap)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),0,lambdap,areal,breal,arg)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
do k1p=0,n_b(1)
|
do k1p=0,n_b(1)
|
||||||
|
array_coefs_B(k1p,1) = binom_func(n_b(1),k1p)*(c(1)-b(1))**(n_b(1)-k1p)
|
||||||
|
enddo
|
||||||
do k2p=0,n_b(2)
|
do k2p=0,n_b(2)
|
||||||
|
array_coefs_B(k2p,2) = binom_func(n_b(2),k2p)*(c(2)-b(2))**(n_b(2)-k2p)
|
||||||
|
enddo
|
||||||
do k3p=0,n_b(3)
|
do k3p=0,n_b(3)
|
||||||
|
array_coefs_B(k3p,3) = binom_func(n_b(3),k3p)*(c(3)-b(3))**(n_b(3)-k3p)
|
||||||
array_coefs_B(k1p,k2p,k3p)=binom_func(n_b(1),k1p)*binom_func(n_b(2),k2p)*binom_func(n_b(3),k3p) &
|
|
||||||
*(c(1)-b(1))**(n_b(1)-k1p)*(c(2)-b(2))**(n_b(2)-k2p)*(c(3)-b(3))**(n_b(3)-k3p)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!=!=!=!=!=!=!=!
|
!=!=!=!=!=!=!=!
|
||||||
@ -460,43 +463,43 @@ else if(ac.eq.0.d0.and.bc.ne.0.d0)then
|
|||||||
|
|
||||||
accu=0.d0
|
accu=0.d0
|
||||||
do l=0,lmax
|
do l=0,lmax
|
||||||
do m=-l,l
|
do m=-l,l
|
||||||
|
|
||||||
do lambdap=0,l+ntotB
|
do k3p=0,n_b(3)
|
||||||
do mup=-lambdap,lambdap
|
do k2p=0,n_b(2)
|
||||||
do k1p=0,n_b(1)
|
|
||||||
do k2p=0,n_b(2)
|
|
||||||
do k3p=0,n_b(3)
|
|
||||||
array_I_B(lambdap,mup,k1p,k2p,k3p)=bigI(lambdap,mup,l,m,k1p,k2p,k3p)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
|
|
||||||
prod=bigI(0,0,l,m,n_a(1),n_a(2),n_a(3))
|
|
||||||
|
|
||||||
do lambdap=0,l+ntotB
|
|
||||||
do mup=-lambdap,lambdap
|
|
||||||
do k1p=0,n_b(1)
|
do k1p=0,n_b(1)
|
||||||
do k2p=0,n_b(2)
|
do lambdap=0,l+ntotB
|
||||||
do k3p=0,n_b(3)
|
do mup=-lambdap,lambdap
|
||||||
|
array_I_B(mup,lambdap,k1p,k2p,k3p)=bigI(lambdap,mup,l,m,k1p,k2p,k3p)
|
||||||
prodp=array_coefs_B(k1p,k2p,k3p)*ylm(lambdap,mup,theta_BC0,phi_BC0)*array_I_B(lambdap,mup,k1p,k2p,k3p)
|
enddo
|
||||||
|
enddo
|
||||||
do k=1,kmax
|
|
||||||
|
|
||||||
ktot=ntotA+k1p+k2p+k3p+n_kl(k,l)
|
|
||||||
accu=accu+prod*prodp*v_kl(k,l)*array_R(ktot,k,l,0,lambdap)
|
|
||||||
|
|
||||||
enddo
|
|
||||||
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
prod=bigI(0,0,l,m,n_a(1),n_a(2),n_a(3))
|
||||||
|
|
||||||
|
do k3p=0,n_b(3)
|
||||||
|
do k2p=0,n_b(2)
|
||||||
|
do k1p=0,n_b(1)
|
||||||
|
do lambdap=0,l+ntotB
|
||||||
|
do mup=-lambdap,lambdap
|
||||||
|
|
||||||
|
prodp=prod*array_coefs_B(k1p,1)*array_coefs_B(k2p,2)*array_coefs_B(k3p,3)*ylm(lambdap,mup,theta_BC0,phi_BC0)*array_I_B(mup,lambdap,k1p,k2p,k3p)
|
||||||
|
|
||||||
|
do k=1,kmax
|
||||||
|
|
||||||
|
ktot=ntotA+k1p+k2p+k3p+n_kl(k,l)
|
||||||
|
accu=accu+prodp*v_kl(k,l)*array_R(k,ktot,l,0,lambdap)
|
||||||
|
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!=!=!=!=!
|
!=!=!=!=!
|
||||||
@ -519,26 +522,24 @@ else if(ac.ne.0.d0.and.bc.eq.0.d0)then
|
|||||||
breal=2.d0*g_b*bc
|
breal=2.d0*g_b*bc
|
||||||
freal=dexp(-g_a*ac**2-g_b*bc**2)
|
freal=dexp(-g_a*ac**2-g_b*bc**2)
|
||||||
|
|
||||||
do ktot=0,ntotA+ntotB+nkl_max
|
do lambda=0,lmax+ntotA
|
||||||
do lambda=0,lmax+ntotA
|
do l=0,lmax
|
||||||
do k=1,kmax
|
do ktot=0,ntotA+ntotB+nkl_max
|
||||||
do l=0,lmax
|
do k=1,kmax
|
||||||
|
array_R(k,ktot,l,lambda,0)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),lambda,0,areal,breal,arg)
|
||||||
array_R(ktot,k,l,lambda,0)= int_prod_bessel(ktot+2,g_a+g_b+dz_kl(k,l),lambda,0,areal,breal,arg)
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
|
array_coefs_A(k1,1) = binom_func(n_a(1),k1)*(c(1)-a(1))**(n_a(1)-k1)
|
||||||
|
enddo
|
||||||
do k2=0,n_a(2)
|
do k2=0,n_a(2)
|
||||||
|
array_coefs_A(k2,2) = binom_func(n_a(2),k2)*(c(2)-a(2))**(n_a(2)-k2)
|
||||||
|
enddo
|
||||||
do k3=0,n_a(3)
|
do k3=0,n_a(3)
|
||||||
|
array_coefs_A(k3,3) = binom_func(n_a(3),k3)*(c(3)-a(3))**(n_a(3)-k3)
|
||||||
array_coefs_A(k1,k2,k3)=binom_func(n_a(1),k1)*binom_func(n_a(2),k2)*binom_func(n_a(3),k3) &
|
|
||||||
*(c(1)-a(1))**(n_a(1)-k1)*(c(2)-a(2))**(n_a(2)-k2)*(c(3)-a(3))**(n_a(3)-k3)
|
|
||||||
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
!=!=!=!=!=!=!=!
|
!=!=!=!=!=!=!=!
|
||||||
@ -549,36 +550,36 @@ else if(ac.ne.0.d0.and.bc.eq.0.d0)then
|
|||||||
do l=0,lmax
|
do l=0,lmax
|
||||||
do m=-l,l
|
do m=-l,l
|
||||||
|
|
||||||
do lambda=0,l+ntotA
|
do k3=0,n_a(3)
|
||||||
do mu=-lambda,lambda
|
do k2=0,n_a(2)
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
do k2=0,n_a(2)
|
do lambda=0,l+ntotA
|
||||||
do k3=0,n_a(3)
|
do mu=-lambda,lambda
|
||||||
array_I_A(lambda,mu,k1,k2,k3)=bigI(lambda,mu,l,m,k1,k2,k3)
|
array_I_A(mu,lambda,k1,k2,k3)=bigI(lambda,mu,l,m,k1,k2,k3)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
do lambda=0,l+ntotA
|
do k3=0,n_a(3)
|
||||||
do mu=-lambda,lambda
|
do k2=0,n_a(2)
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
do k2=0,n_a(2)
|
do lambda=0,l+ntotA
|
||||||
do k3=0,n_a(3)
|
do mu=-lambda,lambda
|
||||||
|
|
||||||
prod=array_coefs_A(k1,k2,k3)*ylm(lambda,mu,theta_AC0,phi_AC0)*array_I_A(lambda,mu,k1,k2,k3)
|
prod=array_coefs_A(k1,1)*array_coefs_A(k2,2)*array_coefs_A(k3,3)*ylm(lambda,mu,theta_AC0,phi_AC0)*array_I_A(mu,lambda,k1,k2,k3)
|
||||||
prodp=bigI(0,0,l,m,n_b(1),n_b(2),n_b(3))
|
prodp=prod*bigI(0,0,l,m,n_b(1),n_b(2),n_b(3))
|
||||||
|
|
||||||
do k=1,kmax
|
do k=1,kmax
|
||||||
ktot=k1+k2+k3+ntotB+n_kl(k,l)
|
ktot=k1+k2+k3+ntotB+n_kl(k,l)
|
||||||
accu=accu+prod*prodp*v_kl(k,l)*array_R(ktot,k,l,lambda,0)
|
accu=accu+prodp*v_kl(k,l)*array_R(k,ktot,l,lambda,0)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
enddo
|
enddo
|
||||||
@ -850,22 +851,27 @@ implicit none
|
|||||||
integer lambda,mu,l,m,k1,k2,k3
|
integer lambda,mu,l,m,k1,k2,k3
|
||||||
integer k,i,kp,ip
|
integer k,i,kp,ip
|
||||||
double precision pi,sum,factor1,factor2,cylm,cylmp,bigA,binom_func,fact,coef_pm
|
double precision pi,sum,factor1,factor2,cylm,cylmp,bigA,binom_func,fact,coef_pm
|
||||||
|
double precision sgn, sgnp
|
||||||
pi=dacos(-1.d0)
|
pi=dacos(-1.d0)
|
||||||
|
|
||||||
if(mu.gt.0.and.m.gt.0)then
|
if(mu.gt.0.and.m.gt.0)then
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
|
sgn = 1.d0
|
||||||
do k=0,mu/2
|
do k=0,mu/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
|
sgnp = 1.d0
|
||||||
do kp=0,m/2
|
do kp=0,m/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm=sgn*factor1*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-2*k+m-2*kp+k1,2*k+2*kp+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(mu-2*k+m-2*kp+k1,2*k+2*kp+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
bigI=sum
|
bigI=sum
|
||||||
return
|
return
|
||||||
@ -888,15 +894,17 @@ endif
|
|||||||
|
|
||||||
if(mu.eq.0.and.m.gt.0)then
|
if(mu.eq.0.and.m.gt.0)then
|
||||||
factor1=dsqrt((2*lambda+1)/(4.d0*pi))
|
factor1=dsqrt((2*lambda+1)/(4.d0*pi))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
do i=0,lambda
|
do i=0,lambda
|
||||||
|
sgnp = 1.d0
|
||||||
do kp=0,m/2
|
do kp=0,m/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=factor1*coef_pm(lambda,i)
|
cylm=factor1*coef_pm(lambda,i)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(m-2*kp+k1,2*kp+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(m-2*kp+k1,2*kp+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
bigI=sum
|
bigI=sum
|
||||||
@ -905,16 +913,18 @@ endif
|
|||||||
|
|
||||||
if(mu.gt.0.and.m.eq.0)then
|
if(mu.gt.0.and.m.eq.0)then
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)/(4.d0*pi))
|
factor2=dsqrt((2*l+1)/(4.d0*pi))
|
||||||
|
sgn = 1.d0
|
||||||
do k=0,mu/2
|
do k=0,mu/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
do ip=0,l
|
do ip=0,l
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm=sgn*factor1*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=factor2*coef_pm(l,ip)
|
cylmp=factor2*coef_pm(l,ip)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-2*k +k1,2*k +k2,i+ip +k3)
|
sum=sum+cylm*cylmp*bigA(mu-2*k +k1,2*k +k2,i+ip +k3)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
bigI=sum
|
bigI=sum
|
||||||
return
|
return
|
||||||
@ -923,19 +933,23 @@ endif
|
|||||||
if(mu.lt.0.and.m.lt.0)then
|
if(mu.lt.0.and.m.lt.0)then
|
||||||
mu=-mu
|
mu=-mu
|
||||||
m=-m
|
m=-m
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
|
sgn = 1.d0
|
||||||
do k=0,(mu-1)/2
|
do k=0,(mu-1)/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
|
sgnp = 1.d0
|
||||||
do kp=0,(m-1)/2
|
do kp=0,(m-1)/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm=sgn*factor1*binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+m-(2*kp+1)+k1,(2*k+1)+(2*kp+1)+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+m-(2*kp+1)+k1,(2*k+1)+(2*kp+1)+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
mu=-mu
|
mu=-mu
|
||||||
m=-m
|
m=-m
|
||||||
@ -946,15 +960,17 @@ endif
|
|||||||
if(mu.eq.0.and.m.lt.0)then
|
if(mu.eq.0.and.m.lt.0)then
|
||||||
m=-m
|
m=-m
|
||||||
factor1=dsqrt((2*lambda+1)/(4.d0*pi))
|
factor1=dsqrt((2*lambda+1)/(4.d0*pi))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
do i=0,lambda
|
do i=0,lambda
|
||||||
|
sgnp = 1.d0
|
||||||
do kp=0,(m-1)/2
|
do kp=0,(m-1)/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=factor1*coef_pm(lambda,i)
|
cylm=factor1*coef_pm(lambda,i)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(m-(2*kp+1)+k1,2*kp+1+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(m-(2*kp+1)+k1,2*kp+1+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
m=-m
|
m=-m
|
||||||
@ -965,16 +981,18 @@ endif
|
|||||||
if(mu.lt.0.and.m.eq.0)then
|
if(mu.lt.0.and.m.eq.0)then
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
mu=-mu
|
mu=-mu
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)/(4.d0*pi))
|
factor2=dsqrt((2*l+1)/(4.d0*pi))
|
||||||
|
sgn = 1.d0
|
||||||
do k=0,(mu-1)/2
|
do k=0,(mu-1)/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
do ip=0,l
|
do ip=0,l
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm=sgn*factor1*binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=factor2*coef_pm(l,ip)
|
cylmp=factor2*coef_pm(l,ip)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+k1,2*k+1+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+k1,2*k+1+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
mu=-mu
|
mu=-mu
|
||||||
bigI=sum
|
bigI=sum
|
||||||
@ -983,19 +1001,23 @@ endif
|
|||||||
|
|
||||||
if(mu.gt.0.and.m.lt.0)then
|
if(mu.gt.0.and.m.lt.0)then
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
m=-m
|
m=-m
|
||||||
|
sgn=1.d0
|
||||||
do k=0,mu/2
|
do k=0,mu/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
|
sgnp=1.d0
|
||||||
do kp=0,(m-1)/2
|
do kp=0,(m-1)/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm =sgn *factor1*binom_func(mu,2*k)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp+1)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-2*k+m-(2*kp+1)+k1,2*k+2*kp+1+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(mu-2*k+m-(2*kp+1)+k1,2*k+2*kp+1+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
m=-m
|
m=-m
|
||||||
bigI=sum
|
bigI=sum
|
||||||
@ -1004,19 +1026,23 @@ endif
|
|||||||
|
|
||||||
if(mu.lt.0.and.m.gt.0)then
|
if(mu.lt.0.and.m.gt.0)then
|
||||||
mu=-mu
|
mu=-mu
|
||||||
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(4.d0*pi*fact(lambda+iabs(mu))))
|
factor1=dsqrt((2*lambda+1)*fact(lambda-iabs(mu))/(2.d0*pi*fact(lambda+iabs(mu))))
|
||||||
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(4.d0*pi*fact(l+iabs(m))))
|
factor2=dsqrt((2*l+1)*fact(l-iabs(m))/(2.d0*pi*fact(l+iabs(m))))
|
||||||
sum=0.d0
|
sum=0.d0
|
||||||
|
sgn = 1.d0
|
||||||
do k=0,(mu-1)/2
|
do k=0,(mu-1)/2
|
||||||
do i=0,lambda-mu
|
do i=0,lambda-mu
|
||||||
|
sgnp = 1.d0
|
||||||
do kp=0,m/2
|
do kp=0,m/2
|
||||||
do ip=0,l-m
|
do ip=0,l-m
|
||||||
cylm=(-1.d0)**k*factor1*dsqrt(2.d0)*binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
cylm=sgn*factor1 *binom_func(mu,2*k+1)*fact(mu+i)/fact(i)*coef_pm(lambda,i+mu)
|
||||||
cylmp=(-1.d0)**kp*factor2*dsqrt(2.d0)*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
cylmp=sgnp*factor2*binom_func(m,2*kp)*fact(m+ip)/fact(ip)*coef_pm(l,ip+m)
|
||||||
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+m-2*kp+k1,2*k+1+2*kp+k2,i+ip+k3)
|
sum=sum+cylm*cylmp*bigA(mu-(2*k+1)+m-2*kp+k1,2*k+1+2*kp+k2,i+ip+k3)
|
||||||
enddo
|
enddo
|
||||||
|
sgnp = -sgnp
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
sgn = -sgn
|
||||||
enddo
|
enddo
|
||||||
bigI=sum
|
bigI=sum
|
||||||
mu=-mu
|
mu=-mu
|
||||||
@ -1128,28 +1154,49 @@ end
|
|||||||
!
|
!
|
||||||
IMPLICIT DOUBLE PRECISION (P,X)
|
IMPLICIT DOUBLE PRECISION (P,X)
|
||||||
DIMENSION PM(0:MM,0:(N+1))
|
DIMENSION PM(0:MM,0:(N+1))
|
||||||
DO 10 I=0,N
|
DOUBLE PRECISION, SAVE :: INVERSE(100) = 0.D0
|
||||||
DO 10 J=0,M
|
DOUBLE PRECISION :: LS, II, JJ
|
||||||
10 PM(J,I)=0.0D0
|
IF (INVERSE(1) == 0.d0) THEN
|
||||||
|
DO I=1,100
|
||||||
|
INVERSE(I) = 1.D0/DBLE(I)
|
||||||
|
ENDDO
|
||||||
|
ENDIF
|
||||||
|
DO I=0,N
|
||||||
|
DO J=0,M
|
||||||
|
PM(J,I)=0.0D0
|
||||||
|
ENDDO
|
||||||
|
ENDDO
|
||||||
PM(0,0)=1.0D0
|
PM(0,0)=1.0D0
|
||||||
IF (DABS(X).EQ.1.0D0) THEN
|
IF (DABS(X).EQ.1.0D0) THEN
|
||||||
DO 15 I=1,N
|
DO I=1,N
|
||||||
15 PM(0,I)=X**I
|
PM(0,I)=X**I
|
||||||
|
ENDDO
|
||||||
RETURN
|
RETURN
|
||||||
ENDIF
|
ENDIF
|
||||||
LS=1
|
LS=1.D0
|
||||||
IF (DABS(X).GT.1.0D0) LS=-1
|
IF (DABS(X).GT.1.0D0) LS=-1.D0
|
||||||
XQ=DSQRT(LS*(1.0D0-X*X))
|
XQ=DSQRT(LS*(1.0D0-X*X))
|
||||||
XS=LS*(1.0D0-X*X)
|
XS=LS*(1.0D0-X*X)
|
||||||
DO 30 I=1,M
|
II = 1.D0
|
||||||
30 PM(I,I)=-LS*(2.0D0*I-1.0D0)*XQ*PM(I-1,I-1)
|
DO I=1,M
|
||||||
DO 35 I=0,M
|
PM(I,I)=-LS*II*XQ*PM(I-1,I-1)
|
||||||
35 PM(I,I+1)=(2.0D0*I+1.0D0)*X*PM(I,I)
|
II = II+2.D0
|
||||||
|
ENDDO
|
||||||
|
II = 1.D0
|
||||||
|
DO I=0,M
|
||||||
|
PM(I,I+1)=II*X*PM(I,I)
|
||||||
|
II = II+2.D0
|
||||||
|
ENDDO
|
||||||
|
|
||||||
DO 40 I=0,M
|
II = 0.D0
|
||||||
DO 40 J=I+2,N
|
DO I=0,M
|
||||||
PM(I,J)=((2.0D0*J-1.0D0)*X*PM(I,J-1)- (I+J-1.0D0)*PM(I,J-2))/(J-I)
|
JJ = II+2.D0
|
||||||
40 CONTINUE
|
DO J=I+2,N
|
||||||
|
PM(I,J)=((2.0D0*JJ-1.0D0)*X*PM(I,J-1)- (II+JJ-1.0D0)*PM(I,J-2))*INVERSE(J-I)
|
||||||
|
JJ = JJ+1.D0
|
||||||
|
ENDDO
|
||||||
|
II = II+1.D0
|
||||||
|
ENDDO
|
||||||
END
|
END
|
||||||
|
|
||||||
|
|
||||||
@ -1703,17 +1750,37 @@ end
|
|||||||
!c
|
!c
|
||||||
double precision function ylm(l,m,theta,phi)
|
double precision function ylm(l,m,theta,phi)
|
||||||
implicit none
|
implicit none
|
||||||
integer l,m
|
integer l,m,i
|
||||||
double precision theta,phi,pm,factor,pi,x,fact,sign
|
double precision theta,phi,pm,factor,twopi,x,fact,sign
|
||||||
DIMENSION PM(0:100,0:100)
|
DIMENSION PM(0:100,0:100)
|
||||||
pi=dacos(-1.d0)
|
twopi=2.d0*dacos(-1.d0)
|
||||||
x=dcos(theta)
|
x=dcos(theta)
|
||||||
sign=(-1.d0)**m
|
if (iand(m,1) == 1) then
|
||||||
|
sign=-1.d0
|
||||||
|
else
|
||||||
|
sign=1.d0
|
||||||
|
endif
|
||||||
CALL LPMN(100,l,l,X,PM)
|
CALL LPMN(100,l,l,X,PM)
|
||||||
factor=dsqrt( (2*l+1)*fact(l-iabs(m)) /(4.d0*pi*fact(l+iabs(m))) )
|
if (m > 0) then
|
||||||
if(m.gt.0)ylm=sign*dsqrt(2.d0)*factor*pm(m,l)*dcos(dfloat(m)*phi)
|
factor=dsqrt((l+l+1)*fact(l-m) /(twopi*fact(l+m)) )
|
||||||
if(m.eq.0)ylm=factor*pm(m,l)
|
! factor = dble(l+m)
|
||||||
if(m.lt.0)ylm=sign*dsqrt(2.d0)*factor*pm(iabs(m),l)*dsin(dfloat(iabs(m))*phi)
|
! do i=-m,m-1
|
||||||
|
! factor = factor * (factor - 1.d0)
|
||||||
|
! enddo
|
||||||
|
! factor=dsqrt(dble(l+l+1)/(twopi*factor) )
|
||||||
|
ylm=sign*factor*pm(m,l)*dcos(dfloat(m)*phi)
|
||||||
|
else if (m == 0) then
|
||||||
|
factor=dsqrt( 0.5d0*(l+l+1) /twopi )
|
||||||
|
ylm=factor*pm(m,l)
|
||||||
|
else if (m < 0) then
|
||||||
|
factor=dsqrt( (l+l+1)*fact(l+m) /(twopi*fact(l-m)) )
|
||||||
|
! factor = dble(l-m)
|
||||||
|
! do i=m,-m-1
|
||||||
|
! factor = factor * (factor - 1.d0)
|
||||||
|
! enddo
|
||||||
|
! factor=dsqrt(dble(l+l+1)/(twopi*factor) )
|
||||||
|
ylm=sign*factor*pm(-m,l)*dsin(dfloat(-m)*phi)
|
||||||
|
endif
|
||||||
end
|
end
|
||||||
|
|
||||||
!c Explicit representation of Legendre polynomials P_n(x)
|
!c Explicit representation of Legendre polynomials P_n(x)
|
||||||
@ -1829,11 +1896,12 @@ end
|
|||||||
double precision function binom_gen(alpha,n)
|
double precision function binom_gen(alpha,n)
|
||||||
implicit none
|
implicit none
|
||||||
integer :: n,k
|
integer :: n,k
|
||||||
double precision :: fact,alpha,prod
|
double precision :: fact,alpha,prod, factn_inv
|
||||||
prod=1.d0
|
prod=1.d0
|
||||||
|
factn_inv = 1.d0/fact(n)
|
||||||
do k=0,n-1
|
do k=0,n-1
|
||||||
prod=prod*(alpha-k)
|
prod=prod*(alpha-k)
|
||||||
binom_gen = prod/(fact(n))
|
binom_gen = prod*factn_inv
|
||||||
enddo
|
enddo
|
||||||
end
|
end
|
||||||
|
|
||||||
@ -1881,6 +1949,7 @@ double precision function int_prod_bessel(l,gam,n,m,a,b,arg)
|
|||||||
double precision :: term_A, term_B, term_rap, expo
|
double precision :: term_A, term_B, term_rap, expo
|
||||||
double precision :: s_q_0, s_q_k, s_0_0, a_over_b_square
|
double precision :: s_q_0, s_q_k, s_0_0, a_over_b_square
|
||||||
double precision :: int_prod_bessel_loc
|
double precision :: int_prod_bessel_loc
|
||||||
|
double precision :: inverses(0:300)
|
||||||
|
|
||||||
logical done
|
logical done
|
||||||
|
|
||||||
@ -1927,18 +1996,32 @@ double precision function int_prod_bessel(l,gam,n,m,a,b,arg)
|
|||||||
! Initialise the first recurence terme for the q loop
|
! Initialise the first recurence terme for the q loop
|
||||||
s_q_0 = s_0_0
|
s_q_0 = s_0_0
|
||||||
|
|
||||||
|
|
||||||
! Loop over q for the convergence of the sequence
|
! Loop over q for the convergence of the sequence
|
||||||
do while (.not.done)
|
do while (.not.done)
|
||||||
|
|
||||||
! Init
|
! Init
|
||||||
sum=0
|
|
||||||
s_q_k=s_q_0
|
s_q_k=s_q_0
|
||||||
|
sum=s_q_0
|
||||||
|
|
||||||
! Iteration of k
|
if (q>300) then
|
||||||
do k=0,q
|
stop 'pseudopot.f90 : q > 300'
|
||||||
|
endif
|
||||||
|
|
||||||
|
do k=0,q-1
|
||||||
|
s_q_k = ( dble(2*(q-k+m)+1)*dble(q-k)*inverses(k) ) * s_q_k
|
||||||
sum=sum+s_q_k
|
sum=sum+s_q_k
|
||||||
s_q_k = a_over_b_square * ( dble(2*(q-k+m)+1)/dble(2*(k+n)+3) ) * ( dble(q-k)/dble(k+1)) * s_q_k
|
|
||||||
enddo
|
enddo
|
||||||
|
inverses(q) = a_over_b_square/(dble(2*(q+n)+3) * dble(q+1))
|
||||||
|
! do k=0,q
|
||||||
|
! sum=sum+s_q_k
|
||||||
|
! s_q_k = a_over_b_square * ( dble(2*(q-k+m)+1)*dble(q-k)/(dble(2*(k+n)+3) * dble(k+1)) ) * s_q_k
|
||||||
|
! enddo
|
||||||
|
! Iteration of k
|
||||||
|
! do k=0,q
|
||||||
|
! sum=sum+s_q_k
|
||||||
|
! s_q_k = a_over_b_square * ( dble(2*(q-k+m)+1)*dble(q-k)/(dble(2*(k+n)+3) * dble(k+1)) ) * s_q_k
|
||||||
|
! enddo
|
||||||
|
|
||||||
int=int+sum
|
int=int+sum
|
||||||
|
|
||||||
@ -2120,15 +2203,15 @@ parameter (ntot_max=14)
|
|||||||
integer l,m
|
integer l,m
|
||||||
double precision a(3),g_a,c(3)
|
double precision a(3),g_a,c(3)
|
||||||
double precision prod,binom_func,accu,bigI,ylm,bessel_mod
|
double precision prod,binom_func,accu,bigI,ylm,bessel_mod
|
||||||
double precision theta_AC0,phi_AC0,ac,factor,fourpi,arg,r,areal
|
double precision theta_AC0,phi_AC0,ac,ac2,factor,fourpi,arg,r,areal
|
||||||
integer ntotA,mu,k1,k2,k3,lambda
|
integer ntotA,mu,k1,k2,k3,lambda
|
||||||
integer n_a(3)
|
integer n_a(3)
|
||||||
double precision &
|
double precision y, f1, f2
|
||||||
array_I_A(0:lmax_max+ntot_max,-(lmax_max+ntot_max):lmax_max+ntot_max,0:ntot_max,0:ntot_max,0:ntot_max)
|
double precision, allocatable :: array_coefs_A(:,:)
|
||||||
double precision array_coefs_A(0:ntot_max,0:ntot_max,0:ntot_max), y
|
|
||||||
|
|
||||||
ac=dsqrt((a(1)-c(1))**2+(a(2)-c(2))**2+(a(3)-c(3))**2)
|
ac2=(a(1)-c(1))**2+(a(2)-c(2))**2+(a(3)-c(3))**2
|
||||||
arg=g_a*(ac**2+r**2)
|
ac=dsqrt(ac2)
|
||||||
|
arg=g_a*(ac2+r*r)
|
||||||
fourpi=4.d0*dacos(-1.d0)
|
fourpi=4.d0*dacos(-1.d0)
|
||||||
factor=fourpi*dexp(-arg)
|
factor=fourpi*dexp(-arg)
|
||||||
areal=2.d0*g_a*ac
|
areal=2.d0*g_a*ac
|
||||||
@ -2144,51 +2227,45 @@ else
|
|||||||
theta_AC0=dacos( (a(3)-c(3))/ac )
|
theta_AC0=dacos( (a(3)-c(3))/ac )
|
||||||
phi_AC0=datan2((a(2)-c(2))/ac,(a(1)-c(1))/ac)
|
phi_AC0=datan2((a(2)-c(2))/ac,(a(1)-c(1))/ac)
|
||||||
|
|
||||||
|
allocate (array_coefs_A(0:ntotA,3))
|
||||||
do k1=0,n_a(1)
|
do k1=0,n_a(1)
|
||||||
do k2=0,n_a(2)
|
array_coefs_A(k1,1) = binom_func(n_a(1),k1)*(c(1)-a(1))**(n_a(1)-k1)*r**(k1)
|
||||||
do k3=0,n_a(3)
|
|
||||||
array_coefs_A(k1,k2,k3)=binom_func(n_a(1),k1)*binom_func(n_a(2),k2)*binom_func(n_a(3),k3) &
|
|
||||||
*(c(1)-a(1))**(n_a(1)-k1)*(c(2)-a(2))**(n_a(2)-k2)*(c(3)-a(3))**(n_a(3)-k3) &
|
|
||||||
*r**(k1+k2+k3)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
do k2=0,n_a(2)
|
||||||
do lambda=0,l+ntotA
|
array_coefs_A(k2,2) = binom_func(n_a(2),k2)*(c(2)-a(2))**(n_a(2)-k2)*r**(k2)
|
||||||
do mu=-lambda,lambda
|
enddo
|
||||||
do k1=0,n_a(1)
|
do k3=0,n_a(3)
|
||||||
do k2=0,n_a(2)
|
array_coefs_A(k3,3) = binom_func(n_a(3),k3)*(c(3)-a(3))**(n_a(3)-k3)*r**(k3)
|
||||||
do k3=0,n_a(3)
|
|
||||||
array_I_A(lambda,mu,k1,k2,k3)=bigI(lambda,mu,l,m,k1,k2,k3)
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
accu=0.d0
|
accu=0.d0
|
||||||
do lambda=0,l+ntotA
|
do lambda=0,l+ntotA
|
||||||
do mu=-lambda,lambda
|
do mu=-lambda,lambda
|
||||||
y = ylm(lambda,mu,theta_AC0,phi_AC0)
|
y = ylm(lambda,mu,theta_AC0,phi_AC0)
|
||||||
if (y == 0.d0) then
|
if (y == 0.d0) then
|
||||||
cycle
|
cycle
|
||||||
endif
|
endif
|
||||||
do k1=0,n_a(1)
|
do k3=0,n_a(3)
|
||||||
do k2=0,n_a(2)
|
f1 = y*array_coefs_A(k3,3)
|
||||||
do k3=0,n_a(3)
|
if (f1 == 0.d0) cycle
|
||||||
prod=y*array_coefs_A(k1,k2,k3)*array_I_A(lambda,mu,k1,k2,k3)
|
do k2=0,n_a(2)
|
||||||
if (prod == 0.d0) then
|
f2 = f1*array_coefs_A(k2,2)
|
||||||
cycle
|
if (f2 == 0.d0) cycle
|
||||||
endif
|
do k1=0,n_a(1)
|
||||||
if (areal*r < 100.d0) then ! overflow!
|
prod=f2*array_coefs_A(k1,1)*bigI(lambda,mu,l,m,k1,k2,k3)
|
||||||
accu=accu+prod*bessel_mod(areal*r,lambda)
|
if (prod == 0.d0) then
|
||||||
endif
|
cycle
|
||||||
|
endif
|
||||||
|
if (areal*r < 100.d0) then ! overflow!
|
||||||
|
accu=accu+prod*bessel_mod(areal*r,lambda)
|
||||||
|
endif
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
|
enddo
|
||||||
enddo
|
enddo
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
|
||||||
enddo
|
enddo
|
||||||
ylm_orb=factor*accu
|
ylm_orb=factor*accu
|
||||||
|
deallocate (array_coefs_A)
|
||||||
return
|
return
|
||||||
endif
|
endif
|
||||||
|
|
||||||
|
Before Width: | Height: | Size: 37 KiB After Width: | Height: | Size: 36 KiB |
Before Width: | Height: | Size: 41 KiB After Width: | Height: | Size: 41 KiB |
Before Width: | Height: | Size: 24 KiB After Width: | Height: | Size: 24 KiB |
Before Width: | Height: | Size: 10 KiB After Width: | Height: | Size: 10 KiB |
@ -51,7 +51,7 @@ size: (nuclei.nucl_num,pseudo.pseudo_kmax,0:pseudo.pseudo_lmax)
|
|||||||
|
|
||||||
[do_pseudo]
|
[do_pseudo]
|
||||||
type: logical
|
type: logical
|
||||||
doc: Using pseudo potential integral of not
|
doc: Using pseudo potential integral or not
|
||||||
interface: ezfio,provider,ocaml
|
interface: ezfio,provider,ocaml
|
||||||
default: False
|
default: False
|
||||||
|
|
||||||
|
@ -29,7 +29,7 @@ Documentation
|
|||||||
|
|
||||||
|
|
||||||
`do_pseudo <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo/ezfio_interface.irp.f#L248>`_
|
`do_pseudo <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo/ezfio_interface.irp.f#L248>`_
|
||||||
Using pseudo potential integral of not
|
Using pseudo potential integral or not
|
||||||
|
|
||||||
|
|
||||||
`pseudo_dz_k <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo/ezfio_interface.irp.f#L204>`_
|
`pseudo_dz_k <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo/ezfio_interface.irp.f#L204>`_
|
||||||
|
Before Width: | Height: | Size: 14 KiB After Width: | Height: | Size: 14 KiB |
Before Width: | Height: | Size: 2.4 KiB After Width: | Height: | Size: 2.3 KiB |