On Tue, May 26, 2009 at 2:15 PM, Matthew Knepley <knepley at gmail.com> wrote: > On Tue, May 26, 2009 at 12:07 PM, Lisandro Dalcin <dalcinl at gmail.com> > wrote: >> >> 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. > > Okay, anyone else not want me to try and do this? >
Could you elaborate a bit more what exactly are you going to try? Are you proposing that rectangular matrices should always have a square diagonal block? >> >> 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? > > No, I? do not think that this should work, because the Krylov methods impose > a more stringent condition on the matrix A, namely that the layout of Ax > should > be the same as the layout of x. This necessitates a symmetric partition, and > is > probably at the root of our problems with rectangular matrices. > > ? Matt > >> >> 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 > > > > -- > 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