From e46ea31b961bfe8db8388f6b8bda8afe34de7862 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 4 Jan 2019 00:03:37 +0100 Subject: [PATCH] Tests OK for CISD --- docs/source/auto_generate.py | 7 +- .../programmers_guide/index_providers.rst | 183 +++++- src/ao_one_e_integrals/pot_ao_erf_ints.irp.f | 617 +++++++++--------- src/ao_one_e_integrals/pot_ao_ints.irp.f | 24 +- .../pot_ao_pseudo_ints.irp.f | 7 - src/cisd/40.cisd.bats | 2 +- src/mo_basis/mos.irp.f | 4 +- 7 files changed, 510 insertions(+), 334 deletions(-) diff --git a/docs/source/auto_generate.py b/docs/source/auto_generate.py index 433cf544..6cf98569 100755 --- a/docs/source/auto_generate.py +++ b/docs/source/auto_generate.py @@ -123,7 +123,12 @@ def generate_providers(abs_module): state = 0 entity = { "decl": [], "doc": [] , "name": name , "module": module } - for line in f.readlines(): + text=f.read() + text_old = None + while text_old != text: + text_old = text + text = text.replace("$"," :math:`",1).replace("$","` ") + for line in text.splitlines(): line = line.rstrip() if line.startswith(".SH Declaration"): state = 1 diff --git a/docs/source/programmers_guide/index_providers.rst b/docs/source/programmers_guide/index_providers.rst index 1888ed4f..eca78fcb 100644 --- a/docs/source/programmers_guide/index_providers.rst +++ b/docs/source/programmers_guide/index_providers.rst @@ -112,6 +112,7 @@ Index of Providers * :c:data:`aos_vx_alpha_pbe_w` * :c:data:`aos_vx_beta_lda_w` * :c:data:`aos_vx_beta_pbe_w` +* :c:data:`apply_exc_to_psi` * :c:data:`barycentric_electronic_energy` * :c:data:`bi_elec_ref_bitmask_energy` * :c:data:`big_array_coulomb_integrals` @@ -133,8 +134,12 @@ Index of Providers * :c:data:`cart_to_sphe_9` * :c:data:`cas_bitmask` * :c:data:`center_of_mass` +* :c:data:`ci_dressed_pt2_new_eigenvectors` +* :c:data:`ci_dressed_pt2_new_eigenvectors_s2` +* :c:data:`ci_dressed_pt2_new_energy` * :c:data:`ci_eigenvectors` * :c:data:`ci_eigenvectors_s2` +* :c:data:`ci_electronic_dressed_pt2_new_energy` * :c:data:`ci_electronic_energy` * :c:data:`ci_energy` * :c:data:`closed_shell_ref_bitmask` @@ -146,6 +151,7 @@ Index of Providers * :c:data:`core_fock_operator_erf` * :c:data:`core_inact_act_bitmask_4` * :c:data:`core_inact_virt_bitmask` +* :c:data:`corr_e_from_1h1p` * :c:data:`correlation_energy_ratio_max` * :c:data:`correlation_functional` * :c:data:`damping_for_rs_dft` @@ -158,6 +164,7 @@ Index of Providers * :c:data:`degree_max_generators` * :c:data:`degree_max_integration_lebedev` * :c:data:`delta_e_per_selector` +* :c:data:`delta_ij_mrpt` * :c:data:`density_for_dft` * :c:data:`det_alpha_norm` * :c:data:`det_beta_norm` @@ -177,6 +184,7 @@ Index of Providers * :c:data:`do_direct_integrals` * :c:data:`do_pseudo` * :c:data:`do_pt2` +* :c:data:`do_third_order_1h1p` * :c:data:`double_exc_bitmask` * :c:data:`double_index_selectors` * :c:data:`dr_radial_integral` @@ -203,6 +211,8 @@ Index of Providers * :c:data:`energy_c` * :c:data:`energy_c_lda` * :c:data:`energy_c_pbe` +* :c:data:`energy_cas_dyall` +* :c:data:`energy_cas_dyall_no_exchange` * :c:data:`energy_iterations` * :c:data:`energy_sr_c_lda` * :c:data:`energy_sr_c_pbe` @@ -224,6 +234,12 @@ Index of Providers * :c:data:`final_grid_points` * :c:data:`final_weight_at_r` * :c:data:`final_weight_at_r_vector` +* :c:data:`final_weight_functions_at_final_grid_points` +* :c:data:`final_weight_functions_at_grid_points` +* :c:data:`fock_core_inactive` +* :c:data:`fock_core_inactive_from_act` +* :c:data:`fock_core_inactive_total` +* :c:data:`fock_core_inactive_total_spin_trace` * :c:data:`fock_matrix_alpha_no_xc_ao` * :c:data:`fock_matrix_ao` * :c:data:`fock_matrix_ao_alpha` @@ -234,7 +250,12 @@ Index of Providers * :c:data:`fock_matrix_mo` * :c:data:`fock_matrix_mo_alpha` * :c:data:`fock_matrix_mo_beta` +* :c:data:`fock_operator_active_from_core_inact` * :c:data:`fock_operator_closed_shell_ref_bitmask` +* :c:data:`fock_virt_from_act` +* :c:data:`fock_virt_from_core_inact` +* :c:data:`fock_virt_total` +* :c:data:`fock_virt_total_spin_trace` * :c:data:`fock_wee_closed_shell` * :c:data:`fps_spf_matrix_ao` * :c:data:`fps_spf_matrix_mo` @@ -242,12 +263,19 @@ Index of Providers * :c:data:`full_ijkl_bitmask_4` * :c:data:`gauleg_t2` * :c:data:`gauleg_w` +* :c:data:`gen_det_ref_idx` +* :c:data:`gen_det_ref_shortcut` +* :c:data:`gen_det_ref_sorted` +* :c:data:`gen_det_ref_version` * :c:data:`general_primitive_integral` * :c:data:`general_primitive_integral_erf` * :c:data:`generators_bitmask` * :c:data:`generators_bitmask_restart` * :c:data:`gga_sr_type_functionals` * :c:data:`gga_type_functionals` +* :c:data:`give_holes_and_particles_in_active_space` +* :c:data:`give_polynom_mult_center_mono_elec_erf` +* :c:data:`give_polynom_mult_center_mono_elec_erf_opt` * :c:data:`grad_aos_dsr_vc_alpha_pbe_w` * :c:data:`grad_aos_dsr_vc_beta_pbe_w` * :c:data:`grad_aos_dsr_vx_alpha_pbe_w` @@ -264,11 +292,15 @@ Index of Providers * :c:data:`h_apply_buffer_lock` * :c:data:`h_matrix_all_dets` * :c:data:`h_matrix_cas` +* :c:data:`h_matrix_ref` * :c:data:`hf_bitmask` * :c:data:`hf_energy` * :c:data:`hf_exchange` * :c:data:`hf_one_electron_energy` * :c:data:`hf_two_electron_energy` +* :c:data:`hmatrix_dressed_pt2_new` +* :c:data:`hmatrix_dressed_pt2_new_symmetrized` +* :c:data:`holes_operators` * :c:data:`i2radix_sort` * :c:data:`i8radix_sort` * :c:data:`i8radix_sort_big` @@ -284,6 +316,10 @@ Index of Providers * :c:data:`i_x2_pol_mult_mono_elec` * :c:data:`idx_cas` * :c:data:`idx_non_cas` +* :c:data:`idx_non_ref` +* :c:data:`idx_non_ref_from_sorted` +* :c:data:`idx_non_ref_rev` +* :c:data:`idx_ref` * :c:data:`inact_bitmask` * :c:data:`inact_virt_bitmask` * :c:data:`index_final_points` @@ -298,7 +334,10 @@ Index of Providers * :c:data:`insert_into_mo_integrals_map` * :c:data:`int_erf_3_index` * :c:data:`int_erf_3_index_exc` +* :c:data:`integral_density_alpha_knowles_becke_per_atom` +* :c:data:`integral_density_beta_knowles_becke_per_atom` * :c:data:`inv_int` +* :c:data:`inv_norm_psi_ref` * :c:data:`inv_selectors_coef_hf` * :c:data:`inv_selectors_coef_hf_squared` * :c:data:`iradix_sort` @@ -403,7 +442,6 @@ Index of Providers * :c:data:`n_det_max` * :c:data:`n_det_max_full` * :c:data:`n_det_non_cas` -* :c:data:`n_det_print_wf` * :c:data:`n_det_non_ref` * :c:data:`n_det_print_wf` * :c:data:`n_det_ref` @@ -434,6 +472,7 @@ Index of Providers * :c:data:`no_oa_or_av_opt` * :c:data:`no_vvv_integrals` * :c:data:`no_vvvv_integrals` +* :c:data:`norm_psi_ref` * :c:data:`nproc` * :c:data:`nthreads_davidson` * :c:data:`nucl_aos` @@ -455,19 +494,31 @@ Index of Providers * :c:data:`nucl_num` * :c:data:`nucl_num_shell_aos` * :c:data:`nuclear_repulsion` +* :c:data:`one_anhil` +* :c:data:`one_anhil_inact` +* :c:data:`one_anhil_one_creat` +* :c:data:`one_anhil_one_creat_inact_virt` +* :c:data:`one_anhil_one_creat_inact_virt_bis` +* :c:data:`one_anhil_one_creat_inact_virt_norm` +* :c:data:`one_body_dm_alpha_and_grad_at_r` * :c:data:`one_body_dm_alpha_ao_for_dft` * :c:data:`one_body_dm_alpha_at_r` * :c:data:`one_body_dm_ao_alpha` * :c:data:`one_body_dm_ao_beta` * :c:data:`one_body_dm_average_mo_for_dft` +* :c:data:`one_body_dm_beta_and_grad_at_r` * :c:data:`one_body_dm_beta_ao_for_dft` * :c:data:`one_body_dm_beta_at_r` * :c:data:`one_body_dm_dagger_mo_spin_index` * :c:data:`one_body_dm_mo` * :c:data:`one_body_dm_mo_alpha` +* :c:data:`one_body_dm_mo_alpha_and_grad_at_grid_points` +* :c:data:`one_body_dm_mo_alpha_at_grid_points` * :c:data:`one_body_dm_mo_alpha_average` * :c:data:`one_body_dm_mo_alpha_for_dft` * :c:data:`one_body_dm_mo_beta` +* :c:data:`one_body_dm_mo_beta_and_grad_at_grid_points` +* :c:data:`one_body_dm_mo_beta_at_grid_points` * :c:data:`one_body_dm_mo_beta_average` * :c:data:`one_body_dm_mo_beta_for_dft` * :c:data:`one_body_dm_mo_diff` @@ -479,6 +530,8 @@ Index of Providers * :c:data:`one_body_single_double_dm_mo_beta` * :c:data:`one_body_spin_density_ao` * :c:data:`one_body_spin_density_mo` +* :c:data:`one_creat` +* :c:data:`one_creat_virt` * :c:data:`one_dm_alpha_in_r` * :c:data:`one_dm_and_grad_alpha_in_r` * :c:data:`one_dm_and_grad_beta_in_r` @@ -488,6 +541,7 @@ Index of Providers * :c:data:`output_cpu_time_0` * :c:data:`output_wall_time_0` * :c:data:`overlap_gaussian_xyz` +* :c:data:`particles_operators` * :c:data:`phi_angular_integration_lebedev` * :c:data:`potential_c_alpha_ao` * :c:data:`potential_c_alpha_ao_lda` @@ -530,6 +584,7 @@ Index of Providers * :c:data:`pseudo_v_k_transp` * :c:data:`pseudo_v_kl` * :c:data:`pseudo_v_kl_transp` +* :c:data:`psi_active` * :c:data:`psi_average_norm_contrib` * :c:data:`psi_average_norm_contrib_sorted` * :c:data:`psi_bilinear_matrix` @@ -583,8 +638,30 @@ Index of Providers * :c:data:`psi_non_cas_coef` * :c:data:`psi_non_cas_coef_sorted_bit` * :c:data:`psi_non_cas_sorted_bit` +* :c:data:`psi_non_ref` +* :c:data:`psi_non_ref_coef` +* :c:data:`psi_non_ref_coef_interm_norm` +* :c:data:`psi_non_ref_coef_restart` +* :c:data:`psi_non_ref_coef_sorted_bit` +* :c:data:`psi_non_ref_coef_transp` +* :c:data:`psi_non_ref_restart` +* :c:data:`psi_non_ref_sorted_bit` * :c:data:`psi_occ_pattern` * :c:data:`psi_occ_pattern_hii` +* :c:data:`psi_ref` +* :c:data:`psi_ref_bis_lock` +* :c:data:`psi_ref_coef` +* :c:data:`psi_ref_coef_diagonalized` +* :c:data:`psi_ref_coef_interm_norm` +* :c:data:`psi_ref_coef_inv` +* :c:data:`psi_ref_coef_normalized` +* :c:data:`psi_ref_coef_restart` +* :c:data:`psi_ref_coef_sorted_bit` +* :c:data:`psi_ref_coef_transp` +* :c:data:`psi_ref_energy` +* :c:data:`psi_ref_energy_diagonalized` +* :c:data:`psi_ref_restart` +* :c:data:`psi_ref_sorted_bit` * :c:data:`psi_selectors` * :c:data:`psi_selectors_coef` * :c:data:`psi_selectors_coef_transp` @@ -625,12 +702,14 @@ Index of Providers * :c:data:`ref_bitmask` * :c:data:`ref_bitmask_energy` * :c:data:`ref_closed_shell_bitmask` +* :c:data:`ref_hamiltonian_matrix` * :c:data:`reunion_of_bitmask` * :c:data:`reunion_of_cas_inact_bitmask` * :c:data:`reunion_of_core_inact_act_bitmask` * :c:data:`reunion_of_core_inact_bitmask` * :c:data:`rs_ks_energy` * :c:data:`s2_eig` +* :c:data:`s2_matrix_all_dets` * :c:data:`s2_values` * :c:data:`s_half` * :c:data:`s_half_inv` @@ -643,6 +722,15 @@ Index of Providers * :c:data:`scf_density_matrix_ao_alpha` * :c:data:`scf_density_matrix_ao_beta` * :c:data:`scf_energy` +* :c:data:`second_order_pt_new` +* :c:data:`second_order_pt_new_1h` +* :c:data:`second_order_pt_new_1h1p` +* :c:data:`second_order_pt_new_1h2p` +* :c:data:`second_order_pt_new_1p` +* :c:data:`second_order_pt_new_2h` +* :c:data:`second_order_pt_new_2h1p` +* :c:data:`second_order_pt_new_2h2p` +* :c:data:`second_order_pt_new_2p` * :c:data:`select_max` * :c:data:`selection_criterion` * :c:data:`selection_criterion_factor` @@ -666,6 +754,8 @@ Index of Providers * :c:data:`state_following` * :c:data:`target_energy` * :c:data:`theta_angular_integration_lebedev` +* :c:data:`three_anhil` +* :c:data:`three_creat` * :c:data:`thresh_scf` * :c:data:`threshold_davidson` * :c:data:`threshold_diis` @@ -677,6 +767,10 @@ Index of Providers * :c:data:`trace_v_hxc` * :c:data:`trace_v_xc` * :c:data:`transpose` +* :c:data:`two_anhil` +* :c:data:`two_anhil_one_creat` +* :c:data:`two_creat` +* :c:data:`two_creat_one_anhil` * :c:data:`two_electron_energy` * :c:data:`unpaired_alpha_electrons` * :c:data:`used_weight` @@ -684,6 +778,7 @@ Index of Providers * :c:data:`virt_bitmask` * :c:data:`virt_bitmask_4` * :c:data:`weight_at_r` +* :c:data:`weight_functions_at_grid_points` * :c:data:`weight_occ_pattern` * :c:data:`weights_angular_integration_lebedev` * :c:data:`weights_angular_points` @@ -777,8 +872,11 @@ Index of Subroutines/Functions * :c:func:`connect_to_taskserver` * :c:func:`connected_to_ref` * :c:func:`connected_to_ref_by_mono` +* :c:func:`contrib_1h2p_dm_based` +* :c:func:`contrib_2h1p_dm_based` * :c:func:`copy_h_apply_buffer_to_wf` * :c:func:`copy_psi_bilinear_to_psi` +* :c:func:`coulomb_value_no_check` * :c:func:`create_guess` * :c:func:`create_microlist` * :c:func:`create_minilist` @@ -817,6 +915,8 @@ Index of Subroutines/Functions * :c:func:`diag_h_mat_elem_au0_h_au0` * :c:func:`diag_h_mat_elem_fock` * :c:func:`diag_h_mat_elem_monoelec` +* :c:func:`diag_h_mat_elem_no_elec_check` +* :c:func:`diag_h_mat_elem_no_elec_check_no_exchange` * :c:func:`diag_s_mat_elem` * :c:func:`diag_wee_mat_elem` * :c:func:`diagonalize_ci` @@ -859,6 +959,7 @@ Index of Subroutines/Functions * :c:func:`example_bitmask` * :c:func:`example_determinants` * :c:func:`example_determinants_psi_det` +* :c:func:`extract_ref` * :c:func:`extrapolate_data` * :c:func:`f_function_becke` * :c:func:`f_integral` @@ -870,6 +971,7 @@ Index of Subroutines/Functions * :c:func:`filter_connected` * :c:func:`filter_connected_i_h_psi0` * :c:func:`filter_not_connected` +* :c:func:`find_connections_previous` * :c:func:`find_rotation` * :c:func:`four_idx_transform` * :c:func:`g0d` @@ -911,6 +1013,8 @@ Index of Subroutines/Functions * :c:func:`get_d0` * :c:func:`get_d1` * :c:func:`get_d2` +* :c:func:`get_delta_e_dyall` +* :c:func:`get_delta_e_dyall_general_mp` * :c:func:`get_double_excitation` * :c:func:`get_double_excitation_spin` * :c:func:`get_excitation` @@ -925,6 +1029,7 @@ Index of Subroutines/Functions * :c:func:`get_index_in_psi_det_alpha_unique` * :c:func:`get_index_in_psi_det_beta_unique` * :c:func:`get_index_in_psi_det_sorted_bit` +* :c:func:`get_index_in_psi_ref_sorted_bit` * :c:func:`get_inverse` * :c:func:`get_m0` * :c:func:`get_m1` @@ -957,6 +1062,19 @@ Index of Subroutines/Functions * :c:func:`get_uj_s2_ui` * :c:func:`getmobiles` * :c:func:`getunitandopen` +* :c:func:`give_1h1p_contrib` +* :c:func:`give_1h1p_only_doubles_spin_cross` +* :c:func:`give_1h1p_sec_order_singles_contrib` +* :c:func:`give_1h2p_contrib` +* :c:func:`give_1h2p_contrib_sec_order` +* :c:func:`give_1h2p_new` +* :c:func:`give_1p_sec_order_singles_contrib` +* :c:func:`give_2h1p_contrib` +* :c:func:`give_2h1p_contrib_sec_order` +* :c:func:`give_2h1p_new` +* :c:func:`give_2h2p` +* :c:func:`give_2p_new` +* :c:func:`give_active_part_determinant` * :c:func:`give_all_aos_and_grad_and_lapl_at_r` * :c:func:`give_all_aos_and_grad_at_r` * :c:func:`give_all_aos_at_r` @@ -965,13 +1083,16 @@ Index of Subroutines/Functions * :c:func:`give_all_mos_and_grad_and_lapl_at_r` * :c:func:`give_all_mos_and_grad_at_r` * :c:func:`give_all_mos_at_r` +* :c:func:`give_core_inactive_part_determinant` * :c:func:`give_explicit_poly_and_gaussian` * :c:func:`give_explicit_poly_and_gaussian_double` * :c:func:`give_explicit_poly_and_gaussian_x` +* :c:func:`give_holes_in_inactive_space` +* :c:func:`give_particles_in_virt_space` * :c:func:`give_polynom_mult_center_mono_elec` -* :c:func:`give_polynom_mult_center_mono_elec_erf` -* :c:func:`give_polynom_mult_center_mono_elec_erf_opt` * :c:func:`give_polynom_mult_center_x` +* :c:func:`give_singles_and_partial_doubles_1h1p_contrib` +* :c:func:`give_virt_part_determinant` * :c:func:`gpw` * :c:func:`grad_rho_ab_to_grad_rho_oc` * :c:func:`gser` @@ -985,6 +1106,51 @@ Index of Subroutines/Functions * :c:func:`h_apply_cisd_diexcorg` * :c:func:`h_apply_cisd_diexcp` * :c:func:`h_apply_cisd_monoexc` +* :c:func:`h_apply_mrpt` +* :c:func:`h_apply_mrpt_1h` +* :c:func:`h_apply_mrpt_1h1p` +* :c:func:`h_apply_mrpt_1h1p_diexc` +* :c:func:`h_apply_mrpt_1h1p_diexcorg` +* :c:func:`h_apply_mrpt_1h1p_diexcp` +* :c:func:`h_apply_mrpt_1h1p_monoexc` +* :c:func:`h_apply_mrpt_1h2p` +* :c:func:`h_apply_mrpt_1h2p_diexc` +* :c:func:`h_apply_mrpt_1h2p_diexcorg` +* :c:func:`h_apply_mrpt_1h2p_diexcp` +* :c:func:`h_apply_mrpt_1h2p_monoexc` +* :c:func:`h_apply_mrpt_1h_diexc` +* :c:func:`h_apply_mrpt_1h_diexcorg` +* :c:func:`h_apply_mrpt_1h_diexcp` +* :c:func:`h_apply_mrpt_1h_monoexc` +* :c:func:`h_apply_mrpt_1p` +* :c:func:`h_apply_mrpt_1p_diexc` +* :c:func:`h_apply_mrpt_1p_diexcorg` +* :c:func:`h_apply_mrpt_1p_diexcp` +* :c:func:`h_apply_mrpt_1p_monoexc` +* :c:func:`h_apply_mrpt_2h` +* :c:func:`h_apply_mrpt_2h1p` +* :c:func:`h_apply_mrpt_2h1p_diexc` +* :c:func:`h_apply_mrpt_2h1p_diexcorg` +* :c:func:`h_apply_mrpt_2h1p_diexcp` +* :c:func:`h_apply_mrpt_2h1p_monoexc` +* :c:func:`h_apply_mrpt_2h2p` +* :c:func:`h_apply_mrpt_2h2p_diexc` +* :c:func:`h_apply_mrpt_2h2p_diexcorg` +* :c:func:`h_apply_mrpt_2h2p_diexcp` +* :c:func:`h_apply_mrpt_2h2p_monoexc` +* :c:func:`h_apply_mrpt_2h_diexc` +* :c:func:`h_apply_mrpt_2h_diexcorg` +* :c:func:`h_apply_mrpt_2h_diexcp` +* :c:func:`h_apply_mrpt_2h_monoexc` +* :c:func:`h_apply_mrpt_2p` +* :c:func:`h_apply_mrpt_2p_diexc` +* :c:func:`h_apply_mrpt_2p_diexcorg` +* :c:func:`h_apply_mrpt_2p_diexcp` +* :c:func:`h_apply_mrpt_2p_monoexc` +* :c:func:`h_apply_mrpt_diexc` +* :c:func:`h_apply_mrpt_diexcorg` +* :c:func:`h_apply_mrpt_diexcp` +* :c:func:`h_apply_mrpt_monoexc` * :c:func:`h_s2_u_0_bielec_nstates_openmp` * :c:func:`h_s2_u_0_bielec_nstates_openmp_work` * :c:func:`h_s2_u_0_bielec_nstates_openmp_work_1` @@ -1023,6 +1189,8 @@ Index of Subroutines/Functions * :c:func:`i_h_j_bielec` * :c:func:`i_h_j_double_alpha_beta` * :c:func:`i_h_j_double_spin` +* :c:func:`i_h_j_dyall` +* :c:func:`i_h_j_dyall_no_exchange` * :c:func:`i_h_j_mono_spin` * :c:func:`i_h_j_mono_spin_monoelec` * :c:func:`i_h_j_monoelec` @@ -1030,6 +1198,7 @@ Index of Subroutines/Functions * :c:func:`i_h_j_verbose` * :c:func:`i_h_psi` * :c:func:`i_h_psi_minilist` +* :c:func:`i_h_psi_pert_new_minilist` * :c:func:`i_s2_psi_minilist` * :c:func:`i_wee_j_mono` * :c:func:`i_x1_pol_mult` @@ -1062,6 +1231,7 @@ Index of Subroutines/Functions * :c:func:`is_connected_to_by_mono` * :c:func:`is_generable_cassd` * :c:func:`is_i_in_virtual` +* :c:func:`is_in_psi_ref` * :c:func:`is_in_wavefunction` * :c:func:`is_spin_flip_possible` * :c:func:`is_the_hole_in_det` @@ -1101,6 +1271,7 @@ Index of Subroutines/Functions * :c:func:`molden` * :c:func:`mono_excitation_wee` * :c:func:`mpi_print` +* :c:func:`mrpt_dress` * :c:func:`multiply_poly` * :c:func:`n_pt_sup` * :c:func:`nai_pol_mult` @@ -1160,7 +1331,6 @@ Index of Subroutines/Functions * :c:func:`print_generators_bitmasks_particles` * :c:func:`print_generators_bitmasks_particles_for_one_generator` * :c:func:`print_memory_usage` -* :c:func:`print_r2` * :c:func:`print_spindet` * :c:func:`print_summary` * :c:func:`print_wf` @@ -1249,6 +1419,7 @@ Index of Subroutines/Functions * :c:func:`set_bit_to_integer` * :c:func:`set_bitmask_hole_as_input` * :c:func:`set_bitmask_particl_as_input` +* :c:func:`set_generators_bitmasks_as_holes_and_particles` * :c:func:`set_natural_mos` * :c:func:`set_order` * :c:func:`set_order_big` @@ -1276,6 +1447,8 @@ Index of Subroutines/Functions * :c:func:`tasks_done_to_taskserver` * :c:func:`testteethbuilding` * :c:func:`total_memory` +* :c:func:`u0_h_dyall_u0` +* :c:func:`u0_h_dyall_u0_no_exchange` * :c:func:`u_0_h_u_0` * :c:func:`u_0_h_u_0_bielec` * :c:func:`u_0_s2_u_0` @@ -1368,4 +1541,4 @@ Index of Subroutines/Functions * :c:func:`zmq_put_psi_det_beta_unique` * :c:func:`zmq_put_psi_det_size` * :c:func:`zmq_selection` -* :c:func:`zmq_set_running` +* :c:func:`zmq_set_running` \ No newline at end of file diff --git a/src/ao_one_e_integrals/pot_ao_erf_ints.irp.f b/src/ao_one_e_integrals/pot_ao_erf_ints.irp.f index 85921dbc..9fafbce3 100644 --- a/src/ao_one_e_integrals/pot_ao_erf_ints.irp.f +++ b/src/ao_one_e_integrals/pot_ao_erf_ints.irp.f @@ -1,368 +1,371 @@ subroutine give_all_erf_kl_ao(integrals_ao,mu_in,C_center) - implicit none - BEGIN_DOC -! subroutine that returs all integrals over r of type erf(mu_in * | r-C_center |)/| r-C_center | - END_DOC - double precision, intent(in) :: mu_in,C_center(3) - double precision, intent(out) :: integrals_ao(ao_num,ao_num) - double precision :: NAI_pol_mult_erf_ao - integer :: i,j,l,k,m - do k = 1, ao_num - do m = 1, ao_num - integrals_ao(m,k) = NAI_pol_mult_erf_ao(m,k,mu_in,C_center) + implicit none + BEGIN_DOC + ! Subroutine that returns all integrals over $r$ of type + ! $\frac{ \erf(\mu * |r-R_C|) }{ |r-R_C| }$ + END_DOC + double precision, intent(in) :: mu_in,C_center(3) + double precision, intent(out) :: integrals_ao(ao_num,ao_num) + double precision :: NAI_pol_mult_erf_ao + integer :: i,j,l,k,m + do k = 1, ao_num + do m = 1, ao_num + integrals_ao(m,k) = NAI_pol_mult_erf_ao(m,k,mu_in,C_center) + enddo enddo - enddo end double precision function NAI_pol_mult_erf_ao(i_ao,j_ao,mu_in,C_center) - implicit none - BEGIN_DOC - ! computes the following integral : - ! int[-infty;+infty] dr AO_i_ao (r) AO_j_ao(r) erf(mu_in * | r-C_center |)/| r-C_center | - END_DOC - integer, intent(in) :: i_ao,j_ao - double precision, intent(in) :: mu_in, C_center(3) - integer :: i,j,num_A,num_B, power_A(3), power_B(3), n_pt_in - double precision :: A_center(3), B_center(3),integral, alpha,beta, NAI_pol_mult_erf - num_A = ao_nucl(i_ao) - power_A(1:3)= ao_power(i_ao,1:3) - A_center(1:3) = nucl_coord(num_A,1:3) - num_B = ao_nucl(j_ao) - power_B(1:3)= ao_power(j_ao,1:3) - B_center(1:3) = nucl_coord(num_B,1:3) - n_pt_in = n_pt_max_integrals - NAI_pol_mult_erf_ao = 0.d0 - do i = 1, ao_prim_num(i_ao) - alpha = ao_expo_ordered_transp(i,i_ao) - do j = 1, ao_prim_num(j_ao) - beta = ao_expo_ordered_transp(j,j_ao) - integral = NAI_pol_mult_erf(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in,mu_in) - NAI_pol_mult_erf_ao += integral * ao_coef_normalized_ordered_transp(j,j_ao)*ao_coef_normalized_ordered_transp(i,i_ao) + implicit none + BEGIN_DOC + ! Computes the following integral : + ! $\int_{-\infty}^{infty} dr \chi_i(r) \chi_j(r) \frac{\erf(\mu |r-R_C|)}{|r-R_C|}$. + END_DOC + integer, intent(in) :: i_ao,j_ao + double precision, intent(in) :: mu_in, C_center(3) + integer :: i,j,num_A,num_B, power_A(3), power_B(3), n_pt_in + double precision :: A_center(3), B_center(3),integral, alpha,beta + double precision :: NAI_pol_mult_erf + num_A = ao_nucl(i_ao) + power_A(1:3)= ao_power(i_ao,1:3) + A_center(1:3) = nucl_coord(num_A,1:3) + num_B = ao_nucl(j_ao) + power_B(1:3)= ao_power(j_ao,1:3) + B_center(1:3) = nucl_coord(num_B,1:3) + n_pt_in = n_pt_max_integrals + NAI_pol_mult_erf_ao = 0.d0 + do i = 1, ao_prim_num(i_ao) + alpha = ao_expo_ordered_transp(i,i_ao) + do j = 1, ao_prim_num(j_ao) + beta = ao_expo_ordered_transp(j,j_ao) + integral = NAI_pol_mult_erf(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in,mu_in) + NAI_pol_mult_erf_ao += integral * ao_coef_normalized_ordered_transp(j,j_ao)*ao_coef_normalized_ordered_transp(i,i_ao) + enddo enddo - enddo - - - + end double precision function NAI_pol_mult_erf(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in,mu_in) -! function that computes the folowing integral : -! int{dr} of (x-A_x)^ax (x-B_X)^bx exp(-alpha (x-A_x)^2 - beta (x-B_x)^2 ) erf(mu_in*(r-R_c))/(r-R_c) - -implicit none -integer, intent(in) :: n_pt_in -double precision,intent(in) :: C_center(3),A_center(3),B_center(3),alpha,beta,mu_in -integer, intent(in) :: power_A(3),power_B(3) -integer :: i,j,k,l,n_pt -double precision :: P_center(3) - -double precision :: d(0:n_pt_in),pouet,coeff,dist,const,pouet_2,factor -double precision :: I_n_special_exact,integrate_bourrin,I_n_bibi -double precision :: V_e_n,const_factor,dist_integral,tmp -double precision :: accu,rint,p_inv,p,rho,p_inv_2 -integer :: n_pt_out,lmax -include 'utils/constants.include.F' + BEGIN_DOC + ! Computes the following integral : + ! $\int dr (x-A_x)^a (x-B_x)^b \exp(-\alpha (x-A_x)^2 - \beta (x-B_x)^2 ) + ! \frac{\erf(\mu |r-R_C|)}{|r-R_c|}$. + END_DOC + + implicit none + integer, intent(in) :: n_pt_in + double precision,intent(in) :: C_center(3),A_center(3),B_center(3),alpha,beta,mu_in + integer, intent(in) :: power_A(3),power_B(3) + integer :: i,j,k,l,n_pt + double precision :: P_center(3) + + double precision :: d(0:n_pt_in),pouet,coeff,dist,const,pouet_2,factor + double precision :: I_n_special_exact,integrate_bourrin,I_n_bibi + double precision :: V_e_n,const_factor,dist_integral,tmp + double precision :: accu,rint,p_inv,p,rho,p_inv_2 + integer :: n_pt_out,lmax + include 'utils/constants.include.F' p = alpha + beta p_inv = 1.d0/p - p_inv_2 = 0.5d0 * p_inv + p_inv_2 = 0.5d0 * p_inv rho = alpha * beta * p_inv - + dist = 0.d0 dist_integral = 0.d0 do i = 1, 3 - P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv - dist += (A_center(i) - B_center(i))*(A_center(i) - B_center(i)) - dist_integral += (P_center(i) - C_center(i))*(P_center(i) - C_center(i)) + P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv + dist += (A_center(i) - B_center(i))*(A_center(i) - B_center(i)) + dist_integral += (P_center(i) - C_center(i))*(P_center(i) - C_center(i)) enddo - const_factor = dist*rho + const_factor = dist*rho if(const_factor > 80.d0)then - NAI_pol_mult_erf = 0.d0 - return + NAI_pol_mult_erf = 0.d0 + return endif - double precision :: p_new + double precision :: p_new p_new = mu_in/dsqrt(p+ mu_in * mu_in) factor = dexp(-const_factor) - coeff = dtwo_pi * factor * p_inv * p_new + coeff = dtwo_pi * factor * p_inv * p_new lmax = 20 - ! print*, "b" + ! print*, "b" do i = 0, n_pt_in d(i) = 0.d0 enddo n_pt = 2 * ( (power_A(1) + power_B(1)) +(power_A(2) + power_B(2)) +(power_A(3) + power_B(3)) ) - const = p * dist_integral * p_new * p_new + const = p * dist_integral * p_new * p_new if (n_pt == 0) then - pouet = rint(0,const) - NAI_pol_mult_erf = coeff * pouet - return + pouet = rint(0,const) + NAI_pol_mult_erf = coeff * pouet + return endif - -! call give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in) - p_new = p_new * p_new + + ! call give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in) + p_new = p_new * p_new call give_polynom_mult_center_mono_elec_erf_opt(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in,p,p_inv,p_inv_2,p_new,P_center) - + if(n_pt_out<0)then - NAI_pol_mult_erf = 0.d0 - return + NAI_pol_mult_erf = 0.d0 + return endif accu = 0.d0 - -! sum of integrals of type : int {t,[0,1]} exp-(rho.(P-Q)^2 * t^2) * t^i + + ! sum of integrals of type : int {t,[0,1]} exp-(rho.(P-Q)^2 * t^2) * t^i do i =0 ,n_pt_out,2 - accu += d(i) * rint(i/2,const) + accu += d(i) * rint(i/2,const) enddo NAI_pol_mult_erf = accu * coeff - + end -subroutine give_polynom_mult_center_mono_elec_erf_opt(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in,p,p_inv,p_inv_2,p_new,P_center) -!!!! subroutine that returns the explicit polynom in term of the "t" variable of the following polynomw :: -!!!! I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q) -!!!! it is for the nuclear electron atraction -implicit none -integer, intent(in) :: n_pt_in -integer,intent(out) :: n_pt_out -double precision, intent(in) :: A_center(3), B_center(3),C_center(3),p,p_inv,p_inv_2,p_new,P_center(3) -double precision, intent(in) :: alpha,beta,mu_in -integer, intent(in) :: power_A(3), power_B(3) -integer :: a_x,b_x,a_y,b_y,a_z,b_z -double precision :: d(0:n_pt_in) -double precision :: d1(0:n_pt_in) -double precision :: d2(0:n_pt_in) -double precision :: d3(0:n_pt_in) -double precision :: accu - accu = 0.d0 -!COMPTEUR irp_rdtsc1 = irp_rdtsc() - ASSERT (n_pt_in > 1) - - double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2) - R1x(0) = (P_center(1) - A_center(1)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(1) - C_center(1))* p_new - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(1) - B_center(1)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(1) - C_center(1))* p_new - !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R2x(0) = p_inv_2 - R2x(1) = 0.d0 - R2x(2) = -p_inv_2* p_new - !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - do i = 0,n_pt_in - d1(i) = 0.d0 - enddo - do i = 0,n_pt_in - d2(i) = 0.d0 - enddo - do i = 0,n_pt_in - d3(i) = 0.d0 - enddo - integer :: n_pt1,n_pt2,n_pt3,dim,i - n_pt1 = n_pt_in - n_pt2 = n_pt_in - n_pt3 = n_pt_in - a_x = power_A(1) - b_x = power_B(1) - call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in) -! print*,'passed the first I_x1' +subroutine give_polynom_mult_center_mono_elec_erf_opt(A_center,B_center,alpha,beta,& + power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in,p,p_inv,p_inv_2,p_new,P_center) + BEGIN_DOC + ! Returns the explicit polynomial in terms of the $t$ variable of the following polynomial: + ! + ! $I_{x1}(a_x, d_x,p,q) \times I_{x1}(a_y, d_y,p,q) \times I_{x1}(a_z, d_z,p,q)$. + END_DOC + implicit none + integer, intent(in) :: n_pt_in + integer,intent(out) :: n_pt_out + double precision, intent(in) :: A_center(3), B_center(3),C_center(3),p,p_inv,p_inv_2,p_new,P_center(3) + double precision, intent(in) :: alpha,beta,mu_in + integer, intent(in) :: power_A(3), power_B(3) + integer :: a_x,b_x,a_y,b_y,a_z,b_z + double precision :: d(0:n_pt_in) + double precision :: d1(0:n_pt_in) + double precision :: d2(0:n_pt_in) + double precision :: d3(0:n_pt_in) + double precision :: accu + accu = 0.d0 + ASSERT (n_pt_in > 1) + + double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2) + R1x(0) = (P_center(1) - A_center(1)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(1) - C_center(1))* p_new + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(1) - B_center(1)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(1) - C_center(1))* p_new + !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R2x(0) = p_inv_2 + R2x(1) = 0.d0 + R2x(2) = -p_inv_2* p_new + !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + do i = 0,n_pt_in + d1(i) = 0.d0 + enddo + do i = 0,n_pt_in + d2(i) = 0.d0 + enddo + do i = 0,n_pt_in + d3(i) = 0.d0 + enddo + integer :: n_pt1,n_pt2,n_pt3,dim,i + n_pt1 = n_pt_in + n_pt2 = n_pt_in + n_pt3 = n_pt_in + a_x = power_A(1) + b_x = power_B(1) + call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in) if(n_pt1<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - - R1x(0) = (P_center(2) - A_center(2)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(2) - C_center(2))* p_new - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(2) - B_center(2)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(2) - C_center(2))* p_new - !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - a_y = power_A(2) - b_y = power_B(2) - call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in) -! print*,'passed the second I_x1' + + R1x(0) = (P_center(2) - A_center(2)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(2) - C_center(2))* p_new + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(2) - B_center(2)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(2) - C_center(2))* p_new + !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + a_y = power_A(2) + b_y = power_B(2) + call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in) if(n_pt2<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - - - R1x(0) = (P_center(3) - A_center(3)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(3) - C_center(3))* p_new - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(3) - B_center(3)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(3) - C_center(3))* p_new - !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 - a_z = power_A(3) - b_z = power_B(3) - -! print*,'a_z = ',a_z -! print*,'b_z = ',b_z - call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in) -! print*,'passed the third I_x1' + + + R1x(0) = (P_center(3) - A_center(3)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(3) - C_center(3))* p_new + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(3) - B_center(3)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(3) - C_center(3))* p_new + !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 + a_z = power_A(3) + b_z = power_B(3) + + call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in) if(n_pt3<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - integer :: n_pt_tmp - n_pt_tmp = 0 - call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp) - do i = 0,n_pt_tmp - d1(i) = 0.d0 - enddo - n_pt_out = 0 - call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out) - do i = 0, n_pt_out - d(i) = d1(i) - enddo - + integer :: n_pt_tmp + n_pt_tmp = 0 + call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp) + do i = 0,n_pt_tmp + d1(i) = 0.d0 + enddo + n_pt_out = 0 + call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out) + do i = 0, n_pt_out + d(i) = d1(i) + enddo + end -subroutine give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in) -!!!! subroutine that returns the explicit polynom in term of the "t" variable of the following polynomw :: -!!!! I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q) -!!!! it is for the nuclear electron atraction -implicit none -integer, intent(in) :: n_pt_in -integer,intent(out) :: n_pt_out -double precision, intent(in) :: A_center(3), B_center(3),C_center(3) -double precision, intent(in) :: alpha,beta,mu_in -integer, intent(in) :: power_A(3), power_B(3) -integer :: a_x,b_x,a_y,b_y,a_z,b_z -double precision :: d(0:n_pt_in) -double precision :: d1(0:n_pt_in) -double precision :: d2(0:n_pt_in) -double precision :: d3(0:n_pt_in) -double precision :: accu, pq_inv, p10_1, p10_2, p01_1, p01_2 -double precision :: p,P_center(3),rho,p_inv,p_inv_2 - accu = 0.d0 -!COMPTEUR irp_rdtsc1 = irp_rdtsc() - ASSERT (n_pt_in > 1) - p = alpha+beta - p_inv = 1.d0/p - p_inv_2 = 0.5d0/p - do i =1, 3 - P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv - enddo - - double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2) - R1x(0) = (P_center(1) - A_center(1)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in) - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(1) - B_center(1)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in) - !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R2x(0) = p_inv_2 - R2x(1) = 0.d0 - R2x(2) = -p_inv_2* mu_in**2 / (p+mu_in*mu_in) - !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - do i = 0,n_pt_in - d1(i) = 0.d0 - enddo - do i = 0,n_pt_in - d2(i) = 0.d0 - enddo - do i = 0,n_pt_in - d3(i) = 0.d0 - enddo - integer :: n_pt1,n_pt2,n_pt3,dim,i - n_pt1 = n_pt_in - n_pt2 = n_pt_in - n_pt3 = n_pt_in - a_x = power_A(1) - b_x = power_B(1) - call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in) -! print*,'passed the first I_x1' +subroutine give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,& + power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in) + BEGIN_DOC + ! Returns the explicit polynomial in terms of the $t$ variable of the following polynomial: + ! + ! $I_{x1}(a_x, d_x,p,q) \times I_{x1}(a_y, d_y,p,q) \times I_{x1}(a_z, d_z,p,q)$. + END_DOC + implicit none + integer, intent(in) :: n_pt_in + integer,intent(out) :: n_pt_out + double precision, intent(in) :: A_center(3), B_center(3),C_center(3) + double precision, intent(in) :: alpha,beta,mu_in + integer, intent(in) :: power_A(3), power_B(3) + integer :: a_x,b_x,a_y,b_y,a_z,b_z + double precision :: d(0:n_pt_in) + double precision :: d1(0:n_pt_in) + double precision :: d2(0:n_pt_in) + double precision :: d3(0:n_pt_in) + double precision :: accu, pq_inv, p10_1, p10_2, p01_1, p01_2 + double precision :: p,P_center(3),rho,p_inv,p_inv_2 + accu = 0.d0 + !COMPTEUR irp_rdtsc1 = irp_rdtsc() + ASSERT (n_pt_in > 1) + p = alpha+beta + p_inv = 1.d0/p + p_inv_2 = 0.5d0/p + do i =1, 3 + P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv + enddo + + double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2) + R1x(0) = (P_center(1) - A_center(1)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in) + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(1) - B_center(1)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in) + !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R2x(0) = p_inv_2 + R2x(1) = 0.d0 + R2x(2) = -p_inv_2* mu_in**2 / (p+mu_in*mu_in) + !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + do i = 0,n_pt_in + d1(i) = 0.d0 + enddo + do i = 0,n_pt_in + d2(i) = 0.d0 + enddo + do i = 0,n_pt_in + d3(i) = 0.d0 + enddo + integer :: n_pt1,n_pt2,n_pt3,dim,i + n_pt1 = n_pt_in + n_pt2 = n_pt_in + n_pt3 = n_pt_in + a_x = power_A(1) + b_x = power_B(1) + call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in) + ! print*,'passed the first I_x1' if(n_pt1<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - - R1x(0) = (P_center(2) - A_center(2)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in) - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(2) - B_center(2)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in) - !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - a_y = power_A(2) - b_y = power_B(2) - call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in) -! print*,'passed the second I_x1' + + R1x(0) = (P_center(2) - A_center(2)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in) + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(2) - B_center(2)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in) + !R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + a_y = power_A(2) + b_y = power_B(2) + call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in) + ! print*,'passed the second I_x1' if(n_pt2<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - - - R1x(0) = (P_center(3) - A_center(3)) - R1x(1) = 0.d0 - R1x(2) = -(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in) - ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 - R1xp(0) = (P_center(3) - B_center(3)) - R1xp(1) = 0.d0 - R1xp(2) =-(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in) - !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 - a_z = power_A(3) - b_z = power_B(3) - -! print*,'a_z = ',a_z -! print*,'b_z = ',b_z - call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in) -! print*,'passed the third I_x1' + + + R1x(0) = (P_center(3) - A_center(3)) + R1x(1) = 0.d0 + R1x(2) = -(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in) + ! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2 + R1xp(0) = (P_center(3) - B_center(3)) + R1xp(1) = 0.d0 + R1xp(2) =-(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in) + !R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2 + a_z = power_A(3) + b_z = power_B(3) + + ! print*,'a_z = ',a_z + ! print*,'b_z = ',b_z + call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in) + ! print*,'passed the third I_x1' if(n_pt3<0)then - n_pt_out = -1 - do i = 0,n_pt_in - d(i) = 0.d0 - enddo - return + n_pt_out = -1 + do i = 0,n_pt_in + d(i) = 0.d0 + enddo + return endif - integer :: n_pt_tmp - n_pt_tmp = 0 - call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp) - do i = 0,n_pt_tmp - d1(i) = 0.d0 - enddo - n_pt_out = 0 - call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out) - do i = 0, n_pt_out - d(i) = d1(i) - enddo - + integer :: n_pt_tmp + n_pt_tmp = 0 + call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp) + do i = 0,n_pt_tmp + d1(i) = 0.d0 + enddo + n_pt_out = 0 + call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out) + do i = 0, n_pt_out + d(i) = d1(i) + enddo + end diff --git a/src/ao_one_e_integrals/pot_ao_ints.irp.f b/src/ao_one_e_integrals/pot_ao_ints.irp.f index 2229d2ec..34e046ae 100644 --- a/src/ao_one_e_integrals/pot_ao_ints.irp.f +++ b/src/ao_one_e_integrals/pot_ao_ints.irp.f @@ -237,9 +237,9 @@ end subroutine give_polynom_mult_center_mono_elec(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out) implicit none BEGIN_DOC -! Returns the explicit polynomial in terms of the "t" variable of the following -! -! :math:`I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q)` + ! Returns the explicit polynomial in terms of the "t" variable of the following + ! + ! $I_{x1}(a_x, d_x,p,q) \times I_{x1}(a_y, d_y,p,q) \times I_{x1}(a_z, d_z,p,q)$. END_DOC integer, intent(in) :: n_pt_in integer,intent(out) :: n_pt_out @@ -481,9 +481,9 @@ double precision function V_e_n(a_x,a_y,a_z,b_x,b_y,b_z,alpha,beta) BEGIN_DOC ! Primitve nuclear attraction between the two primitves centered on the same atom. ! -! primitive_1 = x**(a_x) y**(a_y) z**(a_z) exp(-alpha * r**2) +! $p_1 = x^{a_x} y^{a_y} z^{a_z} \exp(-\alpha r^2)$ ! -! primitive_2 = x**(b_x) y**(b_y) z**(b_z) exp(- beta * r**2) +! $p_2 = x^{b_x} y^{b_y} z^{b_z} \exp(-\beta r^2)$ END_DOC integer :: a_x,a_y,a_z,b_x,b_y,b_z double precision :: alpha,beta @@ -504,7 +504,7 @@ double precision function int_gaus_pol(alpha,n) BEGIN_DOC ! Computes the integral: ! -! :math:`\int_{-\infty}^{\infty} x^n \exp(-\alpha x^2) dx` +! $\int_{-\infty}^{\infty} x^n \exp(-\alpha x^2) dx$. END_DOC double precision :: alpha integer :: n @@ -530,7 +530,7 @@ double precision function V_r(n,alpha) BEGIN_DOC ! Computes the radial part of the nuclear attraction integral: ! - ! :math:`\int_{0}^{\infty} r^n \exp(-\alpha r^2) dr` + ! $\int_{0}^{\infty} r^n \exp(-\alpha r^2) dr$ ! END_DOC double precision :: alpha, fact @@ -547,9 +547,9 @@ end double precision function V_phi(n,m) implicit none BEGIN_DOC - ! Computes the angular "phi" part of the nuclear attraction integral: + ! Computes the angular $\phi$ part of the nuclear attraction integral: ! - ! :math:`\int_{0}^{2 \pi} \cos(\phi)^n \sin(\phi)^m d\phi` + ! $\int_{0}^{2 \pi} \cos(\phi)^n \sin(\phi)^m d\phi$. END_DOC integer :: n,m, i double precision :: prod, Wallis @@ -564,9 +564,9 @@ end double precision function V_theta(n,m) implicit none BEGIN_DOC - ! Computes the angular "theta" part of the nuclear attraction integral: + ! Computes the angular $\theta$ part of the nuclear attraction integral: ! - ! :math:`\int_{0}^{\pi} \cos(\theta)^n \sin(\theta)^m d\theta` + ! $\int_{0}^{\pi} \cos(\theta)^n \sin(\theta)^m d\theta$ END_DOC integer :: n,m,i double precision :: Wallis, prod @@ -585,7 +585,7 @@ double precision function Wallis(n) BEGIN_DOC ! Wallis integral: ! - ! :math:`\int_{0}^{\pi} \cos(\theta)^n d\theta` + ! $\int_{0}^{\pi} \cos(\theta)^n d\theta$. END_DOC double precision :: fact integer :: n,p diff --git a/src/ao_one_e_integrals/pot_ao_pseudo_ints.irp.f b/src/ao_one_e_integrals/pot_ao_pseudo_ints.irp.f index 95703380..382a47b0 100644 --- a/src/ao_one_e_integrals/pot_ao_pseudo_ints.irp.f +++ b/src/ao_one_e_integrals/pot_ao_pseudo_ints.irp.f @@ -104,13 +104,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num,ao_num)] pseudo_dz_k_transp(1,k), & A_center,power_A,alpha,B_center,power_B,beta,C_center) -!if ((k==nucl_num).and.(num_A == nucl_num).and.(num_B == nucl_num)) then -!print *, pseudo_klocmax,pseudo_v_k_transp (1,k),pseudo_n_k_transp (1,k),pseudo_dz_k_transp(1,k) -!print *, A_center(1:3), power_A -!print *, B_center(1:3), power_B -!print *, C_center(1:3) -!print *, c -!endif enddo ao_pseudo_integral_local(i,j) = ao_pseudo_integral_local(i,j) +& ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i)*c diff --git a/src/cisd/40.cisd.bats b/src/cisd/40.cisd.bats index 5279f68a..0d4e57d4 100644 --- a/src/cisd/40.cisd.bats +++ b/src/cisd/40.cisd.bats @@ -70,7 +70,7 @@ function run() { @test "CH4" { qp_set_mo_class ch4.ezfio -core "[1]" -act "[2-30]" -del "[31-59]" - run ch4.ezfio -40.1996180778616 -40.2403950634044 + run ch4.ezfio -40.2403957502264 -39.8433160343413 } @test "CO2" { diff --git a/src/mo_basis/mos.irp.f b/src/mo_basis/mos.irp.f index eceb93b8..4e88e91c 100644 --- a/src/mo_basis/mos.irp.f +++ b/src/mo_basis/mos.irp.f @@ -62,7 +62,9 @@ BEGIN_PROVIDER [ double precision, mo_coef, (ao_num,mo_tot_num) ] implicit none BEGIN_DOC ! Molecular orbital coefficients on AO basis set + ! ! mo_coef(i,j) = coefficient of the ith ao on the jth mo + ! ! mo_label : Label characterizing the MOS (local, canonical, natural, etc) END_DOC integer :: i, j @@ -115,7 +117,7 @@ BEGIN_PROVIDER [ double precision, mo_coef_in_ao_ortho_basis, (ao_num, mo_tot_nu BEGIN_DOC ! MO coefficients in orthogonalized AO basis ! - ! C^(-1).C_mo + ! $C^(-1).C_mo$ END_DOC call dgemm('N','N',ao_num,mo_tot_num,ao_num,1.d0, & ao_ortho_canonical_coef_inv, size(ao_ortho_canonical_coef_inv,1),&