<div dir="ltr"><div dir="ltr"><div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Sun, Dec 16, 2018 at 10:12 PM Chukwudi Chukwudozie via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Barry,</div><div><br></div><div>Thanks for the reply. I have attached two files with the information. </div></div><br><div class="gmail_quote"><div dir="ltr">On Sat, Dec 15, 2018 at 8:58 PM Smith, Barry F. <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><br>
Chuks,<br>
<br>
1) Please run a successful case with -snes_view and send the output so we can see the exact solver options being used.<br>
<br>
2) Run a less successful case with -snes_monitor -ksp_monitor_true_residual -ksp_converged_reason -fieldsplit_ksp_monitor_true_residual -fieldsplit_ksp_converged reason so we can see how the linear solver (and the sub linear solves) are (not) converging.<br></blockquote></div></blockquote><div><br></div><div>The first thing I would do is to put Krylov wrappers around the ILU on the subblocks. ILU could be good enough</div><div>for your simple case, but I bet the blocks are not begin accurately solved in the non-convergent case. You could</div><div>turn on monitor for the blocks to confirm this:</div><div><br></div><div> -P_fieldsplit_0_ksp_type gmres</div><div> -P_fieldsplit_0_ksp_monitor_true_residual</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<br>
Barry<br>
<br>
<br>
> On Dec 15, 2018, at 7:54 PM, Chukwudi Chukwudozie via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>> wrote:<br>
> <br>
> Hi,<br>
> <br>
> I am developing a fully implicit control volume finite element black oil reservoir simulator using petscdmplex and snes. So far, for 2 phase flow problems (3 by 3 block for pressure, saturation and bottomhole pressure), my simulator runs fine for small time step sizes in the order of less than 5 days. However, increasing the time step size, it doesn't converge with -3 (-5 sometimes) diverge reason. Considering that this is a fully implicit but nonlinear problem, I expect unconditional stability for reasonable time step sizes.<br>
> <br>
> The background DM for the 3 block problem is a composite one with the first two being DMPLES's (bag->plexScalNode for pressure and saturation) and the last one, bag->WellRedun, a redundant DM for the wellbore pressures. See below.<br>
> <br>
> ierr = DMRedundantCreate(PETSC_COMM_WORLD,0,bag->numWells,&bag->WellRedun);CHKERRQ(ierr);<br>
> ierr = DMCompositeCreate(PETSC_COMM_WORLD,&bag->MultiPhasePacker);CHKERRQ(ierr);<br>
> ierr = DMCompositeAddDM(bag->MultiPhasePacker,bag->plexScalNode);CHKERRQ(ierr);<br>
> ierr = DMCompositeAddDM(bag->MultiPhasePacker,bag->plexScalNode);CHKERRQ(ierr);<br>
> ierr = DMCompositeAddDM(bag->MultiPhasePacker,bag->WellRedun);CHKERRQ(ierr); <br>
> <br>
> I have used nested matrices as shown below. In addition, I used the fieldsplit (composite multiplicative) preconditioner with FGRES and least square for newton solver. See below.<br>
> <br>
> ierr = MatCreateNest(PETSC_COMM_WORLD,Msize,is,Msize,is,&bK[0],K);CHKERRQ(ierr);<br>
> <br>
> <br>
> ierr = SNESGetKSP(bag->snesP,&kspP);CHKERRQ(ierr);<br>
> ierr = KSPSetTolerances(kspP,1.e-8,1.e-8,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr);<br>
> ierr = KSPSetType(kspP,KSPFGMRES);CHKERRQ(ierr);<br>
> ierr = KSPSetFromOptions(kspP);CHKERRQ(ierr);<br>
> ierr = KSPGetPC(kspP,&pcP);CHKERRQ(ierr);<br>
> ierr = PCSetType(pcP, PCFIELDSPLIT);CHKERRQ(ierr);<br>
> ierr = DMCompositeGetGlobalISs(bag->MultiPhasePacker,&is);CHKERRQ(ierr); <br>
> for(i = 0; i < Msize; i++) ierr = PCFieldSplitSetIS(pcP,NULL,is[i]);CHKERRQ(ierr);<br>
> ierr = PCFieldSplitSetType(pcP,PC_COMPOSITE_MULTIPLICATIVE);CHKERRQ(ierr);<br>
> <br>
> Any ideas on possible precoditioner options/solver options/combinations to improve stability will be appreciated.<br>
> <br>
> Chuks<br>
> <br>
<br>
</blockquote></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>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</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div></div></div>