PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-06-22 05:02:22 +02:00
Code Releases Activity

Merge branch 'dev-lct' into dev-lcpq

Browse Source
This commit is contained in:
Anthony Scemama 2019-05-15 16:17:26 +02:00
commit 9f7b159fe1
22 changed files with 307 additions and 58 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

79
bin/qp_set_frozen_large_core Executable file
View File

@ -0,0 +1,79 @@
#!/usr/bin/env python2
"""
Automatically finds n, the number of core electrons. Calls qp_set_mo_class
setting all MOs as Active, except the n/2 first ones which are set as Core.
If pseudo-potentials are used, all the MOs are set as Active.
Usage:
qp_set_frozen_core [-q|--query] EZFIO_DIR
Options:
-q --query Prints in the standard output the number of frozen MOs
"""
import os
import sys
import os.path
try:
import qp_path
except ImportError:
print "source .quantum_package.rc"
raise
from docopt import docopt
from ezfio import ezfio
def main(arguments):
"""Main function"""
filename = arguments["EZFIO_DIR"]
ezfio.set_filename(filename)
n_frozen = 0
try:
do_pseudo = ezfio.pseudo_do_pseudo
except:
do_pseudo = False
if not do_pseudo:
for charge in ezfio.nuclei_nucl_charge:
if charge <= 2:
pass
elif charge <= 10:
n_frozen += 1
elif charge <= 18:
n_frozen += 5
elif charge <= 36:
n_frozen += 9
elif charge <= 54:
n_frozen += 18
elif charge <= 86:
n_frozen += 27
elif charge <= 118:
n_frozen += 43
mo_num = ezfio.mo_basis_mo_num
if arguments["--query"]:
print n_frozen
sys.exit(0)
if n_frozen == 0:
os.system("""qp_set_mo_class -a "[1-%d]" %s""" %
(mo_num, sys.argv[1]))
else:
os.system("""qp_set_mo_class -c "[1-%d]" -a "[%d-%d]" %s""" %
(n_frozen, n_frozen+1, mo_num, sys.argv[1]))
if __name__ == '__main__':
ARGUMENTS = docopt(__doc__)
main(ARGUMENTS)

2
configure vendored
View File

@ -387,7 +387,7 @@ if [[ ${ZEROMQ} = $(not_found) ]] ; then
fail
fi
F77ZMQ=$(find_lib -lzmq -lf77zmq)
F77ZMQ=$(find_lib -lzmq -lf77zmq -lpthread)
if [[ ${F77ZMQ} = $(not_found) ]] ; then
error "Fortran binding of ZeroMQ (f77zmq) is not installed."
fail

1
etc/.gitignore vendored
View File

@ -1 +1,2 @@
00.qp_root
local.rc

1
etc/local.rc
View File

@ -16,6 +16,7 @@
# export OMP_NUM_THREADS=16
# Name of the network interface to be chosen
# export QP_NIC=lo
# export QP_NIC=ib0

5
ocaml/MO_label.ml
View File

