Extend Symmetrizer to support arbitrary integrals of motion about pomerol HOT 10 CLOSED

@aeantipov
So, what exactly has to be changed? Only Symmetrizer::compute()?
How much effort would it take?

from pomerol.

@krivenko
Thank you for your work, Igor.

Best,
Junya

from pomerol.

General integrals of motion make Pomerol crash.
Further investigation is required.

from pomerol.

Igor, can you upload a failing test?

from pomerol.

I'm now preparing a test for pomerol2triqs (5 orbital atom, no bath) and it segfaults somewhere in libpomerol.so. I admit, this could be a problem of pomerol2triqs, so I need a bit of time to rule this possibility out.

from pomerol.

Okay, here is what I've been able to find so far.

• The crash is triggered by use of a non-abelian integral of motion, such as S^2. Those guys are more tricky than N or S_z, because they partition Hilbert space so that c/c^+ are not necessarily one-to-one mappings between subspaces.

• Eigen reports a failed assertion in HamiltonianPart::prepare().

/usr/include/eigen3/Eigen/src/Core/DenseCoeffsBase.h:337: Eigen::DenseCoeffsBase<Derived, 1>::Scalar& Eigen::DenseCoeffsBase<Derived, 1>::operator()(Eigen::DenseCoeffsBase<Derived, 1>::Index, Eigen::DenseCoeffsBase<Derived, 1>::Index) [with Derived = Eigen::Matrix<std::complex<double>, -1, -1, 1>; Eigen::DenseCoeffsBase<Derived, 1>::Scalar = std::complex<double>; Eigen::DenseCoeffsBase<Derived, 1>::Index = long int]: Assertion `row >= 0 && row < rows() && col >= 0 && col < cols()' failed.

This is from line 45 of src/pomerol/HamiltonianPart.cpp, I guess.

• It's very likely that we have a problem in Operator::actRight() and/or StatesClassification::getInnerState().

from pomerol.

I can now see this being a conceptual problem, rather than a coding issue.

Here is the point: One cannot use non-abelian integrals of motion to partition the Hilbert space. In other words, an integral of motion Q can only be used if all initial basis states (Fock states) are its eigenstates.

Let's say H is our Hamiltonian, Q is an integral of motion we want to use for partitioning, and |\psi_n> is the orthonormal basis of Fock states. We characterize every basis state by an expectation value of Q computed in this state, q_n = <\psi_n|Q|\psi_n>.

By definition of integral of motion,

<\psi_1| [Q, H] |\psi_2> = <\psi_1| Q H |\psi_2> - <\psi_1| H Q |\psi_2> = 0.

In general we have Q|\psi_n> = q_n |\psi_n> + |\tilde\psi_n>, where |\tilde\psi_n> is orthogonal to |\psi_n>.

Let |\psi_1> and |\psi_2> be two different basis states with q_1 != q_2.

q_1 <\psi_1| H |\psi_2> - q_2 <\psi_1| H |\psi_2> +
<\tilde\psi_1| H |\psi_2> - <\psi_1| H |\tilde\psi_2> = 0

From this equation it follows that <\psi_1| H |\psi_2> = 0 (Hamiltonian does not connect subspaces characterized by different quantum numbers), only under an extra assumption |\tilde\psi_1> = |\tilde\psi_2> = 0. And the assumption is equivalent to saying that |\psi_1>, |\psi_2> are eigenvectors of Q.

Now the question is how we detect the inappropriate integrals of motion in Pomerol?

from pomerol.

from pomerol.

Basically what we need to do is to rotate the basis, such as to have the eigenstates of operators as basis states.

Rotation of the basis makes a little sense, as we would effectively need to diagonalize two operators, the Hamiltonian, and the integral of motion in question. And, in fact, even more if we have multiple integrals of motion.

I'd say checkSymmetry() should simply ignore IOM operators that don't commute with all occupation number operators.

from pomerol.

To be reopened if we have any idea how to do better.

from pomerol.