On Thu, Jun 18, 2009 at 10:48 AM, Christian Klettner <span dir="ltr">&lt;<a href="mailto:christian.klettner@ucl.ac.uk">christian.klettner@ucl.ac.uk</a>&gt;</span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
Hi Matt,<br>
Thank you for the reply. I have now gotten MUMPS to work with the<br>
following implementation:<br>
<br>
MatFactorInfoInitialize(&amp;info);<br>
MatGetFactor(A,MAT_SOLVER_MUMPS,MAT_FACTOR_LU,&amp;B);<br>
MatLUFactorSymbolic(B,A, rows_ck_is, cols_ck_is,&amp;info);<br>
MatLUFactorNumeric(B,A,&amp;info);<br>
<br>
for(){         ///temporal loop<br>
/*coding*/<br>
<div class="im">ierr=MatSolve(A,vecr,vecu);CHKERRQ(ierr);<br>
}<br>
<br>
</div>I explicitly use MUMPS to factor the matrix. However, from my<br>
implementation I do not see what package is being used to solve Ax=b.<br>
Could I be using PETSc to solve the system?? Is there any way of knowing<br>
if MUMPS is actually being used to solve the system?<br>
I am not passing in any arguments to use MUMPS explicitly from the command<br>
line, have<br>
not put in any MUMPS header files and I do not use any other MUMPS<br>
functions (e.g. MatCreate_AIJMUMPS()).</blockquote><div><br>Yes, it is determined from the type. However, I think this is completely the wrong way<br>to program. Why hardcode yourself into a corner? It is hard to change packages, to<br>
wrap this in a Krylov solver, compare with other PCs. This is crazy. I recommend<br><br>  KSPCreate(comm, &amp;ksp);<br>  KSPSetOperators(ksp, A, A, ...);<br>  KSPSetFromOptions(ksp);<br>  KSPSolve(ksp, b, x);<br><br>Now with options you can control everything<br>
<br>  -ksp_type preonly -pc_type lu -mat <b>-pc_factor_mat_solver_package mumps</b><br> <br>will reproduce what you have, but you can test it against<br><br> -ksp_type gmres -pc_type asm -sub_pc_type ilu<br><br>To see what is happening you can always use<br>
<br>  -ksp_view<br><br>   Matt<br><br></div><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;"><br>
Thank you for any advice,<br>
Best regards,<br>
<font color="#888888">Christian Klettner<br>
</font><div><div></div><div class="h5"><br>
<br>
&gt; On Tue, Jun 16, 2009 at 7:40 AM, Christian Klettner &lt;<br>
&gt; <a href="mailto:christian.klettner@ucl.ac.uk">christian.klettner@ucl.ac.uk</a>&gt; wrote:<br>
&gt;<br>
&gt;&gt; Hi Matt,<br>
&gt;&gt; I have tried to use MUMPS and came across the following problems. My<br>
&gt;&gt; application solves Ax=b about 10000-100000 times with A remaining<br>
&gt;&gt; constant.<br>
&gt;&gt; When I tried to use it through KSP I was not finding good performance.<br>
&gt;&gt; Could this be because it was refactoring etc. at each time step?<br>
&gt;&gt; With this in mind I have tried to implement the following:<br>
&gt;<br>
&gt;<br>
&gt; 1) I have no idea what you would mean by good performance. Always, always,<br>
&gt; always<br>
&gt; send the -log_summary.<br>
&gt;<br>
&gt; 2) A matrix is never refactored during the KSP iteration<br>
&gt;<br>
&gt;<br>
&gt;&gt; A is a parallel matrix created with MatCreateMPIAIJ().<br>
&gt;&gt; rows is a list of the global row numbers on the process.cols is a list<br>
&gt;&gt; of<br>
&gt;&gt; the global columns that are on the process.<br>
&gt;&gt; I run the program with ./mpiexec -n 4 ./ex4<br>
&gt;&gt; -pc_factor_mat_solver_package<br>
&gt;&gt; mumps<br>
&gt;&gt;<br>
&gt;&gt; (1)     MatFactorInfo *info;<br>
&gt;&gt; (2)     ierr=MatFactorInfoInitialize(info);CHKERRQ(ierr);<br>
&gt;<br>
&gt;<br>
&gt; You have passed a meaningless pointer here.<br>
&gt;<br>
&gt; (1)     MatFactorInfo info;<br>
&gt; (2)     ierr=MatFactorInfoInitialize(&amp;info);CHKERRQ(ierr);<br>
&gt;<br>
&gt;    Matt<br>
&gt;<br>
&gt;<br>
&gt;&gt; (3)<br>
&gt;&gt; ierr=MatGetFactor(A,MAT_SOLVER_MUMPS,MAT_FACTOR_LU,&amp;F);CHKERRQ(ierr);<br>
&gt;&gt; (4)     ierr=MatLUFactorSymbolic(F, A, rows, cols,info);CHKERRQ(ierr);<br>
&gt;&gt; (5)     MatLUFactorNumeric(F, A,info);CHKERRQ(ierr);<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;        for(){  ///TEMPORAL LOOP<br>
&gt;&gt;<br>
&gt;&gt;                /*CODING*/<br>
&gt;&gt;<br>
&gt;&gt;                ierr=MatSolve(A,vecr,vecu);CHKERRQ(ierr);<br>
&gt;&gt;        }<br>
&gt;&gt;<br>
&gt;&gt; I get the following error messages:<br>
&gt;&gt;<br>
&gt;&gt; Line (2) above gives the following error message:<br>
&gt;&gt;<br>
&gt;&gt; [0]PETSC ERROR: --------------------- Error Message<br>
&gt;&gt; ------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: [0]PETSC ERROR: Null argument, when expecting valid<br>
&gt;&gt; pointer!<br>
&gt;&gt; [0]PETSC ERROR: Trying to zero at a null pointer!<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: Petsc Release Version 3.0.0, Patch 4, Fri Mar  6<br>
&gt;&gt; 14:46:08<br>
&gt;&gt; CST 2009<br>
&gt;&gt; [0]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
&gt;&gt; [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>
&gt;&gt; [0]PETSC ERROR: See docs/index.html for manual pages.<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: ./ex4 on a linux-gnu named christian-desktop by<br>
&gt;&gt; christian<br>
&gt;&gt; Tue Jun 16 13:33:33 2009<br>
&gt;&gt; [0]PETSC ERROR: Libraries linked from<br>
&gt;&gt; /home/christian/Desktop/petsc-3.0.0-p4/linux-gnu-c-debug/lib<br>
&gt;&gt; [0]PETSC ERROR: Configure run at Mon Jun 15 17:05:31 2009<br>
&gt;&gt; [0]PETSC ERROR: Configure options --with-cc=&quot;gcc -fPIC&quot;<br>
&gt;&gt; --download-mpich=1<br>
&gt;&gt; --download-f-blas-lapack --download-scalapack --download-blacs<br>
&gt;&gt; --download-mumps --download-parmetis --download-hypre<br>
&gt;&gt; --download-triangle<br>
&gt;&gt; --with-shared=0<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: PetscMemzero() line 189 in src/sys/utils/memc.c<br>
&gt;&gt; [0]PETSC ERROR: MatFactorInfoInitialize() line 7123 in<br>
&gt;&gt; src/mat/interface/matrix.c<br>
&gt;&gt; [0]PETSC ERROR: main() line 1484 in src/dm/ao/examples/tutorials/ex4.c<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 85) - process 0[cli_0]:<br>
&gt;&gt; aborting job:<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 85) - process 0<br>
&gt;&gt; --------------------- Error Message ------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: Null argument, when expecting valid pointer!<br>
&gt;&gt; [2]PETSC ERROR: Trying to zero at a null pointer!<br>
&gt;&gt; [2]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: Petsc Release Version 3.0.0, Patch 4, Fri Mar  6<br>
&gt;&gt; 14:46:08<br>
&gt;&gt; CST 2009<br>
&gt;&gt; [2]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
&gt;&gt; [2]PETSC ERROR: See docs/faq.html[0]0:Return code = 85<br>
&gt;&gt; [0]1:Return code = 0, signaled with Interrupt<br>
&gt;&gt; [0]2:Return code = 0, signaled with Interrupt<br>
&gt;&gt; [0]3:Return code = 0, signaled with Interrupt<br>
&gt;&gt;<br>
&gt;&gt; /////////////////////////////////////////////////////////////////<br>
&gt;&gt;<br>
&gt;&gt; Line (3) gives the following error messages:<br>
&gt;&gt;<br>
&gt;&gt; [0]PETSC ERROR: --------------------- Error Message<br>
&gt;&gt; ------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: Invalid argument!<br>
&gt;&gt; [0]PETSC ERROR: Wrong type of object: Parameter # 2!<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: [2]PETSC ERROR: --------------------- Error Message<br>
&gt;&gt; ------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: Invalid argument!<br>
&gt;&gt; [2]PETSC ERROR: Wrong type of object: Parameter # 2!<br>
&gt;&gt; [2]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: Petsc Release Version 3.0.0, Patch 4, Fri Mar  6<br>
&gt;&gt; 14:46:08<br>
&gt;&gt; CST 2009<br>
&gt;&gt; [2]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
&gt;&gt; [2]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>
&gt;&gt; [2]PETSC ERROR: See docs/index.html for manual pages.<br>
&gt;&gt; [2]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: ./ex4 on a linux-gnu named christian-desktop by<br>
&gt;&gt; christian<br>
&gt;&gt; Tue Jun 16 13:37:36 2009<br>
&gt;&gt; [2]PETSC ERROR: Libraries linked from<br>
&gt;&gt; /home/christian/Desktop/petsc-3.0.0-p4/linux-gnu-c-debug/lib<br>
&gt;&gt; [2]PETSC ERROR: Configure run at Mon Jun 15 17:05:31 2009<br>
&gt;&gt; [2]PETSC ERROR: Configure options --with-cc=&quot;gcc -fPIC&quot;<br>
&gt;&gt; --download-mpich=1<br>
&gt;&gt; --download-f-blas-lapack --download-scalapack --download-blacs<br>
&gt;&gt; --downlPetsc Release Version 3.0.0, Patch 4, Fri Mar  6 14:46:08 CST<br>
&gt;&gt; 2009<br>
&gt;&gt; [0]PETSC ERROR: See docs/changes/index.html for recent updates.<br>
&gt;&gt; [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>
&gt;&gt; [0]PETSC ERROR: See docs/index.html for manual pages.<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: ./ex4 on a linux-gnu named christian-desktop by<br>
&gt;&gt; christian<br>
&gt;&gt; Tue Jun 16 13:37:36 2009<br>
&gt;&gt; [0]PETSC ERROR: Libraries linked from<br>
&gt;&gt; /home/christian/Desktop/petsc-3.0.0-p4/linux-gnu-c-debug/lib<br>
&gt;&gt; [0]PETSC ERROR: Configure run at Mon Jun 15 17:05:31 2009<br>
&gt;&gt; [0]PETSC ERROR: Configure options --with-cc=&quot;gcc -fPIC&quot;<br>
&gt;&gt; --download-mpich=1<br>
&gt;&gt; --download-f-blas-lapack --download-scalapack --download-blacs<br>
&gt;&gt; --download-mumps --download-parmetis --download-hypre<br>
&gt;&gt; --download-triangle<br>
&gt;&gt; --with-shared=0<br>
&gt;&gt; [0]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [0]PETSC ERROR: MatLUFactorSymbolic() line 2311 in<br>
&gt;&gt; src/mat/interface/matrix.c<br>
&gt;&gt; [0]PETSC ERROR: main() line 148oad-mumps --download-parmetis<br>
&gt;&gt; --download-hypre --download-triangle --with-shared=0<br>
&gt;&gt; [2]PETSC ERROR:<br>
&gt;&gt; ------------------------------------------------------------------------<br>
&gt;&gt; [2]PETSC ERROR: MatLUFactorSymbolic() line 2311 in<br>
&gt;&gt; src/mat/interface/matrix.c<br>
&gt;&gt; [2]PETSC ERROR: main() line 1488 in src/dm/ao/examples/tutorials/ex4.c<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 62) - process 2[cli_2]:<br>
&gt;&gt; aborting job:<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 62) - process 2<br>
&gt;&gt; 8 in src/dm/ao/examples/tutorials/ex4.c<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 62) - process 0[cli_0]:<br>
&gt;&gt; aborting job:<br>
&gt;&gt; application called MPI_Abort(MPI_COMM_WORLD, 62) - process 0<br>
&gt;&gt; [0]0:Return code = 62<br>
&gt;&gt; [0]1:Return code = 0, signaled with Interrupt<br>
&gt;&gt; [0]2:Return code = 62<br>
&gt;&gt; [0]3:Return code = 0, signaled with Interrupt<br>
&gt;&gt;<br>
&gt;&gt; Is the matrix (ie parameter 2) in teh wrong state because it&#39;s not a<br>
&gt;&gt; MUMPS<br>
&gt;&gt; matrix?<br>
&gt;&gt;<br>
&gt;&gt; Any help would be greatly appreciated,<br>
&gt;&gt; Best regards,<br>
&gt;&gt; Christian Klettner<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; &gt; The problem is small enough that you might be able to use MUMPS.<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;   Matt<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; On Fri, Jun 12, 2009 at 9:31 AM, Lisandro Dalcin &lt;<a href="mailto:dalcinl@gmail.com">dalcinl@gmail.com</a>&gt;<br>
&gt;&gt; &gt; wrote:<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;&gt; On Fri, Jun 12, 2009 at 11:13 AM, Christian<br>
&gt;&gt; &gt;&gt; Klettner&lt;<a href="mailto:christian.klettner@ucl.ac.uk">christian.klettner@ucl.ac.uk</a>&gt; wrote:<br>
&gt;&gt; &gt;&gt; &gt; Sorry that I sent this twice. No subject in the first one.<br>
&gt;&gt; &gt;&gt; &gt;<br>
&gt;&gt; &gt;&gt; &gt; Dear PETSc Team,<br>
&gt;&gt; &gt;&gt; &gt; I am writing a CFD finite element code in C. From the<br>
&gt;&gt; discretization<br>
&gt;&gt; &gt;&gt; of<br>
&gt;&gt; &gt;&gt; &gt; the governing equations I have to solve a Poisson type equation<br>
&gt;&gt; which<br>
&gt;&gt; &gt;&gt; is<br>
&gt;&gt; &gt;&gt; &gt; really killing my performance. Which solver/preconditioner from<br>
&gt;&gt; PETSc<br>
&gt;&gt; &gt;&gt; or<br>
&gt;&gt; &gt;&gt; &gt; any external packages would you recommend? The size of my problem<br>
&gt;&gt; is<br>
&gt;&gt; &gt;&gt; from<br>
&gt;&gt; &gt;&gt; &gt; ~30000-100000 DOF per core. What kind of performance would I be<br>
&gt;&gt; able<br>
&gt;&gt; &gt;&gt; to<br>
&gt;&gt; &gt;&gt; &gt; expect with this solver/preconditioner?<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; I would suggest KSPCG. As preconditioner I would use ML or<br>
&gt;&gt; &gt;&gt; HYPRE/BoomerAMG (both are external packages)<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; &gt; I am using a 2*quad core 2.3 GHz Opteron. I have decomposed the<br>
&gt;&gt; domain<br>
&gt;&gt; &gt;&gt; &gt; with Parmetis. The mesh is unstructured.<br>
&gt;&gt; &gt;&gt; &gt; Also, I am writing a code which studies free surface phenomena so<br>
&gt;&gt; the<br>
&gt;&gt; &gt;&gt; mesh<br>
&gt;&gt; &gt;&gt; &gt; is continually changing. Does this matter when choosing a<br>
&gt;&gt; &gt;&gt; &gt; solver/preconditioner? My left hand side matrix (A in Ax=b) does<br>
&gt;&gt; not<br>
&gt;&gt; &gt;&gt; &gt; change in time.<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; ML has a faster setup that BoomerAMG, but the convergence is a bit<br>
&gt;&gt; &gt;&gt; slower. If your A matrix do not change, then likely BoomerAMG will be<br>
&gt;&gt; &gt;&gt; better for you. In any case, you can try both: just build PETSc with<br>
&gt;&gt; &gt;&gt; both packages, then you can change the preconditioner by just passing<br>
&gt;&gt; &gt;&gt; a command line option.<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; &gt;<br>
&gt;&gt; &gt;&gt; &gt; Best regards and thank you in advance,<br>
&gt;&gt; &gt;&gt; &gt; Christian Klettner<br>
&gt;&gt; &gt;&gt; &gt;<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; Disclaimer: the convergence of multigrid preconditioners depends a<br>
&gt;&gt; lot<br>
&gt;&gt; &gt;&gt; on your actual problem. What I&#39;ve suggested is just my limited<br>
&gt;&gt; &gt;&gt; experience in a few problems I&#39;ve run solving electric potentials.<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;&gt; --<br>
&gt;&gt; &gt;&gt; Lisandro Dalcín<br>
&gt;&gt; &gt;&gt; ---------------<br>
&gt;&gt; &gt;&gt; Centro Internacional de Métodos Computacionales en Ingeniería (CIMEC)<br>
&gt;&gt; &gt;&gt; Instituto de Desarrollo Tecnológico para la Industria Química (INTEC)<br>
&gt;&gt; &gt;&gt; Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET)<br>
&gt;&gt; &gt;&gt; PTLC - Güemes 3450, (3000) Santa Fe, Argentina<br>
&gt;&gt; &gt;&gt; Tel/Fax: +54-(0)342-451.1594<br>
&gt;&gt; &gt;&gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; --<br>
&gt;&gt; &gt; What most experimenters take for granted before they begin their<br>
&gt;&gt; &gt; experiments<br>
&gt;&gt; &gt; is infinitely more interesting than any results to which their<br>
&gt;&gt; experiments<br>
&gt;&gt; &gt; lead.<br>
&gt;&gt; &gt; -- Norbert Wiener<br>
&gt;&gt; &gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;<br>
&gt;<br>
&gt; --<br>
&gt; What most experimenters take for granted before they begin their<br>
&gt; experiments<br>
&gt; is infinitely more interesting than any results to which their experiments<br>
&gt; lead.<br>
&gt; -- Norbert Wiener<br>
&gt;<br>
<br>
<br>
</div></div></blockquote></div><br><br clear="all"><br>-- <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>