Jed, Barry, thanks for your answers. I've tried to implement Jed's solution as it was looking faster to try. This is what I came up with for two different orders, FENOrder and FWNOrder respectively
283 KSP ksp; 284 PC pc; 285 286 ierr = KSPCreate(PETSC_COMM_WORLD,&ksp); CHKERRQ(ierr); 287 ierr = KSPSetOptionsPrefix(ksp,"mg_"); CHKERRQ(ierr); 288 ierr = KSPSetOperators(ksp,J,J,SAME_NONZERO_PATTERN); CHKERRQ(ierr); 289 ierr = KSPSetFromOptions(ksp); CHKERRQ(ierr); 290 ierr = KSPGetPC(ksp,&pc); CHKERRQ(ierr); 291 ierr = PCFactorSetFill(pc,1); CHKERRQ(ierr); 292 ierr = PCFactorSetMatOrderingType(pc,"FENOrder"); CHKERRQ(ierr); 293 ierr = KSPSolve(ksp,b,x); CHKERRQ(ierr); 294 ierr = KSPDestroy(&ksp); CHKERRQ(ierr); 295 296 ierr = KSPCreate(PETSC_COMM_WORLD,&ksp); CHKERRQ(ierr); 297 ierr = KSPSetOptionsPrefix(ksp,"mg_"); CHKERRQ(ierr); 298 ierr = KSPSetOperators(ksp,J,J,SAME_NONZERO_PATTERN); CHKERRQ(ierr); 299 ierr = KSPSetFromOptions(ksp); CHKERRQ(ierr); 300 ierr = KSPGetPC(ksp,&pc); CHKERRQ(ierr); 301 ierr = PCFactorSetFill(pc,1); CHKERRQ(ierr); 302 ierr = PCFactorSetMatOrderingType(pc,"FWNOrder"); CHKERRQ(ierr); 303 ierr = KSPSolve(ksp,b,x); CHKERRQ(ierr); 304 ierr = KSPDestroy(&ksp); CHKERRQ(ierr); I don't know if this is what you intended... is there a way to make KSP recompute the factorization with the different ordering without destroying and recreating the KSP? Concerning Barry's option, js it possible to use the PETSc ilu factorization function (and eventually the ones provided by external libraries) inside my PCSHELL, passing it the various ordering I will need? If so, what's the best way to access it? Thanks again Gianluca On 20 June 2011 22:20, Barry Smith <bsmith at mcs.anl.gov> wrote: > > ?I would suggest implementing your own PCSHELL that contains the four > factored matrices and orders. The PCSetUp() you write would compute the > orderings and call the Mat factorization routines. The PCApply() you write > would do the appropriate applications of the triangular solves. The PCILU in > PETSc is just for one ordering and factorization so is not a good starting > point for your algorithm. > > ? Barry > > On Jun 20, 2011, at 8:48 AM, Jed Brown wrote: > >> On Mon, Jun 20, 2011 at 16:26, Gianluca Meneghello <gianmail at gmail.com> >> wrote: >> My goal is to perform some relaxation sweeps for each ordering at each >> level of a multigrid process, probably using ksptype PREONLY and >> pctype ILU. >> Is that possible? >> >> Yes, but you would need to perform a separate factorization for each >> ordering. Also, the orderings that are not aligned with the grid in memory >> will have poor performance. This is a general problem with using different >> orderings. >> >> Is there the equivalent of -pc_sor_its with ILU >> (-pc_ilu_its maybe)? >> >> -ksp_type richardson is the same thing, put it inside -pc_type bjacobi for >> -pc_sor_lits. (-pc_sor_its is just a lower overhead way to do that cycling). >> >> >> I also have the problem that in order to build the ordering I would >> need to have access to a structure containing some grid informations, >> and it seems I cannot pass that structure to the YourOrdering function >> you suggested. >> >> You can PetscObjectCompose() your structure to the Mat. You might need >> PetscContainerCreate() to wrap your struct. >> >> >> I guess a solution for me could be to build the IS from an external >> function (not used by petsc) and then attach them directly to the mat >> structure. I also guess that the one to use are the ones at line 36 of >> >> No, those slots are not public. Use PetscObjectCompose(). > > -- "[Je pense que] l'homme est un monde qui vaut des fois les mondes et que les plus ardentes ambitions sont celles qui ont eu l'orgueil de l'Anonymat" -- Non omnibus, sed mihi et tibi Amedeo Modigliani
