<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Sep 16, 2015 at 7:18 AM, Lawrence Mitchell <span dir="ltr"><<a href="mailto:lawrence.mitchell@imperial.ac.uk" target="_blank">lawrence.mitchell@imperial.ac.uk</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">-----BEGIN PGP SIGNED MESSAGE-----<br>
Hash: SHA1<br>
<span class=""><br>
On 16/09/15 12:40, Matthew Knepley wrote:<br>
> On Tue, Sep 15, 2015 at 9:05 PM, Adrian Croucher<br>
</span>> <<a href="mailto:a.croucher@auckland.ac.nz">a.croucher@auckland.ac.nz</a> <mailto:<a href="mailto:a.croucher@auckland.ac.nz">a.croucher@auckland.ac.nz</a>>><br>
<span class="">> wrote:<br>
><br>
> hi<br>
><br>
> I have a test code (attached) that sets up a finite volume mesh<br>
> using DMPlex, with 2 degrees of freedom per cell.<br>
><br>
> I then create a matrix using DMCreateMatrix(), having used<br>
> DMSetMatType() to set the matrix type to MATBAIJ or MATMPIBAIJ, to<br>
> take advantage of the block structure.<br>
><br>
> This works OK and gives me the expected matrix structure when I<br>
> run on > 1 processor, but neither MATBAIJ or MATMPIBAIJ works if I<br>
> run it in serial:<br>
><br>
><br>
> Plex should automatically discover the block size from the Section.<br>
> If not, it uses block size 1. I have to look at the example to see<br>
> why the discovery is not working. Do you have any constraints?<br>
<br>
</span>It looks like block discovery in parallel effectively always<br>
determines a block size of 1.  Running with -mat_view ::ascii_info gives:<br>
<br>
Mat Object: 2 MPI processes<br>
  type: mpibaij<br>
  rows=20, cols=20<br>
  total: nonzeros=112, allocated nonzeros=112<br>
  total number of mallocs used during MatSetValues calls =0<br>
      block size is 1<br>
                    ^^^<br>
<br>
The block size discovery does this:<br>
<br>
        for (p = pStart; p < pEnd; ++p) {<br>
          ierr = PetscSectionGetDof(sectionGlobal, p, &dof);CHKERRQ(ierr);<br>
          ierr = PetscSectionGetConstraintDof(sectionGlobal, p,<br>
&cdof);CHKERRQ(ierr);<br>
          if (dof-cdof) {<br>
            if (bs < 0) {<br>
              bs = dof-cdof;<br>
            } else if (bs != dof-cdof) {<br>
              /* Layout does not admit a pointwise block size */<br>
              bs = 1;<br>
              break;<br>
            }<br>
          }<br>
        }<br>
<br>
In parallel, some of the dofs are remote (encoded as negative).  So we<br>
run through seeing (dof - cdof) == 2, until we hit a "remote" point at<br>
when we see (dof - cdof) != 2 and then we break out and set bs = 1.<br>
<br>
In serial, we correctly determine bs == 2.  But then<br>
DMGetLocalToGlobalMapping always does<br>
<br>
      ierr = ISLocalToGlobalMappingCreate(PETSC_COMM_SELF, 1,size,<br>
ltog, PETSC_OWN_POINTER, &dm->ltogmap);CHKERRQ(ierr);<br>
<br>
<br>
i.e. is set with block size == 1.<br>
<br>
So there are two bugs here.<br>
<br>
1. In parallel, block size discovery in Mat creation has gone wrong<br></blockquote><div><br></div><div>Crap. Hoist on my own petard. Okay I will fix this.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
2. (Always), the lgmap has block size of 1, irrespective of the<br>
discovered block size from the section.<br></blockquote><div><br></div><div>Yep. This can also be fixed. It should work regardless, but would be better with blocking.</div><div><br></div><div>  Thanks</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Lawrence<br>
<br>
<br>
-----BEGIN PGP SIGNATURE-----<br>
Version: GnuPG v1<br>
<br>
iQEcBAEBAgAGBQJV+V4tAAoJECOc1kQ8PEYvqO4IAN4+oIgtBevvDAughPgUVOzq<br>
kESJkb0Bx4a7+y47IPsY/SOuOMjVgErz2SO3tGd7+K/U5fstojJbIC7zZqPIERn0<br>
S68lH3s+y1pVqmcIMFIorz6is+u46M2xIGS6MLe6d0DluslyThaaA7lMhuIvJkIX<br>
FC8QmtqCsIvHv10VuNll/81UJ3pXSZ+E81+Rs6pBhoGMCnRSXdMgFEtfWBBGL2JR<br>
byEOAueTmY+YZXJ5JINxsHG1C1lep5wyYfAERDiRaD9bhCsE4mf9z/yFdvhiv0e3<br>
JC5WB3Q2/t2dRXhKNOSxQJd5mBKGhHgKptTjzmnW8HF0uxCdy7XpKZ6vS4Kt8Gk=<br>
=2JQe<br>
-----END PGP SIGNATURE-----<br>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">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</div>
</div></div>