On Tue, May 26, 2009 at 1:48 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;">
<div class="im">On Tue, May 26, 2009 at 2:15 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br>
> On Tue, May 26, 2009 at 12:07 PM, Lisandro Dalcin <<a href="mailto:dalcinl@gmail.com">dalcinl@gmail.com</a>> wrote:<br>
>><br>
>> 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.<br>
><br>
> Okay, anyone else not want me to try and do this?<br>
><br>
<br>
</div>Could you elaborate a bit more what exactly are you going to try? Are<br>
you proposing that rectangular matrices should always have a square<br>
diagonal block?<br><div><div class="h5"></div></div></blockquote><div><br>Yes, I want to go into the code and every place we query the diagonal block, replace<br>cMap by rMap. This is the only thing that makes sense to me. By my definitions, diagonal<br>
blocks are disjoint, so they must be square. Since we normally own rows, it makes sense<br>to use row blocking to define them.<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">
>><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>
><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<br>
> should<br>
> be the same as the layout of x. This necessitates a symmetric partition, and<br>
> is<br>
> probably at the root of our problems with rectangular matrices.<br>
><br>
> Matt<br>
><br>
>><br>
>> On Tue, May 26, 2009 at 10:35 AM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>><br>
>> wrote:<br>
>> > I think there is an inconsistency in our definitions. If you create a<br>
>> > 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,<br>
>> > 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<br>
>> > weird<br>
>> > shape,<br>
>> > or we can change all queries for the diagonal block to use only the<br>
>> > 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<br>
>> > experiments<br>
>> > is infinitely more interesting than any results to which their<br>
>> > experiments<br>
>> > lead.<br>
>> > -- Norbert Wiener<br>
>> ><br>
>><br>
>><br>
>><br>
>> --<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>
><br>
><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>--<br>
<div><div></div><div class="h5">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>
</div></div></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>