1
$\begingroup$

Apologies if this is the wrong place to ask this kind of question.

I have a simple question about Openfermion. I have two normal ordered FermionOperators A and B in Openfermion which are not equal. However, when I evaluate A == B, it returns True. Usually, the equality function seems to work fine, but in this specific case it doesn't seem to work. Below is the source code to replicate this problem.

My question is: why is this happening? I'm wondering if this is maybe a problem with my installation, or perhaps a bug (although, looking at the source code, I cannot see any problems with the equality function). Or maybe I have some sort of conceptual misunderstanding?

Source code:

from openfermion import FermionOperator, normal_ordered

#string representation of FermionOperators
A_string = '(4.525943261204361+0j) [0^ 0] +\n(-0.3849231244033131+0j) [0^ 1] +\n(0.8547260787716595+0j) [0^ 2] +\n(-1.1299596129151999+0j) [0^ 3] +\n(-0.3849231244033129+0j) [1^ 0] +\n(-4.381884800061626+0j) [1^ 0^ 1 0] +\n(-4.105905711065374+0j) [1^ 0^ 2 0] +\n(-0.944626497232035+0j) [1^ 0^ 2 1] +\n(3.5938875653898616+0j) [1^ 0^ 3 0] +\n(0.5559380952173161+0j) [1^ 0^ 3 1] +\n(-5.184266085346127+0j) [1^ 0^ 3 2] +\n(0.3677762971945129+0j) [1^ 1] +\n(3.6630156211728777+0j) [1^ 2] +\n(-2.5183774065038045+0j) [1^ 3] +\n(0.8547260787716596+0j) [2^ 0] +\n(-4.105905711065371+0j) [2^ 0^ 1 0] +\n(2.0232821070880056+0j) [2^ 0^ 2 0] +\n(0.8540908922866308+0j) [2^ 0^ 2 1] +\n(1.2460804358695565+0j) [2^ 0^ 3 0] +\n(-5.690396818975691+0j) [2^ 0^ 3 1] +\n(3.2577629759360036+0j) [2^ 0^ 3 2] +\n(3.663015621172878+0j) [2^ 1] +\n(-0.360964409584815+0j) [2^ 1^ 0^ 2 1 0] +\n(-0.515120856381232+0j) [2^ 1^ 0^ 3 1 0] +\n(1.2472039494390206+0j) [2^ 1^ 0^ 3 2 0] +\n(1.2551065744601149+0j) [2^ 1^ 0^ 3 2 1] +\n(-0.9446264972320322+0j) [2^ 1^ 1 0] +\n(0.8540908922866324+0j) [2^ 1^ 2 0] +\n(6.764363308219267+0j) [2^ 1^ 2 1] +\n(-0.5061307336295583+0j) [2^ 1^ 3 0] +\n(1.0665005168981174+0j) [2^ 1^ 3 1] +\n(-1.5870659190452852+0j) [2^ 1^ 3 2] +\n(-1.624576120995937+0j) [2^ 2] +\n(2.9075472932497357+0j) [2^ 3] +\n(-1.1299596129151999+0j) [3^ 0] +\n(3.5938875653898643+0j) [3^ 0^ 1 0] +\n(1.2460804358695592+0j) [3^ 0^ 2 0] +\n(-0.5061307336295612+0j) [3^ 0^ 2 1] +\n(4.19902020363552+0j) [3^ 0^ 3 0] +\n(2.44500249921506+0j) [3^ 0^ 3 1] +\n(-1.4836026643174591+0j) [3^ 0^ 3 2] +\n(-2.518377406503804+0j) [3^ 1] +\n(-0.5151208563812295+0j) [3^ 1^ 0^ 2 1 0] +\n(-0.3481302294514168+0j) [3^ 1^ 0^ 3 1 0] +\n(1.6862584012597264+0j) [3^ 1^ 0^ 3 2 0] +\n(3.6148854702404463+0j) [3^ 1^ 0^ 3 2 1] +\n(0.5559380952173167+0j) [3^ 1^ 1 0] +\n(-5.69039681897569+0j) [3^ 1^ 2 0] +\n(1.066500516898118+0j) [3^ 1^ 2 1] +\n(2.445002499215056+0j) [3^ 1^ 3 0] +\n(0.0213735881799324+0j) [3^ 1^ 3 1] +\n(6.043172803857193+0j) [3^ 1^ 3 2] +\n(2.9075472932497357+0j) [3^ 2] +\n(1.2472039494390192+0j) [3^ 2^ 0^ 2 1 0] +\n(1.686258401259725+0j) [3^ 2^ 0^ 3 1 0] +\n(-4.660442362479951+0j) [3^ 2^ 0^ 3 2 0] +\n(-6.470812728417476+0j) [3^ 2^ 0^ 3 2 1] +\n(-5.184266085346124+0j) [3^ 2^ 1 0] +\n(1.2551065744601169+0j) [3^ 2^ 1^ 2 1 0] +\n(3.6148854702404485+0j) [3^ 2^ 1^ 3 1 0] +\n(-6.470812728417465+0j) [3^ 2^ 1^ 3 2 0] +\n(-3.4393374365242755+0j) [3^ 2^ 1^ 3 2 1] +\n(3.2577629759360045+0j) [3^ 2^ 2 0] +\n(-1.5870659190452812+0j) [3^ 2^ 2 1] +\n(-1.4836026643174585+0j) [3^ 2^ 3 0] +\n(6.043172803857192+0j) [3^ 2^ 3 1] +\n(-0.7285308515736753+0j) [3^ 2^ 3 2] +\n(0.996525236164558+0j) [3^ 3]'
B_string = '(4.525943261204361+0j) [0^ 0] +\n(-0.3849231244033139+0j) [0^ 1] +\n(0.8547260787716593+0j) [0^ 2] +\n(-1.1299596129151999+0j) [0^ 3] +\n(-0.3849231244033138+0j) [1^ 0] +\n(-4.381884800061625+0j) [1^ 0^ 1 0] +\n(-4.105905711065372+0j) [1^ 0^ 2 0] +\n(-0.9446264972320364+0j) [1^ 0^ 2 1] +\n(3.5938875653898648+0j) [1^ 0^ 3 0] +\n(0.5559380952173187+0j) [1^ 0^ 3 1] +\n(-5.184266085346124+0j) [1^ 0^ 3 2] +\n(0.36777629719451244+0j) [1^ 1] +\n(3.6630156211728773+0j) [1^ 2] +\n(-2.518377406503804+0j) [1^ 3] +\n(0.8547260787716591+0j) [2^ 0] +\n(-4.105905711065371+0j) [2^ 0^ 1 0] +\n(2.023282107088007+0j) [2^ 0^ 2 0] +\n(0.8540908922866277+0j) [2^ 0^ 2 1] +\n(1.2460804358695592+0j) [2^ 0^ 3 0] +\n(-5.690396818975689+0j) [2^ 0^ 3 1] +\n(3.257762975936007+0j) [2^ 0^ 3 2] +\n(3.6630156211728773+0j) [2^ 1] +\n(-0.36096440958481857+0j) [2^ 1^ 0^ 2 1 0] +\n(-0.5151208563812304+0j) [2^ 1^ 0^ 3 1 0] +\n(1.2472039494390195+0j) [2^ 1^ 0^ 3 2 0] +\n(1.2551065744601182+0j) [2^ 1^ 0^ 3 2 1] +\n(-0.944626497232038+0j) [2^ 1^ 1 0] +\n(0.8540908922866295+0j) [2^ 1^ 2 0] +\n(6.764363308219268+0j) [2^ 1^ 2 1] +\n(-0.5061307336295606+0j) [2^ 1^ 3 0] +\n(1.0665005168981163+0j) [2^ 1^ 3 1] +\n(-1.5870659190452843+0j) [2^ 1^ 3 2] +\n(-1.624576120995937+0j) [2^ 2] +\n(3.9075472932497366+0j) [2^ 3] +\n(-1.1299596129151999+0j) [3^ 0] +\n(3.5938875653898634+0j) [3^ 0^ 1 0] +\n(1.2460804358695587+0j) [3^ 0^ 2 0] +\n(-0.5061307336295595+0j) [3^ 0^ 2 1] +\n(4.199020203635514+0j) [3^ 0^ 3 0] +\n(2.445002499215058+0j) [3^ 0^ 3 1] +\n(-1.483602664317454+0j) [3^ 0^ 3 2] +\n(-2.518377406503804+0j) [3^ 1] +\n(-0.5151208563812286+0j) [3^ 1^ 0^ 2 1 0] +\n(-0.34813022945142125+0j) [3^ 1^ 0^ 3 1 0] +\n(1.6862584012597162+0j) [3^ 1^ 0^ 3 2 0] +\n(3.6148854702404476+0j) [3^ 1^ 0^ 3 2 1] +\n(0.5559380952173192+0j) [3^ 1^ 1 0] +\n(-5.690396818975689+0j) [3^ 1^ 2 0] +\n(1.066500516898117+0j) [3^ 1^ 2 1] +\n(2.4450024992150574+0j) [3^ 1^ 3 0] +\n(0.02137358817993462+0j) [3^ 1^ 3 1] +\n(6.043172803857192+0j) [3^ 1^ 3 2] +\n(2.9075472932497366+0j) [3^ 2] +\n(1.247203949439014+0j) [3^ 2^ 0^ 2 1 0] +\n(1.6862584012597197+0j) [3^ 2^ 0^ 3 1 0] +\n(-4.660442362479957+0j) [3^ 2^ 0^ 3 2 0] +\n(-6.470812728417478+0j) [3^ 2^ 0^ 3 2 1] +\n(-5.184266085346124+0j) [3^ 2^ 1 0] +\n(1.25510657446012+0j) [3^ 2^ 1^ 2 1 0] +\n(3.614885470240445+0j) [3^ 2^ 1^ 3 1 0] +\n(-6.470812728417472+0j) [3^ 2^ 1^ 3 2 0] +\n(-3.439337436524271+0j) [3^ 2^ 1^ 3 2 1] +\n(3.2577629759360067+0j) [3^ 2^ 2 0] +\n(-1.5870659190452825+0j) [3^ 2^ 2 1] +\n(-1.483602664317455+0j) [3^ 2^ 3 0] +\n(6.043172803857192+0j) [3^ 2^ 3 1] +\n(-0.7285308515736739+0j) [3^ 2^ 3 2] +\n(0.996525236164558+0j) [3^ 3]'

