On Tue, May 26, 2009 at 12:07 PM, Lisandro Dalcin <span dir="ltr"><<a href="mailto:dalcinl@gmail.com">dalcinl@gmail.com</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
The only real issue I see here is that the term "diagonal block" could<br>
lead to misinterpretation. But the behaviour you commented is the one<br>
I would expect.</blockquote><div><br>Okay, anyone else not want me to try and do this?<br> </div><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
<br>
The rest of this mail is perhaps off-topic, but it is somewhat related.<br>
<br>
The other issues is that some Mat preallocation routines are lazy<br>
about properly setting the cmap. I had to implement some hackery in<br>
petsc4py to make sure that created Mat's (at least these created from<br>
petsc4py calls) have a cmap consistent with the local/global<br>
row/column sizes.<br>
<br>
Moreover, supose you create a global square matrix with the following<br>
row/col distribution (I know, it is perhaps a nonsense, but IMHO is a<br>
valid distribution)<br>
from petsc4py import PETSc<br>
<br>
size = PETSc.COMM_WORLD.getSize()<br>
assert size == 2<br>
rank = PETSc.COMM_WORLD.getRank()<br>
Ms = [2, 6]<br>
Ns = [6, 2]<br>
<br>
m = Ms[rank]<br>
n = Ns[rank]<br>
<br>
A = PETSc.Mat().createAIJ([(m, None), (n, None)]) # None is<br>
intrepreted as PETSc.DECIDE<br>
for i in range(8):<br>
    A[i,i] = 1.0/(i+1)<br>
A.assemble()<br>
<br>
x, b = A.getVecs()<br>
x.set(1)<br>
A.mult(x,b)<br>
<br>
Up to now, all works as expected, all rmap/cmap are consistent, the<br>
'x' and 'b' vecs have the same global size, though different local<br>
sizes, and the MatMult succeed. All fine... But now, suppose you do:<br>
<br>
ksp = PETSc.KSP().create()<br>
ksp.setOperators(A)<br>
ksp.pc.type = 'none'<br>
ksp.solve(b,x)<br>
<br>
You are lost, this does not work (even if you try with nonzero initial<br>
guess)... But it should, right?<br><div><div class="h5"></div></div></blockquote><div><br>No, I  do not think that this should work, because the Krylov methods impose<br>a more stringent condition on the matrix A, namely that the layout of Ax should<br>
be the same as the layout of x. This necessitates a symmetric partition, and is<br>probably at the root of our problems with rectangular matrices.<br><br>  Matt<br> </div><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
<div><div class="h5">
On Tue, May 26, 2009 at 10:35 AM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br>
> I think there is an inconsistency in our definitions. If you create a matrix<br>
> which<br>
> is 4x8 on 2 procs, then you end up with a row and column map like this:<br>
><br>
>     rowMap           colMap<br>
>   [0..2) [2..4)      [0..4) [4..8)<br>
><br>
> However, we use the column map for queries which then give nonsense. The<br>
> most<br>
> damaging is in MatSetValues(). If I give a value for (2,2) on proc 1, then<br>
> since<br>
> 2 < cstart = 4, this value is considered "off-diagonal" and can cause a<br>
> malloc(). So<br>
> I think we have a choice. We can define "diagonal block" to have this weird<br>
> shape,<br>
> or we can change all queries for the diagonal block to use only the rowMap.<br>
> I vote<br>
> for the latter, but does anyone have an objection?<br>
><br>
>    Matt<br>
><br>
> --<br>
> What most experimenters take for granted before they begin their experiments<br>
> is infinitely more interesting than any results to which their experiments<br>
> lead.<br>
> -- Norbert Wiener<br>
><br>
<br>
<br>
<br>
</div></div><font color="#888888">--<br>
Lisandro Dalcín<br>
---------------<br>
Centro Internacional de Métodos Computacionales en Ingeniería (CIMEC)<br>
Instituto de Desarrollo Tecnológico para la Industria Química (INTEC)<br>
Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET)<br>
PTLC - Güemes 3450, (3000) Santa Fe, Argentina<br>
Tel/Fax: +54-(0)342-451.1594<br>
</font></blockquote></div><br><br clear="all"><br>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener<br>