AssertionError: refs[i] ≈ V' * Cs[i] * V

[jeremypparker]

symmetry_MRE.jl.txt

Here is a script which hopefully reproduces an error I am encountering with symmetries. If you're interested, it's an attempt to apply the methods of Oeri & Goluskin to a 3DoF Hamiltonian system.

I create a group which is generated by four matrices, using OSCAR. This creates a group of order 96. I have written a very simple wrapper around Oscar groups to support GroupsCore. I also use this group to generate bases for an unknown Lyapunov function V and an additional auxiliary function rho, and assemble these into a polynomial P. P is invariant under all of the symmetry actions. When I attempted to use SumOfSquares to constrain P to be nonneg, I get an error from ordered_block_diag

[jeremypparker]

After a lot of playing around I realised I can recreate this AssertionError directly by using the example here by replacing the line
GroupsCore.gens(G::DihedralGroup) = [DihedralElement(G.n, false, 1), DihedralElement(G.n, true, 0)]
with
GroupsCore.gens(G::DihedralGroup) = [DihedralElement(G.n, true, 1), DihedralElement(G.n, true, 0)]

As far as I can tell, these should both give valid generators for the group

[blegat]

Thanks, I could extract MRE's for ordered_block_diag that I added as tests in #384
It remains to corner down the issue ^^

[blegat]

I could fix the issue with dihedral groups but not with your example. In your example, we end up needed to find and orthogonal matrix U such that A = U' * B * U for all (A, B) in zip(As, Bs) where As and B are as defined in:

A1 = [1 0; 0 1]
A2 = [-1 0; 0 -1]
A3 = [1 0; 0 -1]
A4 = [0.5 1
      1 0.5]
A5 = [-0.5 0.866025
      0.866025 -0.5]
As = [A1, A2, A3, A4]
Bs = [A1, A2, -A3, A5]

By computing a schur decomposition of a linear combination of the matrices in As and Bs, we get different eigenvalues which shows that this is not possible. That makes me thing there might be an issue with the input, not with ordered_block_diag here bug I might be wrong.

[jeremypparker]

Currently travelling with a windows system so I can't test my Oscar example, but your fix does seem to solve the problem for my order 96 symmetry group if I assemble it manually without Oscar. So I don't know what to conclude without more playing around.

I agree that A4 and A5 as stated above are not similar, as they have different eigenvalues so there is no such matrix U. Are you suggesting that there is a problem with my test case, or that there is a different problem in SumOfSquares.jl before ordered_block_diag is called? The third possibility is that you made a mistake when extracting the unit test.

[blegat]

Not sure which one of the three is the case, in the meantime, I'll backport the fix to the current release because the master branch cannot be released right now

[jeremypparker]

I tried to cherry pick this fix onto 0.7.3 and apply it to the original code I wanted to work and now I hitting the error
"Could not simultaneously block-diagonalize into $m identical $dx$d blocks" (weddernburn.jl:190), except that this raises the error that dx is undefined! Perhaps I've missed some other relevant bug fixes since 0.7.3?

[blegat]

Can you reproduce it with this branch ?
https://github.com/jump-dev/SumOfSquares.jl/commits/release-0.7/
The tests are failing which prevented me to release. I need to backport some of the fixes of #372

[jeremypparker]

Same thing - "UndefVarError: dx not defined"

But the original script I posted at the top of this thread now seems to work, although I need to check whether I believe the output of the SDP. So clearly in producing my cleaned up version I've stopped hitting this new error. This also suggests that when you extracted the unit test from that script you made a mistake.

So this is good progress. I need to investigate more why I'm hitting this new error "Could not simultaneously block-diagonalize", but the original AssertionError seems to be fixed.

[jeremypparker]

Brief update on this: The original problem is definitely fixed so I will close this issue. However, the other problem I mention persists and annoyingly I can only reproduce it when I have very high-degree polynomials which needs a very large memory cluster. I will do my best to make a minimal example of the bug and open a new issue, I guess

[blegat]

Thanks for the update!