<div><div dir="auto">If none of the suggestions provided is to your taste, why not just build the preconditioner matrix yourself? Seems your have precise requirements and the relevant info of the individual blocks, so you should be able to construct the preconditioner, either using A (original operator) or directly from the discrete problem.</div></div><div dir="auto"><br></div><div dir="auto"><br></div><div><br><div class="gmail_quote"><div dir="ltr">On Mon, 27 Aug 2018 at 05:33, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Sat, Aug 25, 2018 at 2:30 PM Ali Reza Khaz'ali <<a href="mailto:arkhazali@cc.iut.ac.ir" target="_blank">arkhazali@cc.iut.ac.ir</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Dear Barry and Jed,<br>
<br>
Thanks for your great suggestions. I set the options as Jed proposed and <br>
it worked. However, Block Jacobi preconditioner treats the entire <br>
Jacobian matrix as only one giant block. It is good in small systems, <br>
but in larger ones, it requires a lot of memory and runtime, as direct <br>
methods do. In addition, I cannot use -pc_bjacobi_blocks option, <br>
because the block sizes of the Jacobain matrix are variable (In my <br>
sample code, SadeqSize array contain block sizes), and apparently <br>
-pc_bjacobi_blocks option assumes the block sizes are equal.<br>
I have to add an explanation about unusual high tolerances in my code <br>
sample, in fact, I am facing a serious ill-conditioned problem. The <br>
direct methods (like LU) solve the system without problem, but I could <br>
not find any preconditioner+iterative solver that can solve it with the <br>
required precision. Since direct solvers are useless in large systems, <br></blockquote><div><br></div><div>Lets discuss this point a bit further. I assume your system is sparse. Sparse</div><div>direct solvers can solve systems fairly efficiently for hundreds of thousands of</div><div>unknowns. How big do you want? Also, do you plan on having more than 500K</div><div>unknowns per process? If not, why not just use sparse direct solvers on each process?</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div></div></div><div dir="ltr"><div class="gmail_quote"><div><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
I had to design an specific preconditioner for it (I wrote a paper about <br>
it: <a href="https://link.springer.com/article/10.1007/s10596-014-9421-3" rel="noreferrer" target="_blank">https://link.springer.com/article/10.1007/s10596-014-9421-3</a>), but <br>
the designed preconditioner was rather slow (even on GPU), and I hope <br>
that I can find a faster solution in PETSc. So for now I am searching <br>
for a combination that just works for me, and then, I'll refine it.<br>
<br>
Best wishes,<br>
Ali<br>
<br>
On 8/25/2018 4:31 AM, Smith, Barry F. wrote:<br>
> What you would like to do is reasonable but unfortunately the order of the operations in SNES/KSP means that the KSPSetUp() cannot be called before the SNESSolve() because the Jacobian matrix has yet to be provided to KSP. Thus you have to, as Jed suggested, use the options database to set the options. Note that if you want the options "hardwired" into the code you can use PetscOptionsSetValue() from the code.<br>
><br>
> Barry<br>
><br>
><br>
>> On Aug 24, 2018, at 1:48 PM, Ali Reza Khaz'ali <<a href="mailto:arkhazali@cc.iut.ac.ir" target="_blank">arkhazali@cc.iut.ac.ir</a>> wrote:<br>
>><br>
>> Hi,<br>
>><br>
>> I am trying to use block Jacobi preconditioner in SNES (SNESNEWTONLS). However, PCBJacobiGetSubKSP function returns an error stating "Object is in wrong state, Must call KSPSetUp() or PCSetUp() first". When I add KSPSetUp, I got and error from them as: "Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()", and if PCSetUp is added, "Object is in wrong state, Matrix must be set first" error is printed.<br>
>><br>
>> Below is a part of my code. It is run serially. Any help is much appreciated.<br>
>><br>
>> ierr = SNESGetKSP(snes, &Petsc_ksp);<br>
>> CHKERRQ(ierr);<br>
>> ierr = KSPGetPC(Petsc_ksp, &Petsc_pc);<br>
>> CHKERRQ(ierr);<br>
>> ierr = KSPSetTolerances(Petsc_ksp, 1.e-3, 1e-3, PETSC_DEFAULT, 20000);<br>
>> CHKERRQ(ierr);<br>
>> ierr = SNESSetTolerances(snes, 1e-1, 1e-1, 1e-1, 2000, 2000);<br>
>> CHKERRQ(ierr);<br>
>> ierr = SNESSetType(snes, SNESNEWTONLS);<br>
>> CHKERRQ(ierr);<br>
>> ierr = KSPSetType(Petsc_ksp, KSPGMRES);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCSetType(Petsc_pc, PCBJACOBI);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCSetType(Petsc_pc, PCBJACOBI);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCBJacobiSetTotalBlocks(Petsc_pc, 2*Nx*Ny*Nz, SadeqSize);<br>
>> CHKERRQ(ierr);<br>
>><br>
>> SNESSetUp(snes);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCBJacobiGetSubKSP(Petsc_pc, &nLocal, &Firstly, &subKSP);<br>
>> CHKERRQ(ierr);<br>
>><br>
>> for (i = 0; i < nLocal; i++) {<br>
>> ierr = KSPGetPC(subKSP[i], &SubPc);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCSetType(SubPc, PCLU);<br>
>> CHKERRQ(ierr);<br>
>> ierr = PCFactorSetMatSolverPackage(SubPc, "mkl_pardiso");<br>
>> CHKERRQ(ierr);<br>
>> ierr = KSPSetType(subKSP[i], KSPPREONLY);<br>
>> CHKERRQ(ierr);<br>
>> ierr = KSPSetTolerances(subKSP[i], 1.e-6, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT);<br>
>> CHKERRQ(ierr);<br>
>> }<br>
>> ierr = SNESSolve(snes, NULL, Petsc_X);<br>
>> CHKERRQ(ierr);<br>
>><br>
>><br>
>><br>
>> -- <br>
>> Ali Reza Khaz’ali<br>
>> Assistant Professor of Petroleum Engineering,<br>
>> Department of Chemical Engineering<br>
>> Isfahan University of Technology<br>
>> Isfahan, Iran<br>
>><br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="m_2009903902236731568gmail_signature" data-smartmail="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>
</blockquote></div></div>