$p^{\dagger}q$

pdaggerq is a fermionic computer algebra package for bringing strings of creation / annihilation operators to normal order with respect to a true vacuum or the Fermi vacuum. The code can evaluate expressions like projections that arise in coupled cluster theory and can be used to generate working many-body electronic structure codes. In the examples section we provide worked examples that generate equations and Python code for CCSD, lambda-CCSD, CC3, CCSDT, CCSDTQ and more.

Installation

Installing pdaggerq requires cmake is installed on your system. To install, first clone the package

git clone git@github.com:edeprince3/pdaggerq.git

Then, you can install like you would a normal python package. From the package top level directory, run

python -m pip install .

which should compile pdaggerq. This command will produce a build folder that contains the compiled c++ shared library.

Quickstart

The following is an example that generates the energy expression for CCSD.

Python:

import pdaggerq

pq = pdaggerq.pq_helper("fermi")

print('')
print('    < 0 | e(-T) H e(T) | 0> :')
print('')

pq.add_st_operator(1.0,['f'],['t1','t2'])
pq.add_st_operator(1.0,['v'],['t1','t2'])

pq.simplify()
terms = pq.strings()
for term in terms:
    print(term)
    
pq.clear()

Output:

< 0 | e(-T) H e(T) | 0> :

['+1.00', 'f(i,i)']
['+1.00', 'f(i,a)', 't1(a,i)']
['-0.50', '<j,i||j,i>']
['+0.25', '<j,i||a,b>', 't2(a,b,j,i)']
['-0.50', '<j,i||a,b>', 't1(a,i)', 't1(b,j)']

The same result can be generated by explicitly setting commutators, double commutators, etc., that arise in the similarity transformation. For the sake of readability, commutators that evaluate to zero are excluded.

Python:

import pdaggerq

pq = pdaggerq.pq_helper("fermi")

print('')
print('    < 0 | e(-T) H e(T) | 0> :')
print('')

# one-electron part: 

# f
pq.add_operator_product(1.0,['f'])

# [f, T1]
pq.add_commutator(1.0,['f'],['t1'])

# [f, T2]
pq.add_commutator(1.0,['f'],['t2'])

# two-electron part: 

# v
pq.add_operator_product(1.0,['v'])

# [v, T1]
pq.add_commutator(1.0,['v'],['t1'])

# [v, T2]
pq.add_commutator(1.0,['v'],['t2(a,b,i,j)'])

# [[v, T1], T1]]
pq.add_double_commutator(0.5, ['v'],['t1'],['t1'])

pq.simplify()
terms = pq.strings()
for term in terms:
    print(term)
pq.clear()

PQ Graph

The module pq_graph is used to optimize the order of contractions within the expressions and to eliminate common subexpressions. The module can represent the many-body equations as a multidirected graph and can be output in either C++ or Python syntax for evaluation. Refer to the README.md in the pq_graph directory for more information.

How to contribute

We'd love to accept your contributions and patches. All submissions, including submissions by project members, require review. We use GitHub pull requests for this purpose. Consult GitHub Help for more information on using pull requests. Furthermore, please make sure your new code comes with extensive tests! Make sure you adhere to our style guide. Just have a look at our code for clues. We mostly follow PEP 8 and use the corresponding pylint to check. Code should always come with documentation.

We use Github issues for tracking requests and bugs.

How to cite

When using pdaggerq for research projects, please cite:

