<div dir="auto">Thank you for that pointer ! I have on hand a partial SVD of R, so I used that to build the approximate matrix instead.<br>
<br>
It's great that so many nice features of PETSc like STSetPreconditionerMat are accessible through petsc4py !<br>
<br>
Good day,<br>
<br>
Quentin</div><br>
<br>
<br>
On Mon, 17 Jul 2023 at 17:29, Jose E. Roman <<a href="mailto:jroman@dsic.upv.es" target="_blank" rel="noreferrer">jroman@dsic.upv.es</a>> wrote:<br>
><br>
> It is possible to pass a different matrix to build the preconditioner. That is, the shell matrix for B (EPSSetOperators) and an explicit matrix (that approximates B) for the preconditioner. For instance, you can try passing M for building the preconditioner. Since M is an explicit matrix, you can try the default preconditioner (block Jacobi with ILU as local solver) or even a full LU decomposition. The effectiveness of the preconditioner will depend on how the update M+R^H P M P R moves the eigenvalues around.<br>
><br>
> You can do this with STSetSplitPreconditioner() or STSetPreconditionerMat(). In your case any of them will do.<br>
><br>
> Jose<br>
><br>
><br>
> > El 17 jul 2023, a las 15:50, Quentin Chevalier <<a href="mailto:quentin.chevalier@polytechnique.edu" target="_blank" rel="noreferrer">quentin.chevalier@polytechnique.edu</a>> escribió:<br>
> ><br>
> > Thank you for this suggestion, I tried to implement that but it's<br>
> > proven pretty hard to implement MATOP_GET_DIAGONAL without completely<br>
> > tanking performance. After all, B is a shell matrix for a reason : it<br>
> > looks like M+R^H P M P R with R itself a shell matrix.<br>
> ><br>
> > Allow me to point out that I have no shift. My eigenvalue problem is<br>
> > purely about the largest ones out there. Section 8.2 and 3.4.3 led me<br>
> > to think that there was a way to avoid computing (or writing a shell<br>
> > matrix about it) B^-1... But you seem to stress that there's no way<br>
> > around it.<br>
> ><br>
> > Quentin<br>
> ><br>
> ><br>
> ><br>
> > On Mon, 17 Jul 2023 at 11:56, Jose E. Roman <<a href="mailto:jroman@dsic.upv.es" target="_blank" rel="noreferrer">jroman@dsic.upv.es</a>> wrote:<br>
> >><br>
> >> The B-inner product is independent of the ST operator. See Table 3.2. In generalized eigenproblems you always have an inverse.<br>
> >><br>
> >> If your matrix is diagonally dominant, try implementing the MATOP_GET_DIAGONAL operation and using PCJACOBI. Apart from this, you have to build your own preconditioner.<br>
> >><br>
> >> Jose<br>
> >><br>
> >><br>
> >>> El 17 jul 2023, a las 11:48, Quentin Chevalier <<a href="mailto:quentin.chevalier@polytechnique.edu" target="_blank" rel="noreferrer">quentin.chevalier@polytechnique.edu</a>> escribió:<br>
> >>><br>
> >>> Hello Jose,<br>
> >>><br>
> >>> I guess I expected B to not be inverted but instead used as a mass for a problem-specific inner product since I specified GHEP as a problem type. p50 of the same user manual seems to imply that that would indeed be the case. I don't see what problem there would be with using a shell B matrix as a weighting matrix, as long as a mat utility is provided of course.<br>
> >>><br>
> >>> I tried the first approach - I set up my KSP as CG since B is hermitian positive-definite (I made a mistake in my first email), but I'm getting a KSPSolve has not converged, reason DIVERGED_ITS error. I'm letting it run for 1000 iterations already so it seems suspiciously slow for a CG solver.<br>
> >>><br>
> >>> I'm grappling with a shell preconditioner now to try and speed it up, but I'm unsure which one allows for shell matrices.<br>
> >>><br>
> >>> Thank you for your time,<br>
> >>><br>
> >>> Quentin<br>
> >>><br>
> >>><br>
> >>> On Wed, 12 Jul 2023 at 19:24, Jose E. Roman <<a href="mailto:jroman@dsic.upv.es" target="_blank" rel="noreferrer">jroman@dsic.upv.es</a>> wrote:<br>
> >>>><br>
> >>>> By default, it is solving the problem as B^{-1}*A*x=lambda*x (see chapter on Spectral Transformation). That is why A can be a shell matrix without problem. But B needs to be an explicit matrix in order to compute an LU factorization. If B is also a shell matrix then you should set an iterative solver for the associated KSP (see examples in the chapter).<br>
> >>>><br>
> >>>> An alternative is to create a shell matrix M that computes the action of B^{-1}*A, then pass M to the EPS solver as a standard eigenproblem.<br>
> >>>><br>
> >>>> Jose<br>
> >>>><br>
> >>>><br>
> >>>>> El 12 jul 2023, a las 19:04, Quentin Chevalier <<a href="mailto:quentin.chevalier@polytechnique.edu" target="_blank" rel="noreferrer">quentin.chevalier@polytechnique.edu</a>> escribió:<br>
> >>>>><br>
> >>>>> Hello PETSc Users,<br>
> >>>>><br>
> >>>>> I have a generalised eigenvalue problem : Ax= lambda Bx<br>
> >>>>> I used to have only A as a matrix-free method, I used mumps and an LU preconditioner, everything worked fine.<br>
> >>>>><br>
> >>>>> Now B is matrix-free as well, and my solver is returning an error : "MatSolverType mumps does not support matrix type python", which is ironic given it seem to handle A quite fine.<br>
> >>>>><br>
> >>>>> I have read in the user manual here that there some methods may require additional methods to be supplied for B like MATOP_GET_DIAGONAL but it's unclear to me exactly what I should be implementing and what is the best solver for my case.<br>
> >>>>><br>
> >>>>> A is hermitian, B is hermitian positive but not positive-definite or real. Therefore I have specified a GHEP problem type to the EPS object.<br>
> >>>>><br>
> >>>>> I use PETSc in complex mode through the petsc4py bridge.<br>
> >>>>><br>
> >>>>> Any help on how to get EPS to work for a generalised matrix-free case would be welcome. Performance is not a key issue here - I have a tractable high value case on hand.<br>
> >>>>><br>
> >>>>> Thank you for your time,<br>
> >>>>><br>
> >>>>> Quentin<br>
> >>>><br>
> >><br>
><br>