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>
