<div dir="ltr"><div>Hello all and thanks for your great work in bringing this very helpful package to the community !</div><div><br></div><div>That said, I wouldn't need this mailing list if everything was running smoothly. I have a rather involved eigenvalue problem that I've been working on that's been throwing a mysterious error :<br></div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">petsc4py.PETSc.Error: error code 77<br></blockquote><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">[1] EPSSolve() at /usr/local/slepc/src/eps/interface/epssolve.c:149<br>[1] EPSSolve_KrylovSchur_Default() at /usr/local/slepc/src/eps/impls/krylov/krylovschur/krylovschur.c:289<br>[1] EPSGetStartVector() at /usr/local/slepc/src/eps/interface/epssolve.c:824<br>[1] Petsc has generated inconsistent data<br>[1] Initial vector is zero or belongs to the deflation space</blockquote></div><div> </div><div>This problem occurs in parallel with two processors, using the <a href="https://www.mcs.anl.gov/petsc/petsc4py-current/docs/apiref/petsc4py-module.html">petsc4py</a> library using the dolfinx/dolfinx <a href="https://hub.docker.com/r/dolfinx/dolfinx">docker container</a>. I have PETSc version 3.16.0, in complex mode, python 3, and I'm running all of that on a OpenSUSE Leap 15.2 machine (but I think the docker container has a Ubuntu OS).</div><div><br></div><div>I wrote a minimal working example below, but I'm afraid the process for building the matrices is involved, so I decided to directly share the matrices instead : <a href="https://seminaris.polytechnique.fr/share/s/ryJ6L2nR4ketDwP">https://seminaris.polytechnique.fr/share/s/ryJ6L2nR4ketDwP</a></div><div><br></div><div>They are in binary format, but inside the container I hope someone could reproduce my issue. A word on the structure and intent behind these matrices :</div><div><ul><li>QE is a diagonal rectangular real matrix. Think of it as some sort of preconditioner</li><li>L is the least dense of them all, the only one that is complex, and in order to avoid inverting it I'm using two KSPs to compute solve problems on the fly</li><li>Mf is a diagonal square real matrix, its on the right-hand side of the Generalised Hermitian Eigenvalue problem (I'm solving QE^H*L^-1H*L^-1*QE*x=lambda*Mf*x<br></li></ul><div>Full MWE is below :<br></div></div><br><div style="color:rgb(212,212,212);background-color:rgb(30,30,30);font-family:"Droid Sans Mono","monospace",monospace,"Droid Sans Fallback";font-weight:normal;font-size:14px;line-height:19px;white-space:pre"><div><span style="color:rgb(197,134,192)">from</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">petsc4py</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">import</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">PETSc</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">as</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">pet</span></div><div><span style="color:rgb(197,134,192)">from</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">slepc4py</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">import</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">SLEPc</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">as</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">slp</span></div><div><span style="color:rgb(197,134,192)">from</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">mpi4py</span><span style="color:rgb(212,212,212)">.</span><span style="color:rgb(78,201,176)">MPI</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">import</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(79,193,255)">COMM_WORLD</span></div><br><div><span style="color:rgb(106,153,85)"># Global sizes</span></div><div><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">.</span><span style="color:rgb(156,220,254)">rank</span><span style="color:rgb(212,212,212)">*(</span><span style="color:rgb(181,206,168)">490363</span><span style="color:rgb(212,212,212)">-</span><span style="color:rgb(181,206,168)">489780</span><span style="color:rgb(212,212,212)">)+</span><span style="color:rgb(181,206,168)">489780</span></div><div><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">.</span><span style="color:rgb(156,220,254)">rank</span><span style="color:rgb(212,212,212)">*(</span><span style="color:rgb(181,206,168)">452259</span><span style="color:rgb(212,212,212)">-</span><span style="color:rgb(181,206,168)">451743</span><span style="color:rgb(212,212,212)">)+</span><span style="color:rgb(181,206,168)">451743</span></div><div><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(181,206,168)">980143</span></div><div><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(181,206,168)">904002</span></div><br><div><span style="color:rgb(79,193,255)">QE</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Mat().createAIJ([[</span><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">],[</span><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">]])</span></div><div><span style="color:rgb(79,193,255)">L</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Mat().createAIJ([[</span><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">],[</span><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">]])</span></div><div><span style="color:rgb(156,220,254)">Mf</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Mat().createAIJ([[</span><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">],[</span><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">]])</span></div><br><div><span style="color:rgb(156,220,254)">viewerQE</span><span style="color:rgb(212,212,212)"> = </span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Viewer().createMPIIO(</span><span style="color:rgb(206,145,120)">"QE.dat"</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(206,145,120)">'r'</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(79,193,255)">QE</span><span style="color:rgb(212,212,212)">.load(</span><span style="color:rgb(156,220,254)">viewerQE</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(156,220,254)">viewerL</span><span style="color:rgb(212,212,212)"> = </span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Viewer().createMPIIO(</span><span style="color:rgb(206,145,120)">"L.dat"</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(206,145,120)">'r'</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(79,193,255)">L</span><span style="color:rgb(212,212,212)">.load(</span><span style="color:rgb(156,220,254)">viewerL</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(156,220,254)">viewerMf</span><span style="color:rgb(212,212,212)"> = </span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Viewer().createMPIIO(</span><span style="color:rgb(206,145,120)">"Mf.dat"</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(206,145,120)">'r'</span><span style="color:rgb(212,212,212)">, </span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(156,220,254)">Mf</span><span style="color:rgb(212,212,212)">.load(</span><span style="color:rgb(156,220,254)">viewerMf</span><span style="color:rgb(212,212,212)">)</span></div><br><div><span style="color:rgb(79,193,255)">QE</span><span style="color:rgb(212,212,212)">.assemble()</span></div><div><span style="color:rgb(79,193,255)">L</span><span style="color:rgb(212,212,212)">.assemble()</span></div><br><div><span style="color:rgb(156,220,254)">KSPs</span><span style="color:rgb(212,212,212)"> = []</span></div><div><span style="color:rgb(106,153,85)"># Useful solvers (here to put options for computing a smart R)</span></div><div><span style="color:rgb(197,134,192)">for</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(156,220,254)">Mat</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(197,134,192)">in</span><span style="color:rgb(212,212,212)"> [</span><span style="color:rgb(79,193,255)">L</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(79,193,255)">L</span><span style="color:rgb(212,212,212)">.hermitianTranspose()]:</span></div><div><span style="color:rgb(212,212,212)">    </span><span style="color:rgb(79,193,255)">KSP</span><span style="color:rgb(212,212,212)"> = </span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.KSP().create()</span></div><div><span style="color:rgb(212,212,212)">    </span><span style="color:rgb(79,193,255)">KSP</span><span style="color:rgb(212,212,212)">.setOperators(</span><span style="color:rgb(156,220,254)">Mat</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">    </span><span style="color:rgb(79,193,255)">KSP</span><span style="color:rgb(212,212,212)">.setFromOptions()</span></div><div><span style="color:rgb(212,212,212)">    </span><span style="color:rgb(156,220,254)">KSPs</span><span style="color:rgb(212,212,212)">.</span><span style="color:rgb(220,220,170)">append</span><span style="color:rgb(212,212,212)">(</span><span style="color:rgb(79,193,255)">KSP</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">    </span></div><div><span style="color:rgb(86,156,214)">class</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(78,201,176)">LHS_class</span><span style="color:rgb(212,212,212)">:            </span></div><div><span style="color:rgb(212,212,212)">    </span><span style="color:rgb(86,156,214)">def</span><span style="color:rgb(212,212,212)"> </span><span style="color:rgb(220,220,170)">mult</span><span style="color:rgb(212,212,212)">(</span><span style="color:rgb(156,220,254)">self</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">A</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">x</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">y</span><span style="color:rgb(212,212,212)">):</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(156,220,254)">w</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Vec().createMPI([</span><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">],</span><span style="color:rgb(156,220,254)">comm</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(156,220,254)">z</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Vec().createMPI([</span><span style="color:rgb(156,220,254)">m_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">m</span><span style="color:rgb(212,212,212)">],</span><span style="color:rgb(156,220,254)">comm</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(79,193,255)">QE</span><span style="color:rgb(212,212,212)">.mult(</span><span style="color:rgb(156,220,254)">x</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">w</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(156,220,254)">KSPs</span><span style="color:rgb(212,212,212)">[</span><span style="color:rgb(181,206,168)">0</span><span style="color:rgb(212,212,212)">].solve(</span><span style="color:rgb(156,220,254)">w</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">z</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(156,220,254)">KSPs</span><span style="color:rgb(212,212,212)">[</span><span style="color:rgb(181,206,168)">1</span><span style="color:rgb(212,212,212)">].solve(</span><span style="color:rgb(156,220,254)">z</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">w</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(212,212,212)">        </span><span style="color:rgb(79,193,255)">QE</span><span style="color:rgb(212,212,212)">.multTranspose(</span><span style="color:rgb(156,220,254)">w</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">y</span><span style="color:rgb(212,212,212)">)</span></div><br><div><span style="color:rgb(106,153,85)"># Matrix free operator</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Mat()</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">.create(</span><span style="color:rgb(156,220,254)">comm</span><span style="color:rgb(212,212,212)">=</span><span style="color:rgb(79,193,255)">COMM_WORLD</span><span style="color:rgb(212,212,212)">)</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">.setSizes([[</span><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">],[</span><span style="color:rgb(156,220,254)">n_local</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">n</span><span style="color:rgb(212,212,212)">]])</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">.setType(</span><span style="color:rgb(78,201,176)">pet</span><span style="color:rgb(212,212,212)">.Mat.Type.PYTHON)</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">.setPythonContext(</span><span style="color:rgb(78,201,176)">LHS_class</span><span style="color:rgb(212,212,212)">())</span></div><div><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">.setUp()</span></div><br><div><span style="color:rgb(106,153,85)"># Eigensolver</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)"> = </span><span style="color:rgb(78,201,176)">slp</span><span style="color:rgb(212,212,212)">.EPS(); </span><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.create()</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.setOperators(</span><span style="color:rgb(79,193,255)">LHS</span><span style="color:rgb(212,212,212)">,</span><span style="color:rgb(156,220,254)">Mf</span><span style="color:rgb(212,212,212)">) </span><span style="color:rgb(106,153,85)"># Solve QE^T*L^-1H*L^-1*QE*x=lambda*Mf*x</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.setProblemType(</span><span style="color:rgb(78,201,176)">slp</span><span style="color:rgb(212,212,212)">.EPS.ProblemType.GHEP) </span><span style="color:rgb(106,153,85)"># Specify that A is hermitian (by construction), and B is semi-positive</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.setWhichEigenpairs(</span><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.Which.LARGEST_MAGNITUDE) </span><span style="color:rgb(106,153,85)"># Find largest eigenvalues</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.setFromOptions()</span></div><div><span style="color:rgb(79,193,255)">EPS</span><span style="color:rgb(212,212,212)">.solve()</span></div></div><br>Quentin CHEVALIER – IA parcours recherche<br><br>LadHyX - Ecole polytechnique<br><br>__________</div>