#create FermionOperators
A = normal_ordered(FermionOperator(A_string))
B = normal_ordered(FermionOperator(B_string))

#if A and B are equal, A - B should be zero up to numerical precision
print("A - B = {}".format(A - B)) 

#since A and B are not equal, A == B should return False
print("A == B: {}".format(A == B))

The output is

A - B = (-1.0000000000000009+0j) [2^ 3]
A == B: True
Improve this question
$\endgroup$

1 Answer 1

Reset to default
1
$\begingroup$

This looks like a bug$^1$. Look at following lines of code from the v1.3.0 release branch for the implementation of SymbolicOperator.isclose (which implements FermionOperator.__eq__):

if not (isinstance(a, sympy.Expr) or isinstance(b, sympy.Expr)):
    tol *= max(1, abs(a), abs(b))

This occurs inside of a for loop comparing A and B componentwise. tol is provided to isclose once at call time, then gets modified by this line to grow larger every time a term of A or B has a (non-sympy) coefficient greater than 1. If the coefficients of each term in A or B are all greater than 1, this will eventually corrupt the tolerance after some number of terms and the comparison becomes vacuous.

Here's a minimal reproducing example

x = FermionOperator("0^ 0")
y = FermionOperator("0^ 0")
print(EQ_TOLERANCE) # should be 1e-08 for 

