The attached example exposes a bug in some implementations of MatGetSubMatrices():<div><div>The �extracted submatrices have row such as these:�</div><div><div>row 0: (0, 4) �(3, -1) �(1, -1)�</div><div>row 1: (0, -1) �(1, 4) �(4, -1) �(2, -1)�</div>
</div><div>Because of there are inversions among the column indices, binary search will fail for XXXAIJ matrices.</div><div>This is the result of using unsorted column ISs in MatGetSubMatrices().</div><div>According to our documentation, it's legal. �Still, some (but not all) implementations (e.g., MPIDENSE, SEQAIJ)�</div>
<div>will reject unsorted column indices. �</div><div><br></div><div><br></div><div>Note that this subroutine�is used by PCASM, where the user's subdomain ISs are�</div><div>used to extract the subdomain�submatrices. �It looks like the problem never got exposed because</div>
<div>by default the subdomain ISs are sorted by ASM internally, unless the user explicitly sets a flag�</div><div>(and, incidentally, it doesn't look like the flag can be set on the command line).</div>
<div><br></div><div>The problem seems to exist for the following matrix types:</div><div>MPIAIJ, MPIBAIJ, MPISBAIJ</div><div><br></div><div><div>Unless someone objects I'm going to push a fix to petsc-dev that will check that�</div>
<div>all column ISs are sorted in the three implementations above and will throw an error otherwise.�</div></div><div><br></div><div>Dmitry.</div><div><br></div><div>PS For those that want to run this test:</div><div><div>
make runsubmat_test NP=2 ARGS="-mat_type mpiaij -sort_rows -show_inversions"</div><div><br></div><div>NP must be 1 if using a sequential -mat_type.</div><div><br></div><div>-show_inversions requires MatGetRow(), which doesn't work with symmetric types.</div>
</div><div><br></div><div>-sort_rows will sort the row ISs, but not the column ISs (that's how one gets past the first line of defense�</div><div>for matrix types like MPIDENSE or SEQAIJ, amazingly, which bail as soon as they see unsorted row ISs);�</div>
<div>symmetric matrix type will not like this, since they need identical row and column ISs</div><div><br></div>
<div>Here's the output of make runsubmat_test NP=2 ARGS="-mat_type mpiaij -sort_rows -show_inversions" (I can't figure out what's causing the memory leak):</div><div><div>Original matrix</div><div>Matrix Object: 1 MPI processes</div>
<div>� type: seqaij</div><div>row 0: (0, 4) �(1, -1) �(4, -1)�</div><div>row 1: (0, -1) �(1, 4) �(2, -1) �(5, -1)�</div><div>row 2: (1, -1) �(2, 4) �(3, -1) �(6, -1)�</div><div>row 3: (2, -1) �(3, 4) �(7, -1)�</div><div>
row 4: (0, -1) �(4, 4) �(5, -1) �(8, -1)�</div>
<div>row 5: (1, -1) �(4, -1) �(5, 4) �(6, -1) �(9, -1)�</div><div>row 6: (2, -1) �(5, -1) �(6, 4) �(7, -1) �(10, -1)�</div><div>row 7: (3, -1) �(6, -1) �(7, 4) �(11, -1)�</div><div>row 8: (4, -1) �(8, 4) �(9, -1) �(12, -1)�</div>
<div>row 9: (5, -1) �(8, -1) �(9, 4) �(10, -1) �(13, -1)�</div><div>row 10: (6, -1) �(9, -1) �(10, 4) �(11, -1) �(14, -1)�</div><div>row 11: (7, -1) �(10, -1) �(11, 4) �(15, -1)�</div><div>row 12: (8, -1) �(12, 4) �(13, -1)�</div>
<div>row 13: (9, -1) �(12, -1) �(13, 4) �(14, -1)�</div><div>row 14: (10, -1) �(13, -1) �(14, 4) �(15, -1)�</div><div>row 15: (11, -1) �(14, -1) �(15, 4)�</div><div>[0:1]: Number of subdomains: 2:</div><div>[0:1]: Subdomain row IS 0:</div>
<div>Number of indices in set 12</div><div>0 0</div><div>1 1</div><div>2 2</div><div>3 3</div><div>4 4</div><div>5 5</div><div>6 6</div><div>7 7</div><div>8 8</div><div>9 9</div><div>10 10</div><div>11 11</div><div>[0:1]: Subdomain col IS 0:</div>
<div>Number of indices in set 12</div><div>0 0</div><div>1 4</div><div>2 8</div><div>3 1</div><div>4 5</div><div>5 9</div><div>6 2</div><div>7 6</div><div>8 10</div><div>9 3</div><div>10 7</div><div>11 11</div><div>[0:1]: Submatrix 0:</div>
<div>Matrix Object: 1 MPI processes</div><div>� type: seqaij</div><div>row 0: (0, 4) �(3, -1) �(1, -1)�</div><div>row 1: (0, -1) �(3, 4) �(6, -1) �(4, -1)�</div><div>row 2: (3, -1) �(6, 4) �(9, -1) �(7, -1)�</div><div>row 3: (6, -1) �(9, 4) �(10, -1)�</div>
<div>row 4: (0, -1) �(1, 4) �(4, -1) �(2, -1)�</div><div>row 5: (3, -1) �(1, -1) �(4, 4) �(7, -1) �(5, -1)�</div><div>row 6: (6, -1) �(4, -1) �(7, 4) �(10, -1) �(8, -1)�</div><div>row 7: (9, -1) �(7, -1) �(10, 4) �(11, -1)�</div>
<div>row 8: (1, -1) �(2, 4) �(5, -1)�</div><div>row 9: (4, -1) �(2, -1) �(5, 4) �(8, -1)�</div><div>row 10: (7, -1) �(5, -1) �(8, 4) �(11, -1)�</div><div>row 11: (10, -1) �(8, -1) �(11, 4)�</div><div>***Inversion in row 0: col[2] = 1 < 3 = col[1]</div>
<div>***Inversion in row 1: col[3] = 4 < 6 = col[2]</div><div>***Inversion in row 2: col[3] = 7 < 9 = col[2]</div><div>***Inversion in row 4: col[3] = 2 < 4 = col[2]</div><div>***Inversion in row 5: col[1] = 1 < 3 = col[0]</div>
<div>***Inversion in row 5: col[4] = 5 < 7 = col[3]</div><div>***Inversion in row 6: col[1] = 4 < 6 = col[0]</div><div>***Inversion in row 6: col[4] = 8 < 10 = col[3]</div><div>***Inversion in row 7: col[1] = 7 < 9 = col[0]</div>
<div>***Inversion in row 9: col[1] = 2 < 4 = col[0]</div><div>***Inversion in row 10: col[1] = 5 < 7 = col[0]</div><div>***Inversion in row 11: col[1] = 8 < 10 = col[0]</div><div>[0:1]: Subdomain row IS 1:</div>
<div>
Number of indices in set 12</div><div>0 4</div><div>1 5</div><div>2 6</div><div>3 7</div><div>4 8</div><div>5 9</div><div>6 10</div><div>7 11</div><div>8 12</div><div>9 13</div><div>10 14</div><div>11 15</div><div>[0:1]: Subdomain col IS 1:</div>
<div>Number of indices in set 12</div><div>0 4</div><div>1 8</div><div>2 12</div><div>3 5</div><div>4 9</div><div>5 13</div><div>6 6</div><div>7 10</div><div>8 14</div><div>9 7</div><div>10 11</div><div>11 15</div><div>[0:1]: Submatrix 1:</div>
<div>Matrix Object: 1 MPI processes</div><div>� type: seqaij</div><div>row 0: (0, 4) �(3, -1) �(1, -1)�</div><div>row 1: (0, -1) �(3, 4) �(6, -1) �(4, -1)�</div><div>row 2: (3, -1) �(6, 4) �(9, -1) �(7, -1)�</div><div>row 3: (6, -1) �(9, 4) �(10, -1)�</div>
<div>row 4: (0, -1) �(1, 4) �(4, -1) �(2, -1)�</div><div>row 5: (3, -1) �(1, -1) �(4, 4) �(7, -1) �(5, -1)�</div><div>row 6: (6, -1) �(4, -1) �(7, 4) �(10, -1) �(8, -1)�</div><div>row 7: (9, -1) �(7, -1) �(10, 4) �(11, -1)�</div>
<div>row 8: (1, -1) �(2, 4) �(5, -1)�</div><div>row 9: (4, -1) �(2, -1) �(5, 4) �(8, -1)�</div><div>row 10: (7, -1) �(5, -1) �(8, 4) �(11, -1)�</div><div>row 11: (10, -1) �(8, -1) �(11, 4)�</div><div>***Inversion in row 0: col[2] = 1 < 3 = col[1]</div>
<div>***Inversion in row 1: col[3] = 4 < 6 = col[2]</div><div>***Inversion in row 2: col[3] = 7 < 9 = col[2]</div><div>***Inversion in row 4: col[3] = 2 < 4 = col[2]</div><div>***Inversion in row 5: col[1] = 1 < 3 = col[0]</div>
<div>***Inversion in row 5: col[4] = 5 < 7 = col[3]</div><div>***Inversion in row 6: col[1] = 4 < 6 = col[0]</div><div>***Inversion in row 6: col[4] = 8 < 10 = col[3]</div><div>***Inversion in row 7: col[1] = 7 < 9 = col[0]</div>
<div>***Inversion in row 9: col[1] = 2 < 4 = col[0]</div><div>***Inversion in row 10: col[1] = 5 < 7 = col[0]</div><div>***Inversion in row 11: col[1] = 8 < 10 = col[0]</div><div>[0]Total space allocated 16 bytes</div>
<div>[ 0]16 bytes main() line 99 in submat_test.c</div></div><div><br></div>
</div>