On Wed, Nov 30, 2011 at 12:41 AM, Dave Nystrom <span dir="ltr"><<a href="mailto:dnystrom1@comcast.net">dnystrom1@comcast.net</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
I have a linear system in a code that I have interfaced to petsc that is<br>
taking about 80 percent of the run time per timestep.  This linear system is<br>
a symmetric block banded matrix where the blocks are 2x2.  The matrix looks<br>
as follows:<br>
<br>
 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0<br>
1X X                     Y Y Y<br>
2X X X                     Y Y Y<br>
3  X X X                     Y Y Y<br>
4    X X X                     Y Y Y<br>
5      X X X                     Y Y Y<br>
6        X X X                     Y Y Y<br>
7          X X X                     Y Y Y<br>
8            X X X                     Y Y Y<br>
9              X X X                     Y Y Y<br>
0                X X X                     Y Y Y<br>
1                  X X X                     Y Y Y<br>
2                    X X X                     Y Y Y<br>
3Z                     X X X                     Y Y Y<br>
4Z Z                     X X X                     Y Y Y<br>
5Z Z Z                     X X X                     Y Y Y<br>
6  Z Z Z                     X X X                     Y Y Y<br>
7    Z Z Z                     X X X                     Y Y Y<br>
8      Z Z Z                     X X X                     Y Y Y<br>
9        Z Z Z                     X X X                     Y Y Y<br>
0          Z Z Z                     X X X                     Y Y Y<br>
<br>
So in my diagram above, X, Y and Z are 2x2 blocks.  The symmetry of the<br>
matrix requires that X_ij = transpose(X_ji) and Y_ij = transpose(Z_ji).  So<br>
far, I have just input this matrix to petsc without indicating that it was<br>
block banded with 2x2 blocks.  I have also not told petsc that the matrix is<br>
symmetric.  And I have allowed petsc to decide the best way to store the<br>
matrix.<br>
<br>
I can solve this linear system over the course of a run using -ksp_type<br>
preonly -pc_type lu.  But that will not scale very well to larger problems<br>
that I want to solve.  I can also solve this system over the course of a run<br>
using -ksp_type cg -pc_type jacobi -vec_type cusp -mat_type aijcusp.<br>
However, over the course of a run, the iteration count ranges from 771 to<br>
47300.  I have also tried sacusp, ainvcusp, sacusppoly, ilu(k) and icc(k)<br>
with k=0.  The sacusppoly preconditioner fails because of a thrust error<br>
related to an invalid device pointer, if I am remembering correctly.  I<br>
reported this problem to petsc-maint a while back and have also reported it<br>
for the cusp bugtracker.  But it does not appear that anyone has really<br>
looked into the bug.  For the other preconditioners of sacusp, ilu(k) and<br>
icc(k), they do not result in convergence to a solution and the runs fail.<br></blockquote><div><br></div><div>All preconditioners are custom. Have you done a literature search for PCs</div><div>known to work for this problem? Can yu say anything about the spectrum of the</div>
<div>operator? conditioning? what is the principal symbol (if its a PDE)? The pattern</div><div>is not enough to recommend a PC.</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;">

I'm wondering if there are suggestions of other preconditioners in petsc that<br>
I should try.  The only third party package that I have tried is the<br>
txpetscgpu package.  I have not tried hypre or any of the multigrid<br>
preconditioners yet.  I'm not sure how difficult it is to try those<br>
packages.  Anyway, so far I have not found a preconditioner available in<br>
petsc that provides a robust solution to this problem and would be interested<br>
in any suggestions that anyone might have of things to try.<br>
<br>
I'd be happy to provide additional info and am planning on packaging up a<br>
couple of examples of the matrix and rhs for people I am interacting with at<br>
Tech-X and EMPhotonics.  So I'd be happy to provide the matrix examples for<br>
this forum as well if anyone wants a copy.<br>
<br>
Thanks,<br>
<br>
Dave<br>
<font color="#888888"><br>
--<br>
Dave Nystrom<br>
<br>
phone: <a href="tel:505-661-9943" value="+15056619943">505-661-9943</a> (home office)<br>
       <a href="tel:505-662-6893" value="+15056626893">505-662-6893</a> (home)<br>
skype: dave.nystrom76<br>
email: <a href="mailto:dnystrom1@comcast.net">dnystrom1@comcast.net</a><br>
smail: 219 Loma del Escolar<br>
       Los Alamos, NM 87544<br>
</font></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<br>