<div dir="ltr"><div dir="ltr">On Sun, May 3, 2020 at 12:34 AM Zhuo Chen <<a href="mailto:chenzhuotj@gmail.com">chenzhuotj@gmail.com</a>> wrote:<br></div><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"><div dir="ltr">Dear Petsc users,<div><br></div><div>My name is Zhuo Chen. I am very new to Petsc.</div><div><br></div><div>I have encountered a very weird problem. I try to use Petsc to solve a 2D diffusion problem (<a href="https://en.wikipedia.org/wiki/Heat_equation" target="_blank">https://en.wikipedia.org/wiki/Heat_equation</a>). When I use PCLU, the solution looks very good. Then, I change to PCJACOBI because I plan to use MPI in the future. If I treat the whole matrix as a single block and set the preconditioner of the "sub-matrix" with PCLU, the solution is exactly the same as the PCLU setting. However, if I increase the number of blocks, the result becomes worse. I attached the solution of a 2D diffusion simulation I have done. The red dots are the analytic solution and the blue circles are the solution from Petsc with 8 blocks. The size of the computational domain is 300*8. The result with PCLU only is very close to the analytic solution so I do not include it here.</div><div><br></div><div>Then I implemented the same problem with Matlab with x=gmres(mat,rhs,10,tol,maxit,p), where mat is the same as the one in my Petsc code, and p is the block Jacobi matrix (8 blocks). The result from Matlab is very good, though it is slow. I attach the Matlab result to this email as well.</div></div></blockquote><div><br></div><div>The default tolerance for KSP is 1e-5. If you make it lower, you will get what you want,</div><div><br></div><div> -ksp_rtol 1e-10</div><div><br></div><div>However, BJACOBI is a terrible way to solve the diffusion equation. You should use multigrid, for example</div><div><br></div><div> -pc_type gamg</div><div><br></div><div>or Hypre/BoomerAMG or ML.</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 dir="ltr"><div>I would like to paste a part of my Petsc code here as well.</div><div><br></div><div>179 call KSPSetOperators(ksp,A,A,ierr);CHKERRQ(ierr)<br>180 call KSPGetPC(ksp,pc,ierr);CHKERRQ(ierr)<br>181 tol=1d-6<br>182 call KSPSetTolerances(ksp,tol,PETSC_DEFAULT_REAL,PETSC_DEFAULT_REAL,PETSC_DEFAULT_INTEGER,ierr);CHKERRQ(ierr)<br>183 call PCSetType(pc,PCBJACOBI,ierr);CHKERRQ(ierr)<br>184 nblks=8<br>185 allocate(blks(nblks))<br>186 blks=Ntot/nblks<br>187 call PCBJacobiSetTotalBlocks(pc,nblks,blks,ierr);CHKERRQ(ierr)<br>188 deallocate(blks)<br>189 call KSPsetup(ksp,ierr)<br>190 call PCBJacobiGetSubKSP(pc,nlocal,first,PETSC_NULL_KSP,ierr)<br>191 allocate(subksp(nlocal))<br>192 call PCBJacobiGetSubKSP(pc,nlocal,first,subksp,ierr)<br>193 do i=0,nlocal-1<br>194 call KSPGetPC(subksp(i+1),subpc,ierr)<br>195 call PCSetType(subpc,PCLU,ierr); CHKERRA(ierr)<br>196 call KSPSetType(subksp(i+1),KSPGMRES,ierr); CHKERRA(ierr)<br>197 tol=1d-6<br>198 call KSPSetTolerances(subksp(i+1),tol,PETSC_DEFAULT_REAL,PETSC_DEFAULT_REAL,PETSC_DEFAULT_INTEGER,ierr);CHKERRQ(ierr)<br>199 end do<br>200 deallocate(subksp)<br></div><div><br></div><div>It would be great if anyone can help me with this issue. Many thanks!</div><div><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr"><div dir="ltr"><div dir="ltr"><font color="#0000ff" style="font-family:arial,sans-serif;font-size:12.8px">Zhuo Chen</font><div style="font-family:arial,sans-serif;font-size:12.8px"><font color="#0000ff">Department of Physics</font></div><div style="font-family:arial,sans-serif;font-size:12.8px"><font color="#0000ff">University of Alberta</font></div><div style="font-family:arial,sans-serif;font-size:12.8px"><font color="#0000ff">Edmonton Alberta, Canada T6G 2E1</font></div><div><font face="'comic sans ms', sans-serif"><a href="http://www.pas.rochester.edu/~zchen25/" style="background-color:rgb(255,255,255)" target="_blank"><font color="#000000">http://www.pas.rochester.edu/~zchen25/</font></a><br></font></div></div></div></div></div></div></div></div></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>