@ -6,6 +6,7 @@ type t =
| Natural
| Localized
| Orthonormalized
| MCSCF
| None
[@@deriving sexp]
@ -16,6 +17,7 @@ let to_string = function
| Orthonormalized -> "Orthonormalized"
| Natural -> "Natural"
| Localized -> "Localized"
| MCSCF -> "MCSCF"
| None -> "None"
;;
@ -26,7 +28,8 @@ let of_string s =
| "natural" -> Natural
| "localized" -> Localized
| "orthonormalized" -> Orthonormalized
| "mcscf" -> MCSCF
| "none" -> None
| _ -> (print_endline s ; failwith "MO_label should be one of:
Guess | Orthonormalized | Canonical | Natural | Localized | None.")
Guess | Orthonormalized | Canonical | Natural | Localized | MCSCF | None.")
;;

1
ocaml/MO_label.mli
View File

@ -4,6 +4,7 @@ type t =
| Natural
| Localized
| Orthonormalized
| MCSCF
| None
[@@deriving sexp]

48
scripts/generate_h_apply.py
View File

@ -5,6 +5,8 @@ import os
keywords = """
check_double_excitation
copy_buffer
filter_only_connected_to_hf_single
filter_only_connected_to_hf_double
declarations
decls_main
deinit_thread
@ -205,84 +207,84 @@ class H_apply(object):
def filter_only_2h(self):
self["only_2h_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_2h(hole).eqv. .False.) cycle
if (.not.is_a_2h(hole)) cycle
"""
self["only_2h_double"] = """
! ! DIR$ FORCEINLINE
if ( is_a_2h(key).eqv. .False. )cycle
if (.not.is_a_2h(key))cycle
"""
def filter_only_1h(self):
self["only_1h_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h(hole) .eqv. .False.) cycle
if (.not.is_a_1h(hole)) cycle
"""
self["only_1h_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h(key) .eqv. .False.) cycle
if (.not.is_a_1h(key) ) cycle
"""
def filter_only_1p(self):
self["only_1p_single"] = """
! ! DIR$ FORCEINLINE
if ( is_a_1p(hole) .eqv. .False.) cycle
if (.not. is_a_1p(hole) ) cycle
"""
self["only_1p_double"] = """
! ! DIR$ FORCEINLINE
if ( is_a_1p(key) .eqv. .False.) cycle
if (.not. is_a_1p(key) ) cycle
"""
def filter_only_2h1p(self):
self["only_2h1p_single"] = """
! ! DIR$ FORCEINLINE
if ( is_a_2h1p(hole) .eqv. .False.) cycle
if (.not. is_a_2h1p(hole) ) cycle
"""
self["only_2h1p_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_2h1p(key) .eqv. .False.) cycle
if (.not.is_a_2h1p(key) ) cycle
"""
def filter_only_2p(self):
self["only_2p_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_2p(hole).eqv. .False.) cycle
if (.not.is_a_2p(hole)) cycle
"""
self["only_2p_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_2p(key).eqv. .False.) cycle
if (.not.is_a_2p(key)) cycle
"""
def filter_only_1h1p(self):
self["filter_only_1h1p_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h1p(hole).eqv..False.) cycle
if (.not.is_a_1h1p(hole)) cycle
"""
self["filter_only_1h1p_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h1p(key).eqv..False.) cycle
if (.not.is_a_1h1p(key)) cycle
"""
def filter_only_2h2p(self):
self["filter_only_2h2p_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_two_holes_two_particles(hole).eqv..False.) cycle
if (.not.is_a_two_holes_two_particles(hole)) cycle
"""
self["filter_only_2h2p_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_two_holes_two_particles(key).eqv..False.) cycle
if (.not.is_a_two_holes_two_particles(key)) cycle
"""
def filter_only_1h2p(self):
self["filter_only_1h2p_single"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h2p(hole).eqv..False.) cycle
if (.not.is_a_1h2p(hole)) cycle
"""
self["filter_only_1h2p_double"] = """
! ! DIR$ FORCEINLINE
if (is_a_1h2p(key).eqv..False.) cycle
if (.not.is_a_1h2p(key)) cycle
"""
@ -294,6 +296,16 @@ class H_apply(object):
if (is_a_two_holes_two_particles(hole)) cycle
"""
def filter_only_connected_to_hf(self):
self["filter_only_connected_to_hf_single"] = """
call connected_to_hf(hole,yes_no)
if (.not.yes_no) cycle
"""
self["filter_only_connected_to_hf_double"] = """
call connected_to_hf(key,yes_no)
if (.not.yes_no) cycle
"""
def set_perturbation(self,pert):
if self.perturbation is not None:

1
src/ao_two_e_ints/map_integrals.irp.f
View File

@ -213,6 +213,7 @@ subroutine get_ao_two_e_integrals(j,k,l,sze,out_val)
BEGIN_DOC
! Gets multiple AO bi-electronic integral from the AO map .
! All i are retrieved for j,k,l fixed.
! physicist convention : <ij|kl>
END_DOC
implicit none
integer, intent(in) :: j,k,l, sze

14
src/aux_quantities/EZFIO.cfg
View File

@ -24,3 +24,17 @@ type: double precision
size: (mo_basis.mo_num,mo_basis.mo_num,determinants.n_states)
[data_one_e_dm_alpha_ao]
interface: ezfio, provider
doc: Alpha one body density matrix on the |AO| basis computed with the wave function
type: double precision
size: (ao_basis.ao_num,ao_basis.ao_num,determinants.n_states)
[data_one_e_dm_beta_ao]
interface: ezfio, provider
doc: Beta one body density matrix on the |AO| basis computed with the wave function
type: double precision
size: (ao_basis.ao_num,ao_basis.ao_num,determinants.n_states)

26
src/bitmask/bitmasks.irp.f
View File

@ -466,35 +466,17 @@ END_PROVIDER
BEGIN_PROVIDER [integer(bit_kind), reunion_of_core_inact_act_bitmask, (N_int,2)]
&BEGIN_PROVIDER [ integer, n_core_inact_act_orb ]
implicit none
BEGIN_DOC
! Reunion of the core, inactive and active bitmasks
END_DOC
integer :: i,j
n_core_inact_act_orb = 0
do i = 1, N_int
reunion_of_core_inact_act_bitmask(i,1) = ior(reunion_of_core_inact_bitmask(i,1),cas_bitmask(i,1,1))
reunion_of_core_inact_act_bitmask(i,2) = ior(reunion_of_core_inact_bitmask(i,2),cas_bitmask(i,2,1))
n_core_inact_act_orb +=popcnt(reunion_of_core_inact_act_bitmask(i,1))
reunion_of_core_inact_act_bitmask(i,1) = ior(reunion_of_core_inact_bitmask(i,1),act_bitmask(i,1))
reunion_of_core_inact_act_bitmask(i,2) = ior(reunion_of_core_inact_bitmask(i,2),act_bitmask(i,2))
enddo
END_PROVIDER
BEGIN_PROVIDER [ integer, list_core_inact_act, (n_core_inact_act_orb)]
&BEGIN_PROVIDER [ integer, list_core_inact_act_reverse, (mo_num)]
implicit none
integer :: occ_inact(N_int*bit_kind_size)
integer :: itest,i
occ_inact = 0
call bitstring_to_list(reunion_of_core_inact_act_bitmask(1,1), occ_inact(1), itest, N_int)
list_inact_reverse = 0
do i = 1, n_core_inact_act_orb
list_core_inact_act(i) = occ_inact(i)
list_core_inact_act_reverse(occ_inact(i)) = i
enddo
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), reunion_of_bitmask, (N_int,2)]
@ -563,8 +545,8 @@ END_PROVIDER
END_DOC
integer :: i,j
do i = 1, N_int
reunion_of_cas_inact_bitmask(i,1) = ior(cas_bitmask(i,1,1),inact_bitmask(i,1))
reunion_of_cas_inact_bitmask(i,2) = ior(cas_bitmask(i,2,1),inact_bitmask(i,2))
reunion_of_cas_inact_bitmask(i,1) = ior(act_bitmask(i,1),inact_bitmask(i,1))
reunion_of_cas_inact_bitmask(i,2) = ior(act_bitmask(i,2),inact_bitmask(i,2))
enddo
END_PROVIDER

50
src/bitmask/core_inact_act_virt.irp.f
View File

@ -194,3 +194,53 @@ END_PROVIDER
END_PROVIDER
BEGIN_PROVIDER [integer, n_inact_act_orb ]
implicit none
n_inact_act_orb = (n_inact_orb+n_act_orb)
END_PROVIDER
BEGIN_PROVIDER [integer, list_inact_act, (n_inact_act_orb)]
integer :: i,itmp
itmp = 0
do i = 1, n_inact_orb
itmp += 1
list_inact_act(itmp) = list_inact(i)
enddo
do i = 1, n_act_orb
itmp += 1
list_inact_act(itmp) = list_act(i)
enddo
END_PROVIDER
BEGIN_PROVIDER [integer, n_core_inact_act_orb ]
implicit none
n_core_inact_act_orb = (n_core_orb + n_inact_orb + n_act_orb)
END_PROVIDER
BEGIN_PROVIDER [integer, list_core_inact_act, (n_core_inact_act_orb)]
&BEGIN_PROVIDER [ integer, list_core_inact_act_reverse, (n_core_inact_act_orb)]
integer :: i,itmp
itmp = 0
do i = 1, n_core_orb
itmp += 1
list_core_inact_act(itmp) = list_core(i)
enddo
do i = 1, n_inact_orb
itmp += 1
list_core_inact_act(itmp) = list_inact(i)
enddo
do i = 1, n_act_orb
itmp += 1
list_core_inact_act(itmp) = list_act(i)
enddo
integer :: occ_inact(N_int*bit_kind_size)
occ_inact = 0
call bitstring_to_list(reunion_of_core_inact_act_bitmask(1,1), occ_inact(1), itest, N_int)
list_inact_reverse = 0
do i = 1, n_core_inact_act_orb
list_core_inact_act_reverse(occ_inact(i)) = i
enddo
END_PROVIDER

5
src/cipsi/selection.irp.f
View File

@ -627,6 +627,11 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
sum_e_pert = sum_e_pert + e_pert * selection_weight(istate)
! endif
end do
if(pseudo_sym)then
if(dabs(mat(1, p1, p2)).lt.thresh_sym)then
sum_e_pert = 10.d0
endif
endif
if(sum_e_pert <= buf%mini) then
call add_to_selection_buffer(buf, det, sum_e_pert)

6
src/cis/cis.irp.f
View File

@ -57,7 +57,11 @@ subroutine run
implicit none
integer :: i
call H_apply_cis
if(pseudo_sym)then
call H_apply_cis_sym
else
call H_apply_cis
endif
print *, 'N_det = ', N_det
print*,'******************************'
print *, 'Energies of the states:'

5
src/cis/h_apply.irp.f
View File

@ -5,5 +5,10 @@ BEGIN_SHELL [ /usr/bin/env python2 ]
from generate_h_apply import H_apply
H = H_apply("cis",do_double_exc=False)
print H
H = H_apply("cis_sym",do_double_exc=False)
H.filter_only_connected_to_hf()
print H
END_SHELL

6
src/cisd/cisd.irp.f
View File

@ -53,7 +53,11 @@ subroutine run
implicit none
integer :: i
call H_apply_cisd
if(pseudo_sym)then
call H_apply_cisd_sym
else
call H_apply_cisd
endif
print *, 'N_det = ', N_det
print*,'******************************'
print *, 'Energies of the states:'

4
src/cisd/h_apply.irp.f
View File

@ -5,5 +5,9 @@ BEGIN_SHELL [ /usr/bin/env python2 ]
from generate_h_apply import H_apply
H = H_apply("cisd",do_double_exc=True)
print H
H = H_apply("cisd_sym",do_double_exc=True)
H.filter_only_connected_to_hf()
print H
END_SHELL

30
src/density_for_dft/density_for_dft.irp.f
View File

@ -9,6 +9,8 @@ BEGIN_PROVIDER [double precision, one_e_dm_mo_alpha_for_dft, (mo_num,mo_num, N_s
one_e_dm_mo_alpha_for_dft = data_one_e_dm_alpha_mo + damping_for_rs_dft * delta_alpha
else if (density_for_dft .EQ. "input_density")then
one_e_dm_mo_alpha_for_dft = data_one_e_dm_alpha_mo
else if (density_for_dft .EQ. "input_density_ao")then
call ao_to_mo(data_one_e_dm_alpha_mo,size(data_one_e_dm_alpha_mo,1),one_e_dm_mo_alpha_for_dft,size(one_e_dm_mo_alpha_for_dft,1))
else if (density_for_dft .EQ. "WFT")then
provide mo_coef
one_e_dm_mo_alpha_for_dft = one_e_dm_mo_alpha
@ -58,6 +60,8 @@ BEGIN_PROVIDER [double precision, one_e_dm_mo_beta_for_dft, (mo_num,mo_num, N_st
one_e_dm_mo_beta_for_dft = data_one_e_dm_beta_mo + damping_for_rs_dft * delta_beta
else if (density_for_dft .EQ. "input_density")then
one_e_dm_mo_beta_for_dft = data_one_e_dm_beta_mo
else if (density_for_dft .EQ. "input_density_ao")then
call ao_to_mo(data_one_e_dm_beta_mo,size(data_one_e_dm_beta_mo,1),one_e_dm_mo_beta_for_dft,size(one_e_dm_mo_beta_for_dft,1))
else if (density_for_dft .EQ. "WFT")then
provide mo_coef
one_e_dm_mo_beta_for_dft = one_e_dm_mo_beta
@ -119,16 +123,22 @@ END_PROVIDER
one_e_dm_alpha_ao_for_dft = 0.d0
one_e_dm_beta_ao_for_dft = 0.d0
do istate = 1, N_states
call mo_to_ao_no_overlap( one_e_dm_mo_alpha_for_dft(1,1,istate), &
size(one_e_dm_mo_alpha_for_dft,1), &
one_e_dm_alpha_ao_for_dft(1,1,istate), &
size(one_e_dm_alpha_ao_for_dft,1) )
call mo_to_ao_no_overlap( one_e_dm_mo_beta_for_dft(1,1,istate), &
size(one_e_dm_mo_beta_for_dft,1), &
one_e_dm_beta_ao_for_dft(1,1,istate), &
size(one_e_dm_beta_ao_for_dft,1) )
enddo
if (density_for_dft .EQ. "input_density_ao")then
one_e_dm_alpha_ao_for_dft = data_one_e_dm_alpha_ao
one_e_dm_beta_ao_for_dft = data_one_e_dm_beta_ao
else
do istate = 1, N_states
call mo_to_ao_no_overlap( one_e_dm_mo_alpha_for_dft(1,1,istate), &
size(one_e_dm_mo_alpha_for_dft,1), &
one_e_dm_alpha_ao_for_dft(1,1,istate), &
size(one_e_dm_alpha_ao_for_dft,1) )
call mo_to_ao_no_overlap( one_e_dm_mo_beta_for_dft(1,1,istate), &
size(one_e_dm_mo_beta_for_dft,1), &
one_e_dm_beta_ao_for_dft(1,1,istate), &
size(one_e_dm_beta_ao_for_dft,1) )
enddo
endif
END_PROVIDER

13
src/determinants/EZFIO.cfg
View File

@ -89,3 +89,16 @@ doc: Weight of the states in state-average calculations.
interface: ezfio
size: (determinants.n_states)
[thresh_sym]
type: Threshold
doc: Thresholds to check if a determinant is connected with HF
interface: ezfio,provider,ocaml
default: 1.e-15
[pseudo_sym]
type: logical
doc: If |true|, discard any Slater determinants with an interaction smaller than thresh_sym with HF.
interface: ezfio,provider,ocaml
default: False

5
src/determinants/h_apply.template.f
View File

@ -150,6 +150,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
logical :: is_a_2h1p
logical :: is_a_2h
logical :: b_cycle
logical :: yes_no
check_double_excitation = .True.
iproc = iproc_in
@ -284,6 +285,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
$only_1h_double
$only_1p_double
$only_2h1p_double
$filter_only_connected_to_hf_double
key_idx += 1
do k=1,N_int
keys_out(k,1,key_idx) = key(k,1)
@ -339,6 +341,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
$only_1h_double
$only_1p_double
$only_2h1p_double
$filter_only_connected_to_hf_double
key_idx += 1
do k=1,N_int
keys_out(k,1,key_idx) = key(k,1)
@ -412,6 +415,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,fock_diag_tmp,i_generato
logical :: is_a_1h
logical :: is_a_1p
logical :: is_a_2p
logical :: yes_no
do k=1,N_int
key_mask(k,1) = 0_bit_kind
@ -493,6 +497,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,fock_diag_tmp,i_generato
$filter_only_1h1p_single
$filter_only_1h2p_single
$filter_only_2h2p_single
$filter_only_connected_to_hf_single
key_idx += 1
do k=1,N_int
keys_out(k,1,key_idx) = hole(k,1)

35
src/determinants/slater_rules.irp.f
View File

@ -2257,3 +2257,38 @@ subroutine i_H_j_double_alpha_beta(key_i,key_j,Nint,hij)
end
subroutine connected_to_hf(key_i,yes_no)
implicit none
use bitmasks
integer(bit_kind), intent(in) :: key_i(N_int,2)
logical , intent(out) :: yes_no
double precision :: hij,phase
integer :: exc(0:2,2,2)
integer :: degree
integer :: m,p
yes_no = .True.
call get_excitation_degree(ref_bitmask,key_i,degree,N_int)
if(degree == 2)then
call i_H_j(ref_bitmask,key_i,N_int,hij)
if(dabs(hij) .lt. thresh_sym)then
yes_no = .False.
endif
else if(degree == 1)then
call get_single_excitation(ref_bitmask,key_i,exc,phase,N_int)
! Single alpha
if (exc(0,1,1) == 1) then
m = exc(1,1,1)
p = exc(1,2,1)
! Single beta
else
m = exc(1,1,2)
p = exc(1,2,2)
endif
hij = mo_one_e_integrals(m,p)
if(dabs(hij) .lt. thresh_sym)then
yes_no = .False.
endif
else if(degree == 0)then
yes_no = .True.
endif
end

18
src/tools/rotate_mos.irp.f Normal file
View File

@ -0,0 +1,18 @@
program rotate_mos
implicit none
integer :: iorb,jorb
read(5,*)iorb,jorb
double precision, allocatable :: mo_coef_tmp(:,:)
allocate(mo_coef_tmp(ao_num,mo_num))
mo_coef_tmp = mo_coef
integer :: i,j
double precision :: dsqrt2_inv
dsqrt2_inv = 1.d0/dsqrt(2.d0)
do i = 1, ao_num
mo_coef(i,iorb) = dsqrt2_inv * ( mo_coef_tmp(i,iorb) + mo_coef_tmp(i,jorb) )
mo_coef(i,jorb) = dsqrt2_inv * ( mo_coef_tmp(i,iorb) - mo_coef_tmp(i,jorb) )
enddo
touch mo_coef
call save_mos
end

10
src/tools/save_one_e_dm.irp.f
View File

@ -1,15 +1,15 @@
program save_one_e_dm
implicit none
BEGIN_DOC
! Program that computes the one body density on the |MO| basis
! Program that computes the one body density on the |MO| and |AO| basis
! for $\alpha$ and $\beta$ electrons from the wave function
! stored in the |EZFIO| directory, and then saves it into the
! :ref:`module_aux_quantities`.
!
! Then, the global variable :option:`aux_quantities data_one_e_dm_alpha_mo`
! and :option:`aux_quantities data_one_e_dm_beta_mo` will automatically
! read this density in the next calculation. This can be used to perform
! damping on the density in |RSDFT| calculations (see
! and :option:`aux_quantities data_one_e_dm_beta_mo` (and the corresponding for |AO|)
! will automatically ! read this density in the next calculation.
! This can be used to perform damping on the density in |RSDFT| calculations (see
! :ref:`module_density_for_dft`).
END_DOC
read_wf = .True.
@ -25,4 +25,6 @@ subroutine routine_save_one_e_dm
END_DOC
call ezfio_set_aux_quantities_data_one_e_dm_alpha_mo(one_e_dm_mo_alpha)
call ezfio_set_aux_quantities_data_one_e_dm_beta_mo(one_e_dm_mo_beta)
call ezfio_set_aux_quantities_data_one_e_dm_alpha_ao(one_e_dm_ao_alpha)
call ezfio_set_aux_quantities_data_one_e_dm_beta_ao(one_e_dm_ao_beta)
end