Sympy tensors for symbolic manipulation of quantum chemistry integrals

367 views
Skip to first unread message

Rick Muller

unread,
May 4, 2017, 10:51:42 AM5/4/17
to sympy
I was hoping that someone could give me some help getting started with the sympy tensor objects. I'd like to define symbolic objects to represent one- and two-electron integrals in quantum chemistry with the proper index permutation symmetries. These are real-valued integrals, so commutation relations aren't a problem (and, when they are, can be handled by the physics.secondquant module.

The one-electron integrals are symmetric, i.e. I1[i,j] = I1[j,i], which I assume should be straightforward.

The two-electron integrals are a little trickier, for I2[i,j,k,l] the integral is symmetric when i,j are permuted, and/or k,l are permuted, and/or i,j is permuted with k,l. I've never been able to derive a symbolic object that captures this, and it would be really convenient, for example, to derive equations for orbital optimization for different MC-SCF wave functions.

I'm familiar with techniques to compute the orbitals numerically, e.g., https://github.com/rpmuller/pyquante2. What I'm interested here is to derive and simplify equations for the symbolic manipulations of equations containing these terms. Has anyone done any work on this?

Thanks in advance,

Rick

Rick Muller

unread,
May 4, 2017, 12:55:38 PM5/4/17
to sympy
Here's what I did: I copied what Brian G did in the sympy.physics.secondquant module:

class TensorSymbol(sympy.Expr):
    is_commutative=True

class SymmetricTensor(TensorSymbol):
    def __new__(cls,symbol,inds1,inds2):
        symbol = sympy.sympify(symbol)
        inds1 = list(sympy.ordered(list(inds1)))
        inds2 = list(sympy.ordered(list(inds2)))
        inds1,inds2 = list(sympy.ordered([inds1,inds2]))
        return TensorSymbol.__new__(cls,symbol,inds1,inds2)

    def __str__(self):
        return "%s(%s,%s)" % self.args
        

>>> a,b,c,d = sympy.symbols('a b c d')

>>> SymmetricTensor('I2',(a,b),(c,d))
I2([a,b],[c,d])
>>> SymmetricTensor('I2',(c,d),(a,b))
I2([a,b],[c,d])

>>> SymmetricTensor('I1',(a,),(b,))
I1([a],[b])
>>> SymmetricTensor('I1',(b,),(a,))
I1([a],[b])


This is pretty ugly, but it works. Can anyone suggest improving it, for example using something more formal with the sympy.tensor module?

Ondřej Čertík

unread,
May 5, 2017, 12:31:59 AM5/5/17
to sympy
Hi Rick,
Can you give an example of an equation that contains it?

I.e., are you interested in symbolically calculating the two-particle
integrals of, say, Gaussians? Or rather just representing the
two-particle integral, and let SymPy take care of the freedom due to
the symmetries, as in the example you posted:

>>> SymmetricTensor('I2',(a,b),(c,d))
I2([a,b],[c,d])
>>> SymmetricTensor('I2',(c,d),(a,b))
I2([a,b],[c,d])

So you already have some proof-of-concept solution. How do you
actually want to use it? What is the application?

Ondrej

Rick Muller

unread,
May 5, 2017, 9:37:47 AM5/5/17
to sympy
Ondrej,

(Apologies for the length...)

Here's my larger plan. It's strictly a night-and-weekends thing, but I'd be interested in working with you (or anyone else here) who thinks this is fun.

I'd like to get better at generating expressions for correlated wave functions. Both to evaluate and optimize their energy using programs like PyQuante, Pyscf, Psi4, or whatever, and to generate quantum circuits (i.e. sequences of Clifford group operators in a quantum circuit) to run on quantum emulators like ProjectQ and others.

Hirata's tensor contraction engine is really impressive, and it's written in Python, although it doesn't make use of a lot of tools like Sympy that I think could substantially simplify it. I'm not as interested in generating optimized, parallel code (yet), mostly just with generating correct code right now.

There are two sympy paths to this: first, using the sympy.physics.secondquant operators to generate expressions. I'm getting up to speed on this as quickly as I can.

My post in this thread is for the second path, which is more of a first quantization approach, where we create operators that are terms of the energy expression. A lot of my thesis work took simple MC-SCF energy expressions:

E = sum_i f(i)I1(i,i) + sum_ij a(i,j)I2(i,i,j,j) + b(i,j)I2(i,j,i,j)

And creating expressions to optimize the corresponding set of equations, essentially using variational calculus (i.e. j <- j + dj). Sympy should be able to do this automatically, provided we can generate I1 and I2 operators that can simplify expressions in the right way.

I think the little toy code I posted will work for the above expressions. But it doesn't use any of the substantial power of the sympy tensor operators, which I'm only just now learning about. I'd be interested in seeing whether using something more powerful than my little sympy objects that simplify in the right way could open up new avenues, again, with the goal of generating something like TCE in python.

Sorry for the long post, trying to give you context. I think there's a real opportunity here for sympy, and I'd be willing to get together for an offline discussion with you or anyone else who would be interested in this.

Ondřej Čertík

unread,
May 5, 2017, 2:18:39 PM5/5/17
to sympy
Hi Rick,

On Fri, May 5, 2017 at 7:37 AM, Rick Muller <rpmu...@gmail.com> wrote:
> Ondrej,
>
> (Apologies for the length...)
>
> Here's my larger plan. It's strictly a night-and-weekends thing, but I'd be interested in working with you (or anyone else here) who thinks this is fun.
>
> I'd like to get better at generating expressions for correlated wave functions. Both to evaluate and optimize their energy using programs like PyQuante, Pyscf, Psi4, or whatever, and to generate quantum circuits (i.e. sequences of Clifford group operators in a quantum circuit) to run on quantum emulators like ProjectQ and others.

Is it like the Hylleraas’ method for He?

>
> Hirata's tensor contraction engine is really impressive, and it's written in Python, although it doesn't make use of a lot of tools like Sympy that I think could substantially simplify it. I'm not as interested in generating optimized, parallel code (yet), mostly just with generating correct code right now.
>
> There are two sympy paths to this: first, using the sympy.physics.secondquant operators to generate expressions. I'm getting up to speed on this as quickly as I can.
>
> My post in this thread is for the second path, which is more of a first quantization approach, where we create operators that are terms of the energy expression. A lot of my thesis work took simple MC-SCF energy expressions:
>
> E = sum_i f(i)I1(i,i) + sum_ij a(i,j)I2(i,i,j,j) + b(i,j)I2(i,j,i,j)
>
> And creating expressions to optimize the corresponding set of equations, essentially using variational calculus (i.e. j <- j + dj). Sympy should be able to do this automatically, provided we can generate I1 and I2 operators that can simplify expressions in the right way.

Would you mind sending me your thesis? I know you can write expression
for the energy E like that, but I thought that's the end, then you
just have to implement it. I've written code generation for MBPT,
where you have lots of such expressions too and one has to keep track
of them, but I generated the expressions directly from the graph, and
then I wrote code to generate all graphs.

>
> I think the little toy code I posted will work for the above expressions. But it doesn't use any of the substantial power of the sympy tensor operators, which I'm only just now learning about. I'd be interested in seeing whether using something more powerful than my little sympy objects that simplify in the right way could open up new avenues, again, with the goal of generating something like TCE in python.
>
> Sorry for the long post, trying to give you context. I think there's a real opportunity here for sympy, and I'd be willing to get together for an offline discussion with you or anyone else who would be interested in this.

I'd like to understand more what your idea is.

Ondrej

Rick Muller

unread,
May 6, 2017, 7:22:12 AM5/6/17
to sympy
Ondrej,

The best source is some course notes I wrote up while a staff member at Caltech. I just put them on github in http://github.com/rpmuller/PracticalQuantumChemistry. The PDF is at:

Still embarrassingly incomplete, but Sect 2.6-2.7, and Ch 9 have some of what you want. I can expand at length anytime. In fact, it's sometimes hard to stop me from talking about it ;-).

I'm talking about something along the lines of eq 2.33. Given an energy expression made of one- and two-electron terms that depend upon orbitals constructed of linear combinations of atomic basis functions, we can write a general expression for varying an orbital by some small amount \delta. These result in a set of orbital optimization equations that you can solve to find the optimal orbitals. For closed-shell HF, this results in the well-known Fock equations, and they're pretty simple, by which I mean they can be derived in a lot of different ways. General MC-SCF energy expressions like eq 9.12 are a little more complicated, but it's the same idea. However, rather than just diagonalizing a set of fock equations, you also have to do some occupied orbital rotations, and solve a CI equation.

I think this is well within the capabilities of sympy.

Apologies if this isn't clear. Again, this is one of my loves, so I'm happy to go into exhaustive detail on this. I'm only not doing so to spare you a lot of stuff, in case you're not interested.


On Thursday, May 4, 2017 at 8:51:42 AM UTC-6, Rick Muller wrote:
Reply all
Reply to author
Forward
0 new messages