[petsc-users] How to set up -mat_mumps_icntl13 1 in the code?

[Matthew Knepley]

On Mon, May 14, 2012 at 12:01 PM, Giacomo Mulas
<gmulas at oa-cagliari.inaf.it>wrote:

> Hello.
>
> I am relatively new to petsc and slepc. I am trying to use slepc to solve a
> sparse eigenvalue problem, in particular getting all eigenvalues in an
> interval (spectral slicing). To run in parallel, the documentation states
> that I have to specify an external linear solver, e.g. mumps, and explains
> that I can use the command line options
>
> -st_ksp_type preonly -st_pc_type cholesky \
> -st_pc_factor_mat_solver_**package mumps -mat_mumps_icntl_13 1
>
> to set everything up for it. I was indeed able to get my code to work with
> these command line options (even if I will need some additional work to
> make
> it run faster in parallel, but that's another story). However, I would
> rather like to be able to put these options in my code itself, as defaults.
>

The problems below are exactly why we do not advocate this strategy.
However, you
can get this effect using
http://www.mcs.anl.gov/petsc/petsc-dev/docs/manualpages/Sys/PetscOptionsSetValue.html

   Matt


> So far, I was able to do so for all command line options except for
> -mat_mumps_icntl_13 1
>
> How should I do that? Here goes a snippet of the relevant code:
>
> (... omit creation/population of H matrix ...)
>      ierr = EPSCreate(mixpars->slepc_comm,**&eps);CHKERRQ(ierr);
>      ierr = EPSSetOperators(eps,H,PETSC_**NULL);CHKERRQ(ierr);
>      ierr = EPSSetProblemType(eps,EPS_HEP)**;CHKERRQ(ierr);
>      ierr = EPSSetTolerances(eps, tol, PETSC_DECIDE);
>      if (statesinlist>500)
>        ierr = EPSSetDimensions(eps, 500, PETSC_IGNORE, PETSC_IGNORE);
>      else
>        ierr = EPSSetDimensions(eps, statesinlist, PETSC_IGNORE,
>                                PETSC_IGNORE);
>      CHKERRQ(ierr);
>      if (statesinlist<50) ierr = EPSSetType(eps, EPSLAPACK);
>      else {
>        ierr = EPSGetST(eps,&st);CHKERRQ(**ierr);
>        ierr = STSetType(st,STSINVERT);**CHKERRQ(ierr);
>        ierr = STGetKSP(st,&ksp);CHKERRQ(**ierr);
>        ierr = KSPSetType(ksp,KSPPREONLY);
>        ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
>        ierr = PCSetType(pc, PCCHOLESKY);CHKERRQ(ierr);
>        if (mixpars->slepc_size>1) {
>          ierr = PCFactorSetMatSolverPackage(**pc,MATSOLVERMUMPS);CHKERRQ(*
> *ierr);
>          /*
>          Mat F;
>          ierr = PCFactorGetMatrix(pc,&F);**CHKERRQ(ierr);
>          PetscInt icntl = 13;
>          PetscInt ival = 1;
>          ierr = MatMumpsSetIcntl(F,icntl,ival)**; CHKERRQ(ierr);
>          */
>        }
>      }
>      ierr = EPSSetFromOptions(eps); CHKERRQ(ierr);
>      ierr = EPSSetInterval(eps,
>                            (PetscScalar) currsubclass->lowtargetlimit[**
> j],
>                            (PetscScalar) currsubclass->uptargetlimit[j]**
> );
>      CHKERRQ(ierr);
>      ierr = EPSSetWhichEigenpairs(eps, EPS_ALL); CHKERRQ(ierr);
>      ierr = EPSSolve(eps);CHKERRQ(ierr);
>
> The part which is commented out is what I tried to use to try to obtain the
> same effect as the -mat_mumps_icntl_13 1 command line option. However, if I
> decomment it, I get the following error:
>
> [0]PETSC ERROR: --------------------- Error Message
> ------------------------------**------
> [0]PETSC ERROR: Operation done in wrong order!
> [0]PETSC ERROR: Matrix not yet factored; call after KSPSetUp() or
> PCSetUp()!
> [0]PETSC ERROR:
> ------------------------------**------------------------------**
> ------------
> [0]PETSC ERROR: Petsc Release Version 3.2.0, Patch 5, Sat Oct 29 13:45:54
> CDT 2011 [0]PETSC ERROR: See docs/changes/index.html for recent updates.
> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
> [0]PETSC ERROR: See docs/index.html for manual pages.
>
> The code works fine if I comment out the lines from Mat F;
> to
> ierr = MatMumpsSetIcntl(F,icntl,ival)**; CHKERRQ(ierr);
> and use the command line option -mat_mumps_icntl_13 1
>
> Is it indeed possible to set this as a default in the code, and if yes, is
> it documented somewhere? Can someone give me a hint?
>
> A side issue: should I explicitly free the PC, KSP and ST I get using the
> EPSGetST(), STGetKSP() and KSPGetPC() calls respectively? I only use them to
> set up default options (instead of giving them on the command line). Or can
> I trust that EPSDestroy() will do it?
>
> Thanks in advance, bye
> Giacomo Mulas
>
> --
> ______________________________**______________________________**_____
>
> Giacomo Mulas <gmulas at oa-cagliari.inaf.it>
> ______________________________**______________________________**_____
>
> OSSERVATORIO ASTRONOMICO DI CAGLIARI
> Str. 54, Loc. Poggio dei Pini * 09012 Capoterra (CA)
>
> Tel. (OAC): +39 070 71180 248     Fax : +39 070 71180 222
> Tel. (UNICA): +39 070 675 4916
> ______________________________**______________________________**_____
>
> "When the storms are raging around you, stay right where you are"
>                         (Freddy Mercury)
> ______________________________**______________________________**_____
>



-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120514/d5b10427/attachment.htm>