# construct two identical operators up to some number of terms
num_terms_before_ineq = 30
for i in range(num_terms_before_ineq):
    x += FermionOperator(f" (10+0j) [0^ {i}]")
    y += FermionOperator(f" (10+0j) [0^ {i}]")
    
# add a final term that is equal within tol but gets missed by the isclose check
xfinal = FermionOperator(f" (1+0j) [0^ {num_terms_before_ineq + 1}]")
yfinal = FermionOperator(f" (2+0j) [0^ {num_terms_before_ineq + 1}]")
# these two terms are not equal within tol..
print(xfinal == yfinal) # >>> False
print(xfinal - yfinal) # (-1+0j) [0^ 31]

# ...but these two terms will be, because the `tol` argument to `isclose` was corrupted
x += xfinal
y += yfinal
print(x == y) # >>> True
print(x - y) # (-1+0j) [0^ 31]

If you make num_terms_before_ineq a bit smaller then this issue doesn't come up. But in your case there's a ton of terms, so this seems to be exactly what's happening.

$^\dagger$ Their intention here almost looks like comparing the terms based on a relative difference compared to $\lVert \mathbf{a} \rVert_2$, i.e. check if each term agrees within a relative tolerance $(a_i - b_i) / \text{max}(\lVert \mathbf{a} \rVert_2, \lVert \mathbf{b} \rVert_2) \leq \epsilon$. But in that case the docstring for equality checking is wrong,

Comparison is done for each term individually. Return True
if the difference between each term in self and other is
less than EQ_TOLERANCE
Improve this answer
$\endgroup$
2
  • $\begingroup$ in the spirit of open source i went ahead and filed an issue on their github $\endgroup$
    – forky40
    Commented Jan 23, 2022 at 20:36
  • 1
    $\begingroup$ Thanks a lot for your help and filing the issue on Github. My current work around is using A - B == 0 instead of A == B to check equality, which seems to avoid this issue a lot better. $\endgroup$
    – Solarflare0
    Commented Jan 23, 2022 at 21:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.