<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><span class=""><br>
<br>
</span>  Yeah, gamg tries to put all the rows of the matrix on one process and leave the other processes empty, </blockquote><div><br></div><div>Actually, GAMG does not do that anymore.  But there are two parameters that govern this and the defaults result in one active process on the coarse grid:  -pc_gamg_process_eq_limit <a> and -pc_gamg_coarse_eq_limit <b> .</div><div><br></div><div>Basically if b < a then you get one active process on the coarse grid.  If you are using SuperLU then you have a large coarse grid (b).  </div><div><br></div><div>You can also set the number of MG levels to get a large coarse grid.</div><div><br></div><div>Mark</div><div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">I suspect superlu_dist doesn't like empty matrices. Just use -mg_coarse_mg_coarse_sub_pc_type lu it is PETSc's native sequential sparse solver and will run faster than superlu_dist anyways except for really large coarse problems which you don't want anyway.<br>
<br>
  How many iterations is the multigrid using?<br>
<br>
  Instead of using -mg_coarse_pc_type some of the PETSc guys are fans of -mg_coarse_pc_type  chebyshev which may decrease the iterations slightly.<br>
<br>
  You can run with -log_summary and send the output, this shows where the time is being spent.<br>
<br>
  Barry<br>
<span class=""><br>
<br>
><br>
> but I got an error:<br>
><br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> ****** Error in MC64A/AD. INFO(1) = -2<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
> symbfact() error returns 0<br>
><br>
><br>
> Thank you,<br>
> Michele<br>
><br>
><br>
> On Thu, 2015-07-16 at 18:07 -0500, Barry Smith wrote:<br>
</span><span class="">>> > On Jul 16, 2015, at 5:42 PM, Michele Rosso <<a href="mailto:mrosso@uci.edu">mrosso@uci.edu</a>> wrote:<br>
>> ><br>
>> > Barry,<br>
>> ><br>
>> > thanks for your reply. So if I want it fixed, I will have to use the master branch, correct?<br>
>><br>
>><br>
>>   Yes, or edit mg.c and remove the offending lines of code (easy enough).<br>
>><br>
>> ><br>
>> > On a side note, what I am trying to achieve is to be able to use how many levels of MG I want, despite the limitation imposed by the local number of grid nodes.<br>
>><br>
>><br>
>>    I assume you are talking about with DMDA? There is no generic limitation for PETSc's multigrid, it is only with the way the DMDA code figures out the interpolation that causes a restriction.<br>
>><br>
>><br>
>> > So far I am using a borrowed code that implements a PC that creates a sub communicator and perform MG on it.<br>
>> > While reading the documentation I found out that PCMGSetLevels takes in an optional array of communicators. How does this work?<br>
>><br>
>><br>
>>    It doesn't work. It was an idea that never got pursued.<br>
>><br>
>><br>
>> > Can I can simply define my matrix and rhs on the fine grid as I would do normally ( I do not use kspsetoperators and kspsetrhs ) and KSP would take care of it by using the correct communicator for each level?<br>
>><br>
>><br>
>>    No.<br>
>><br>
>>    You can use the PCMG geometric multigrid with DMDA for as many levels as it works and then use PCGAMG as the coarse grid solver. PCGAMG automatically uses fewer processes for the coarse level matrices and vectors. You could do this all from the command line without writing code.<br>
>><br>
>>    For example if your code uses a DMDA and calls KSPSetDM() use for example -da_refine 3 -pc_type mg -pc_mg_galerkin -mg_coarse_pc_type gamg  -ksp_view<br>
>><br>
>><br>
>><br>
>>   Barry<br>
>><br>
>><br>
>><br>
>> ><br>
>> > Thanks,<br>
>> > Michele<br>
>> ><br>
>> ><br>
>> ><br>
>> ><br>
>> > On Thu, 2015-07-16 at 17:30 -0500, Barry Smith wrote:<br>
>> >>    Michel,<br>
>> >><br>
>> >>     This is a very annoying feature that has been fixed in master<br>
>> >> <a href="http://www.mcs.anl.gov/petsc/developers/index.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/developers/index.html</a><br>
</span><span class="">>> >>   I would like to have changed it in maint but Jed would have a shit-fit :-) since it changes behavior.<br>
>> >><br>
>> >>   Barry<br>
>> >><br>
>> >><br>
</span><span class="">>> >> > On Jul 16, 2015, at 4:53 PM, Michele Rosso <<a href="mailto:mrosso@uci.edu">mrosso@uci.edu</a>> wrote:<br>
>> >> ><br>
>> >> > Hi,<br>
>> >> ><br>
>> >> > I am performing a series of solves inside a loop. The matrix for each solve changes but not enough to justify a rebuilt of the PC at each solve.<br>
>> >> > Therefore I am using  KSPSetReusePreconditioner to avoid rebuilding unless necessary. The solver is CG + MG with a custom  PC at the coarse level.<br>
>> >> > If KSP is not updated each time, everything works as it is supposed to.<br>
>> >> > When instead I allow the default PETSc  behavior, i.e. updating PC every time the matrix changes, the coarse level KSP , initially set to PREONLY, is changed into GMRES<br>
>> >> > after the first solve. I am not sure where the problem lies (my PC or PETSc), so I would like to have your opinion on this.<br>
>> >> > I attached the ksp_view for the 2 successive solve and the options stack.<br>
>> >> ><br>
>> >> > Thanks for your help,<br>
>> >> > Michel<br>
>> >> ><br>
>> >> ><br>
>> >> ><br>
>> >> > <ksp_view.txt><petsc_options.txt><br>
>> >><br>
>> >><br>
>> >><br>
>> ><br>
>><br>
>><br>
>><br>
><br>
</span>> <ksp_view.txt><br>
<br>
</blockquote></div><br></div></div>