Problem with rectangular matrices

[Lisandro Dalcin]

The only real issue I see here is that the term "diagonal block" could
lead to misinterpretation. But the behaviour you commented is the one
I would expect.

The rest of this mail is perhaps off-topic, but it is somewhat related.

The other issues is that some Mat preallocation routines are lazy
about properly setting the cmap. I had to implement some hackery in
petsc4py to make sure that created Mat's (at least these created from
petsc4py calls) have a cmap consistent with the local/global
row/column sizes.

Moreover, supose you create a global square matrix with the following
row/col distribution (I know, it is perhaps a nonsense, but IMHO is a
valid distribution)
from petsc4py import PETSc

size = PETSc.COMM_WORLD.getSize()
assert size == 2
rank = PETSc.COMM_WORLD.getRank()
Ms = [2, 6]
Ns = [6, 2]

m = Ms[rank]
n = Ns[rank]

A = PETSc.Mat().createAIJ([(m, None), (n, None)]) # None is
intrepreted as PETSc.DECIDE
for i in range(8):
    A[i,i] = 1.0/(i+1)
A.assemble()

x, b = A.getVecs()
x.set(1)
A.mult(x,b)

Up to now, all works as expected, all rmap/cmap are consistent, the
'x' and 'b' vecs have the same global size, though different local
sizes, and the MatMult succeed. All fine... But now, suppose you do:

ksp = PETSc.KSP().create()
ksp.setOperators(A)
ksp.pc.type = 'none'
ksp.solve(b,x)

You are lost, this does not work (even if you try with nonzero initial
guess)... But it should, right?



On Tue, May 26, 2009 at 10:35 AM, Matthew Knepley <knepley at gmail.com> wrote:
> I think there is an inconsistency in our definitions. If you create a matrix
> which
> is 4x8 on 2 procs, then you end up with a row and column map like this:
>
>     rowMap           colMap
>   [0..2) [2..4)      [0..4) [4..8)
>
> However, we use the column map for queries which then give nonsense. The
> most
> damaging is in MatSetValues(). If I give a value for (2,2) on proc 1, then
> since
> 2 < cstart = 4, this value is considered "off-diagonal" and can cause a
> malloc(). So
> I think we have a choice. We can define "diagonal block" to have this weird
> shape,
> or we can change all queries for the diagonal block to use only the rowMap.
> I vote
> for the latter, but does anyone have an objection?
>
>    Matt
>
> --
> What most experimenters take for granted before they begin their experiments
> is infinitely more interesting than any results to which their experiments
> lead.
> -- Norbert Wiener
>



-- 
Lisandro Dalcín
---------------
Centro Internacional de Métodos Computacionales en Ingeniería (CIMEC)
Instituto de Desarrollo Tecnológico para la Industria Química (INTEC)
Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET)
PTLC - Güemes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594