@misc{pdaggerq_2021,
    author       = {A. Eugene DePrince and Nicholas C. Rubin},
    title        = {{pdaggerq}: https://github.com/edeprince3/pdaggerq},
    month        = {June},
    year         = {2021},
    url          = {https://github.com/edeprince3/pdaggerq} 
}

Methods and Functionality

Normal Order

Normal order may be defined relative to the true vacuum or the fermi vacuum. This selection is made when creating the pq_helper class:

# true vacuum
pq = pdaggerq.pq_helper("")

or

# true vacuum
pq = pdaggerq.pq_helper("true")

or

# fermi vacuum
pq = pdaggerq.pq_helper("fermi")

Label convention

We follow the usual convention for labeling orbitals: i, j, k, l, m, and n represent occupied orbitals and a, b, c, d, e, and f represent virtual orbitals. Additionally, any label starting with i or a will be considered occupied or virtual, respectively (e.g., i_1 or a2). All other labels are considered general labels. When normal order is defined relative to the fermi vacuum, sums involving general labels are split into sums involving occupied and virtual orbitals using internal labels o1, o2, o3, and o4 (occupied) or v1, v2, v3, and v4 (virtual). So, we recommend avoiding using these labels when specifying any other components of your strings.

Orbital labels refer to spin orbitals. There is functionality to integrate out spin in final expressions, see below.

Operators

Currently supported fermionic operators include

a general one-body operator

'h' 

a general antisymmetrized two-body operator

'g' 

the fock operator

'f' 

the fluctuation potental operator

'v' 

fermionic creation / annihilation operators. e.g., $\hat{a}^{\dagger}_i$, $\hat{a}_j$, etc.

a*(i)
a(j)

up to four-body transition operators, e.g., p*q, p*q*rs, etc.

'e1(p,q)' 
'e2(p,q,r,s)' 

singles, doubles, triples, and quadruples t-amplitudes

't1'
't2' 
't3'
't4'

reference, singles, doubles, triples, and quadruples left-hand amplitudes

'l0'
'l1'  
'l2'  
'l3'
'l4'

reference, singles, doubles, triples, and quadruples right-hand amplitudes

'r0'
'r1'  
'r2'
'r3'
'r4'

Note that all factors such as the 1/4 associated with t2, l2, and r2 are handled internally.

Currently supported bosonic operators include

boson creation / annihilation operators, $\hat{b}^{\dagger}$ and $\hat{b}$

b+
b-

cluster operators that include both n-body electron excitations and m-body photon creation ($T_{n,m}$, with $n=0,1,2,3,4$ and $m=0,1,2$...)

't0,1'
't1,1'
etc.

left-hand operators that include both n-body electron de-excitations and m-body photon annihilation ($L_{n,m}$, with $n=0,1,2,3,4$ and $m=0,1,2$...)

'l2,2'
etc.

right-hand operators that include both n-body electron excitations and m-body photon creation ($R_{n,m}$, with $n=0,1,2,3,4$ and $m=0,1,2$...)

'r1,2'
etc.

a one-body fermionic operator times a boson creation operator, $\sum_{pq} d_{pq} \hat{a}^{\dagger}_p \hat{a}_q \hat{b}^{\dagger}$

'd+' 

a one-body fermionic operator times a boson annihilation operator, $\sum_{pq} d_{pq} \hat{a}^{\dagger}_p \hat{a}_q \hat{b}$

'd-' 

a diagonal boson operator, $\omega_0 \hat{b}^{\dagger}\hat{b}$

'w0'

Other operators

the unit operator

'1'

Methods

Strings are defined in Python using the pq_helper class, which has the following functions:

add_operator_product:

set strings corresponding to a product of operators.

add_operator_product( 0.5, ['h','t1','t2'])

add_commutator:

set strings corresponding to a commutator of operators. Note that the arguments are lists to allow for commutators of products of operators.

add_commutator(1.0, ['f'], ['t2'])

add_double_commutator:

set strings corresponding to a double commutator involving three operators. Note that the arguments are lists to allow for commutators of products of operators.

add_double_commutator(1.0/2.0, ['f'], ['t2'], ['t1'])

add_triple_commutator:

set strings corresponding to a triple commutator involving four operators. Note that the arguments are lists to allow for commutators of products of operators.

add_triple_commutator(1.0/6.0, ['f'], ['t2'], ['t1'], ['t1'])

add_quadruple_commutator:

set strings corresponding to a quadruple commutator involving five operators. Note that the arguments are lists to allow for commutators of products of operators.

add_quadruple_commutator(1.0/24.0, ['f'], ['t2'], ['t1'], ['t1'], ['t1'])

add_st_operator:

set strings corresponding to a similarity transformed operator commutator involving five operators. The first argument after the numerical value is a list of operators; the product of these operators will be similarity transformed. The next argument is a list of operators appearing as a sum the exponential function. The similarity transformation is performed by applying the BCH expansion with four nested commutators.

add_st_operator(1.0, ['v'],['t1','t2'])

set_print_level:

Control the amount of output. Any value greater than the default value of 0 will cause the code to print starting strings.

set_print_level(0)

set_left_operators:

set a sum (outer list) of products (inner lists) of operators that define the bra state

set_left_operators([['1'],['l1'],['l2']])

set_left_operators_type:

set the type of operator defined by 'l1', 'l2', etc., for different flavors of EOM-CC methods. Valid options are 'EE' (excitation energy), 'IP' (ionization potential), 'DIP' (double ionization potential, 'EA' (electron attachment), and 'DEA' (double electron attachment).

set_right_operators:

set a sum (outer list) of products (inner lists) of operators that define the ket state

set_right_operators([['r0'],['r1'],['r2']])

set_right_operators_type:

set the type of operator defined by 'r1', 'r2', etc., for different flavors of EOM-CC methods. Valid options are 'EE' (excitation energy), 'IP' (ionization potential), 'DIP' (double ionization potential, 'EA' (electron attachment), and 'DEA' (double electron attachment).

simplify:

consolidate/cancel terms and zero any delta functions that involve occupied / virtual combinations.

simplify()

strings:

returns the current list of strings. If normal order is defined with respect to the Fermi vacuum, only fully-contracted strings (those containing no creation / annihilation operators) will be returned. If normal order is defined with respect to the true vacuum, all strings are returned.

strings()

Strings could be spin integrated to eliminate non-spin-conserving terms by passing a dictionary of spins for non-summed labels to this function. If there are no non-summed labels, then an empty dictionary would suffice.

spin_labels = {
    'e' : 'a',
    'f' : 'b',
    'm' : 'a',
    'n' : 'b'
}
strings(spin_labels = spin_labels)

Active-space methods [e.g., CCSDt/CCSDtq, J. Chem. Phys. 110, 6103–6122 (1999)] could be realized by passing a dictionary of label ranges that specifies orbital spaces over which the amplitudes are defined.

label_ranges = {
    't3' : ['act', 'act', 'all', 'act', 'act', 'all'],
    't2' : ['all', 'all', 'all', 'all'],
    't1' : ['all', 'all'],
    'a' : ['act'],
    'i' : ['act']
}
strings(label_ranges = label_ranges)

Note that spin labels and label ranges cannot currently be specified simultaneously.

clear:

clear the current set of strings. Note that this function will not reset operator types specified using set_right_operators_type and set_left_operators_type.

clear()