<html><head><meta http-equiv="Content-Type" content="text/html; charset=us-ascii"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><div class=""> Marcel,</div><div class=""><br class=""></div> Would you like to do the following? Compute <div class=""><br class=""></div><div class=""> Q A P where Q is a row permutation, P a column permutation and then apply LSQR on QAP? </div><div class=""> </div><div class=""><br class=""></div><div class=""> From the manual page: </div><div class=""><br class=""></div><div class="">In exact arithmetic the LSQR method (with no preconditioning) is identical to the KSPCG algorithm applied to the normal equations.</div><div class=""><br class=""></div><div class=""> [Q A P]' [Q A P] = P' A' A P = P'(A'A) P the Q drops out because permutation matrices' transposes are their inverse</div><div class=""><br class=""></div><div class=""> Note that P is a small square matrix. </div><div class=""><br class=""></div><div class=""> So my conclusion is that any column permutation of A is also a symmetric permutation of A'A so you can just try using regular reorderings of A'A if </div><div class="">you want to "concentrate" the "important" parts of A'A into your "block diagonal" preconditioner (and throw away the other parts)</div><div class=""><br class=""></div><div class=""> I don't know what it will do to the convergence. I've never had much luck generically trying to symmetrically reorder matrices to improve preconditioners but </div><div class="">for certain situation maybe it might help. For example if the matrix is [0 1; 1 0] and you permute it you get the [1 0; 0 1] which looks better.</div><div class=""><br class=""></div><div class=""> There is this <a href="https://epubs.siam.org/doi/10.1137/S1064827599361308" class="">https://epubs.siam.org/doi/10.1137/S1064827599361308</a> but it is for non-symmetric permutations and in your case if you use a non symmetric permeation you can no longer use LSQR. </div><div class=""><br class=""></div><div class=""> Barry</div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><div class=""><br class=""></div><div><br class=""><blockquote type="cite" class=""><div class="">On Oct 22, 2020, at 4:55 AM, Matthew Knepley <<a href="mailto:knepley@gmail.com" class="">knepley@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class=""><div dir="ltr" class="">On Thu, Oct 22, 2020 at 4:24 AM Marcel Huysegoms <<a href="mailto:m.huysegoms@fz-juelich.de" class="">m.huysegoms@fz-juelich.de</a>> wrote:<br class=""></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Hi all,<br class="">
<br class="">
I'm currently implementing a Gauss-Newton approach for minimizing a<br class="">
non-linear cost function using PETSc4py.<br class="">
The (rectangular) linear systems I am trying to solve have dimensions of<br class="">
about (5N, N), where N is in the range of several hundred millions.<br class="">
<br class="">
Due to its size and because it's an over-determined system, I use LSQR<br class="">
in conjunction with a preconditioner (which operates on A^T x A, e.g.<br class="">
BJacobi).<br class="">
Depending on the ordering of the unknowns the algorithm only converges<br class="">
for special cases. When I use a direct LR solver (as preconditioner) it<br class="">
consistently converges, but consumes too much memory. I have read in the<br class="">
manual that the LR solver internally also applies a matrix reordering<br class="">
beforehand.<br class="">
<br class="">
My question would be:<br class="">
How can I improve the ordering of the unknowns for a rectangular matrix<br class="">
(in order to converge also with iterative preconditioners)? If I use<br class="">
MatGetOrdering(), it only works for square matrices. Is there a way to<br class="">
achieve this from within PETSc4py?<br class="">
ParMETIS seems to be a promising framework for that task. Is it possible<br class="">
to apply its reordering algorithm to a rectangular PETSc-matrix?<br class="">
<br class="">
I would be thankful for every bit of advice that might help.<br class=""></blockquote><div class=""><br class=""></div><div class="">We do not have any rectangular reordering algorithms. I think your first step is to</div><div class="">find something in the literature that you think will work.</div><div class=""><br class=""></div><div class=""> Thanks,</div><div class=""><br class=""></div><div class=""> Matt</div><div class=""> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Best regards,<br class="">
Marcel<br class="">
<br class="">
<br class="">
------------------------------------------------------------------------------------------------<br class="">
------------------------------------------------------------------------------------------------<br class="">
Forschungszentrum Juelich GmbH<br class="">
52425 Juelich<br class="">
Sitz der Gesellschaft: Juelich<br class="">
Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498<br class="">
Vorsitzender des Aufsichtsrats: MinDir Volker Rieke<br class="">
Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender),<br class="">
Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt<br class="">
------------------------------------------------------------------------------------------------<br class="">
------------------------------------------------------------------------------------------------<br class="">
<br class="">
</blockquote></div><br clear="all" class=""><div class=""><br class=""></div>-- <br class=""><div dir="ltr" class="gmail_signature"><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div 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</div><div class=""><br class=""></div><div class=""><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class=""></div></div></div></div></div></div></div></div>
</div></blockquote></div><br class=""></div></body></html>