<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Sat, Dec 28, 2013 at 9:01 AM, Stephan Kramer <span dir="ltr"><<a href="mailto:s.kramer@imperial.ac.uk" target="_blank">s.kramer@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">On 19/12/13 21:03, Matthew Knepley wrote:<div class="im"><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
On Sat, Nov 30, 2013 at 12:13 PM, Stephan Kramer <<a href="mailto:s.kramer@imperial.ac.uk" target="_blank">s.kramer@imperial.ac.uk</a> <mailto:<a href="mailto:s.kramer@imperial.ac.uk" target="_blank">s.kramer@imperial.ac.<u></u>uk</a>>> wrote:<br>
<br>
Dear petsc-devs<br>
<br>
Trying to write a test case for MatZeroRowColumns (as requested by Matt) in order to extend it to MPIBAIJ, I came across a bug in the current implementation for MPIAIJ. The attached test code<br>
assembles a simple 2D Laplacian (with the standard 4,-1,-1,-1 stencil) on the following grid:<br>
<br>
0 1 2 3<br>
4 5 6 7<br>
8 9 10 11<br>
<br>
with boundary conditions applied to nodes 0,1,2,3,4 and 8 with values (stored in Vec x) 0., 1., 2., 3., 4. and 8. For simplicity I've kept a zero rhs (before the call to MatZeroRowsColumns) After<br>
the call to MatZeroRowColumns I therefore expect the following rhs:<br>
<br>
0. 1. 2. 3. 4. 5. 2. 3. 8. 8. 0. 0.<br>
<br>
The entries at 0,1,2,3,4 and 8 correspond to the boundary values (multiplied by a diagonal value of 1.0) and the values at 5,6,7 and 9 to the the lifted values: 1.+4.=5.,2.,3. and 8. However, if<br>
you run the attached example that's not what you get, the lifted values at 6 and 7 are wrong. This works fine in serial (you have to change nlocal=4 to get the same grid).<br>
<br>
I think this is caused by the way that entries in x are communicated with other processes that have off-diagonal entries referencing those. In line 944 of mpiaij.c it uses a VecScatter with<br>
ADD_VALUES to do this, but it uses the local work vector l->lvec without zeroing it first. This means it will still have values in it of whatever it was used for before. In this case that would be<br>
the MatMult of line 125 of test.c (in fact you can "fix" the example by adding the line 136-139 where it does a spurious MatMult with a zero vector).<br>
<br>
The attached code could be easily extended to test block matrices (simply increasing bs). The commented out lines in the assembly were basically to make the problem a bit more challenging (have an<br>
asymmetric stencil and the setting of off-process boundary nodes). The expected results for MPIBAIJ can be obtained by comparing with the MPIAIJ results (once that works correctly).<br>
<br>
<br>
Okay, I have put this in knepley/feature-mat-<u></u>zerorowscols-baij and pushed to next. I also fixed the bug that you found. It should<br>
not take much over the holidays to get the BAIJ functionality you want.<br>
<br>
Thanks,<br>
<br>
Matt<br>
</blockquote>
<br></div>
Great. Let me know if there's anything else I can do to help. Just a small heads up: the new src/mat/examples/tests/ex18.c has an uninitialized bool nonlocalBC<br></blockquote><div><br></div><div>Okay, I have checked in preliminary code. It passes the tests for Mat ex18. I have merged it to next. Let me know how</div>
<div>it works for you.</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">
Cheers<span class="HOEnZb"><font color="#888888"><br>
Stephan<br>
</font></span></blockquote></div><br><br clear="all"><div><br></div>-- <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
</div></div>