On Sun, May 3, 2020 at 12:34 AM Zhuo Chen <chenzhu...@gmail.com> wrote:
> Dear Petsc users, > > My name is Zhuo Chen. I am very new to Petsc. > > I have encountered a very weird problem. I try to use Petsc to solve a 2D > diffusion problem (https://en.wikipedia.org/wiki/Heat_equation). 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. > > 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. > The default tolerance for KSP is 1e-5. If you make it lower, you will get what you want, -ksp_rtol 1e-10 However, BJACOBI is a terrible way to solve the diffusion equation. You should use multigrid, for example -pc_type gamg or Hypre/BoomerAMG or ML. Thanks, Matt > I would like to paste a part of my Petsc code here as well. > > 179 call KSPSetOperators(ksp,A,A,ierr);CHKERRQ(ierr) > 180 call KSPGetPC(ksp,pc,ierr);CHKERRQ(ierr) > 181 tol=1d-6 > 182 call > KSPSetTolerances(ksp,tol,PETSC_DEFAULT_REAL,PETSC_DEFAULT_REAL,PETSC_DEFAULT_INTEGER,ierr);CHKERRQ(ierr) > 183 call PCSetType(pc,PCBJACOBI,ierr);CHKERRQ(ierr) > 184 nblks=8 > 185 allocate(blks(nblks)) > 186 blks=Ntot/nblks > 187 call PCBJacobiSetTotalBlocks(pc,nblks,blks,ierr);CHKERRQ(ierr) > 188 deallocate(blks) > 189 call KSPsetup(ksp,ierr) > 190 call PCBJacobiGetSubKSP(pc,nlocal,first,PETSC_NULL_KSP,ierr) > 191 allocate(subksp(nlocal)) > 192 call PCBJacobiGetSubKSP(pc,nlocal,first,subksp,ierr) > 193 do i=0,nlocal-1 > 194 call KSPGetPC(subksp(i+1),subpc,ierr) > 195 call PCSetType(subpc,PCLU,ierr); CHKERRA(ierr) > 196 call KSPSetType(subksp(i+1),KSPGMRES,ierr); CHKERRA(ierr) > 197 tol=1d-6 > 198 call > KSPSetTolerances(subksp(i+1),tol,PETSC_DEFAULT_REAL,PETSC_DEFAULT_REAL,PETSC_DEFAULT_INTEGER,ierr);CHKERRQ(ierr) > 199 end do > 200 deallocate(subksp) > > It would be great if anyone can help me with this issue. Many thanks! > > -- > Zhuo Chen > Department of Physics > University of Alberta > Edmonton Alberta, Canada T6G 2E1 > http://www.pas.rochester.edu/~zchen25/ > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>