<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Feb 12, 2021, at 2:32 AM, Florian Bruckner <<a href="mailto:e0425375@gmail.com" class="">e0425375@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class=""><div class="">Dear Jose, Dear Matt, <br class=""></div><div class=""><br class=""></div><div class="">I needed some time to think about your answers. <br class=""></div><div class="">If I understand correctly, the eigenmode solver internally uses A0^{-1}*B0, which is normally handled by the ST object, which creates a KSP solver and a corresponding preconditioner.</div><div class="">What I would need is an interface to provide not only the system Matrix A0 (which is an operator), but also a preconditioning matrix (sparse approximation of the operator).</div><div class="">Unfortunately this interface is not available, right?</div></div></div></blockquote><div><br class=""></div> If SLEPc does not provide this directly it is still intended to be trivial to provide the "preconditioner matrix" (that is matrix from which the preconditioner is built). Just get the KSP from the ST object and use KSPSetOperators() to provide the "preconditioner matrix" . </div><div><br class=""></div><div> Barry</div><div><br class=""></div><div><blockquote type="cite" class=""><div class=""><div dir="ltr" class=""><div class=""><br class=""></div><div class="">Matt directly creates A0^{-1}*B0 as a matshell operator. The operator uses a KSP with a proper PC internally. SLEPc would directly get A0^{-1}*B0 and solve a standard eigenvalue problem with this modified operator. Did I understand this correctly?</div><div class=""><br class=""></div><div class="">I have two further points, which I did not mention yet: the matrix B0 is Hermitian, but it is (purely) imaginary (B0.real=0). Right now, I am using Firedrake to set up the PETSc system matrices A0, i*B0 (which is real). Then I convert them into ScipyLinearOperators and use scipy.sparse.eigsh(B0, b=A0, Minv=Minv) to calculate the eigenvalues. Minv=A0^-1 is also solving within scipy using a preconditioned gmres. Advantage of this setup is that the imaginary B0 can be handled efficiently and also the post-processing of the eigenvectors (which requires complex arithmetics) is simplified. <br class=""></div><div class=""><br class=""></div><div class="">Nevertheless I think that the mixing of PETSc and Scipy looks too complicated and is not very flexible. <br class=""></div><div class="">If I would use Matt's approach, could I then simply switch between multiple standard eigenvalue methods (e.g. LOBPCG)? or is it limited due to the use of matshell?</div><div class="">Is there a solution for the imaginary B0, or do I have to use the non-hermitian methods? Is this a large performance drawback?</div><div class=""><br class=""></div><div class=""></div><div class="">thanks again,</div><div class="">and best wishes</div><div class="">Florian<br class=""></div></div><br class=""><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Feb 8, 2021 at 3:37 PM Jose E. Roman <<a href="mailto:jroman@dsic.upv.es" class="">jroman@dsic.upv.es</a>> wrote:<br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">The problem can be written as A0*v=omega*B0*v and you want the eigenvalues omega closest to zero. If the matrices were explicitly available, you would do shift-and-invert with target=0, that is<br class="">
<br class="">
(A0-sigma*B0)^{-1}*B0*v=theta*v for sigma=0, that is<br class="">
<br class="">
A0^{-1}*B0*v=theta*v<br class="">
<br class="">
and you compute EPS_LARGEST_MAGNITUDE eigenvalues theta=1/omega.<br class="">
<br class="">
Matt: I guess you should have EPS_LARGEST_MAGNITUDE instead of EPS_SMALLEST_REAL in your code. Are you getting the eigenvalues you need? EPS_SMALLEST_REAL will give slow convergence.<br class="">
<br class="">
Florian: I would not recommend setting the KSP matrices directly, it may produce strange side-effects. We should have an interface function to pass this matrix. Currently there is STPrecondSetMatForPC() but it has two problems: (1) it is intended for STPRECOND, so cannot be used with Krylov-Schur, and (2) it is not currently available in the python interface.<br class="">
<br class="">
The approach used by Matt is a workaround that does not use ST, so you can handle linear solves with a KSP of your own.<br class="">
<br class="">
As an alternative, since your problem is symmetric, you could try LOBPCG, assuming that the leftmost eigenvalues are those that you want (e.g. if all eigenvalues are non-negative). In that case you could use STPrecondSetMatForPC(), but the remaining issue is calling it from python.<br class="">
<br class="">
If you are using the git repo, I could add the relevant code.<br class="">
<br class="">
Jose<br class="">
<br class="">
<br class="">
<br class="">
> El 8 feb 2021, a las 14:22, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank" class="">knepley@gmail.com</a>> escribió:<br class="">
> <br class="">
> On Mon, Feb 8, 2021 at 7:04 AM Florian Bruckner <<a href="mailto:e0425375@gmail.com" target="_blank" class="">e0425375@gmail.com</a>> wrote:<br class="">
> Dear PETSc / SLEPc Users,<br class="">
> <br class="">
> my question is very similar to the one posted here: <br class="">
> <a href="https://lists.mcs.anl.gov/pipermail/petsc-users/2018-August/035878.html" rel="noreferrer" target="_blank" class="">https://lists.mcs.anl.gov/pipermail/petsc-users/2018-August/035878.html</a><br class="">
> <br class="">
> The eigensystem I would like to solve looks like:<br class="">
> B0 v = 1/omega A0 v<br class="">
> B0 and A0 are both hermitian, A0 is positive definite, but only given as a linear operator (matshell). I am looking for the largest eigenvalues (=smallest omega). <br class="">
> <br class="">
> I also have a sparse approximation P0 of the A0 operator, which i would like to use as precondtioner, using something like this:<br class="">
> <br class="">
> es = SLEPc.EPS().create(comm=fd.COMM_WORLD)<br class="">
> st = es.getST()<br class="">
> ksp = st.getKSP()<br class="">
> ksp.setOperators(self.A0, self.P0)<br class="">
> <br class="">
> Unfortunately PETSc still complains that it cannot create a preconditioner for a type 'python' matrix although P0.type == 'seqaij' (but A0.type == 'python'). <br class="">
> By the way, should P0 be an approximation of A0 or does it have to include B0?<br class="">
> <br class="">
> Right now I am using the krylov-schur method. Are there any alternatives if A0 is only given as an operator?<br class="">
> <br class="">
> Jose can correct me if I say something wrong.<br class="">
> <br class="">
> When I did this, I made a shell operator for the action of A0^{-1} B0 which has a KSPSolve() in it, so you can use your P0 preconditioning matrix, and<br class="">
> then handed that to EPS. You can see me do it here:<br class="">
> <br class="">
> <a href="https://gitlab.com/knepley/bamg/-/blob/master/src/coarse/bamgCoarseSpace.c#L123" rel="noreferrer" target="_blank" class="">https://gitlab.com/knepley/bamg/-/blob/master/src/coarse/bamgCoarseSpace.c#L123</a><br class="">
> <br class="">
> I had a hard time getting the embedded solver to work the way I wanted, but maybe that is the better way.<br class="">
> <br class="">
> Thanks,<br class="">
> <br class="">
> Matt<br class="">
> <br class="">
> thanks for any advice<br class="">
> best wishes<br class="">
> Florian<br class="">
> <br class="">
> <br class="">
> -- <br class="">
> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br class="">
> -- Norbert Wiener<br class="">
> <br class="">
> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class="">
<br class="">
</blockquote></div>
</div></blockquote></div><br class=""></body></html>