Changing the matrix types worked. When solving AX=B, X and B have to be dense and A can be a sparse format? I'd like to avoid storing the identity as a dense format, if possible. By the way, the man page for MATAIJ says it's a sparse format <a href="http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/docs/manualpages/Mat/MATAIJ.html">http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/docs/manualpages/Mat/MATAIJ.html</a>. <div>
<br></div><div>I must ask, is there somewhere in the documentation I should be looking to figure these sorts of questions out? I feel like a pest, but I'm not sure where else to look for answers.</div><div><br></div><div>
<br><div class="gmail_quote">On Thu, Jul 7, 2011 at 12:20 AM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov">jedbrown@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div class="im"><div class="gmail_quote">On Wed, Jul 6, 2011 at 09:53, Adam Byrd <span dir="ltr"><<a href="mailto:adam1.byrd@gmail.com" target="_blank">adam1.byrd@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
I do not want the transpose, but that is what is coming out as the solution. Both my tests show the correct original matrix, and both yield the transpose of the solution verified against Mathematica and the website from which I pulled the 4x4 matrix in test.cpp. Could it be related to the factoring?</blockquote>
</div><br></div><div>The problem is that you are trying to solve with an AIJ matrix instead of a DENSE matrix. MatMatSolve() calls MatGetArray(). Since this is serial and you preallocated enough for a dense matrix, it goes through with no memory errors. The storage format for MATAIJ (that happens to be dense) and MATDENSE is transposed, so you see the transpose. If you make the solution matrix dense, then it will work correctly.</div>
<div><br></div><div><br></div><div><div>--- a/test.cpp</div><div>+++ b/test.cpp</div><div>@@ -48,7 +48,9 @@ int main(int argc,char **args)</div><div> ierr = MatAssemblyEnd(identityMat, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div>
<div> </div><div> //Setup inverseMat</div><div>- ierr = MatDuplicate(testMat, MAT_DO_NOT_COPY_VALUES, &inverseMat);CHKERRQ(ierr);</div><div>+ ierr = MatCreate(PETSC_COMM_WORLD, &inverseMat);CHKERRQ(ierr);</div>
<div>+ ierr = MatSetSizes(inverseMat, PETSC_DECIDE, PETSC_DECIDE, rows, columns);CHKERRQ(ierr);</div><div>+ ierr = MatSetType(inverseMat, MATDENSE);CHKERRQ(ierr);</div><div> ierr = MatAssemblyBegin(inverseMat, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div>
<div> ierr = MatAssemblyEnd(inverseMat, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);</div></div><div><br></div><div><br></div><div><br></div><div>PETSc should be more careful about matching types here.</div>
</blockquote></div><br></div>