[petsc-users] a strange bug on using DMDA_STENCIL_BOX

[Barry Smith]

On Jul 10, 2014, at 10:07 PM, Sun, Hui <hus003 at ucsd.edu> wrote:

> I have also noticed those unnecessary zeros, in both cases DMDA_STENTIL_BOX and DMDA_STENTIL_STAR there are many unnecessary zeros. Do you have any idea how to get rid of those entries? 

   These zeros are for all the locations allowed by the stencil (assuming a full stencil for every variable). The only way to get rid of them is to use a custom matrix allocation.

   The box stencil will have many extra of these zeros since it includes a lot of “coupling” in the matrix that doesn’t exist.

  Barry

> 
> 
> ________________________________________
> From: Barry Smith [bsmith at mcs.anl.gov]
> Sent: Thursday, July 10, 2014 7:46 PM
> To: Sun, Hui
> Cc: petsc-users at mcs.anl.gov
> Subject: Re: [petsc-users] a strange bug on using DMDA_STENCIL_BOX
> 
>   Hmm. Because the box stencil puts extra zeros in the original matrix (due to the “extra” stencil locations) ilu(0) will converge differently with the box stencil than the star stencil with the linear solver.  Thus one will actually get a different convergence history. With Jacobi (which doesn’t depend on the non-zero structure of the matrix) the convergence won’t change. I checked with Jacobi and -da_refine 1 and they converge the same. Are you sure there is an issue with box stencil? You wrote originally
> 
>>>> However, if I change it to DMDA_STENTIL_BOX, it does not converge, the residue \|Ax-b\|_2 remains large.
> 
> can you explain to me how to reproduce this “does not converge, the residue \|Ax-b\|_2 remains large.”  business.
> 
>   Barry
> 
> On Jul 10, 2014, at 8:46 PM, Sun, Hui <hus003 at ucsd.edu> wrote:
> 
>> Thanks for the hint. That is a possible bug. However, I still haven't caught the bug after rethinking about it. I've got three degrees of freedom on each node, MAC grid, and std 5 points laplacian, central differencing for derivatives. Here I've also attached my code.
>> 
>> 
>> ________________________________________
>> From: Barry Smith [bsmith at mcs.anl.gov]
>> Sent: Thursday, July 10, 2014 6:14 PM
>> To: Sun, Hui
>> Cc: petsc-users at mcs.anl.gov
>> Subject: Re: [petsc-users] a strange bug on using DMDA_STENCIL_BOX
>> 
>> Perhaps the code is not correct for the star stencil? That is perhaps you are actually setting values in the “corner” locations that are not being communicated with stencil but are communicated with star?  If your finite difference stencil involves the “corner” locations then you need to always use the box stencil type.
>> 
>>  What is exactly your finite difference stencil? How many degrees of freedom per point?
>> 
>> 
>> Barry
>> 
>> On Jul 10, 2014, at 8:07 PM, Sun, Hui <hus003 at ucsd.edu> wrote:
>> 
>>> If I run:
>>> ./test_29 -da_grid_x 6 -da_grid_y 6 -ksp_type bcgs -ksp_monitor_short -ksp_converged_reason
>>> 
>>> For DMDA_STENCIL_BOX, I get:
>>> 
>>> 0 KSP Residual norm 17.5431
>>> 1 KSP Residual norm 3.23193
>>> 2 KSP Residual norm 0.399307
>>> 3 KSP Residual norm 0.0103842
>>> 4 KSP Residual norm 0.000124391
>>> Linear solve converged due to CONVERGED_RTOL iterations 4
>>> 
>>> For DMDA_STENCIL_STAR, I get:
>>> 
>>> 0 KSP Residual norm 60.7996
>>> 1 KSP Residual norm 13.7279
>>> 2 KSP Residual norm 4.15583
>>> 3 KSP Residual norm 0.430568
>>> 4 KSP Residual norm 0.0999957
>>> 5 KSP Residual norm 0.00260724
>>> 6 KSP Residual norm 0.00191392
>>> 7 KSP Residual norm 5.4681e-05
>>> Linear solve converged due to CONVERGED_RTOL iterations 7
>>> 
>>> The history seems different.
>>> 
>>> 
>>> ________________________________________
>>> From: Barry Smith [bsmith at mcs.anl.gov]
>>> Sent: Thursday, July 10, 2014 5:54 PM
>>> To: Sun, Hui
>>> Cc: petsc-users at mcs.anl.gov
>>> Subject: Re: [petsc-users] a strange bug on using DMDA_STENCIL_BOX
>>> 
>>> Your interpretation is correct. If DMDA_STENTIL_STAR works then one absolutely expects that DMDA_STENTIL_BOX would produce the same answer.  If you are running with linux could you run with valgrind to see if there is some strange memory corruption or lack of initialization http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind
>>> 
>>> What happens with -ksp_monitor on both cases? Same convergence history or different?
>>> 
>>> Barry
>>> 
>>> 
>>> 
>>> On Jul 10, 2014, at 7:49 PM, Sun, Hui <hus003 at ucsd.edu> wrote:
>>> 
>>>> I'm using KSP for 2D Navier Stokes, finite difference. I have DM object created with DMDA_STENTIL_STAR. I'm only using one single process to run the code. It works fine, the linear system converges to the correct solution with correct order, as expected. However, if I change it to DMDA_STENTIL_BOX, it does not converge, the residue \|Ax-b\|_2 remains large.
>>>> 
>>>> I have checked that the rhs for both cases are exactly the same, and the matrices also seem to be the same because I have checked the output of A*y, where y is the analytic solution, and both cases give exactly the same answer. If the linear system is exactly the same, then it seems very strange to me that DMDA_STENTIL_STAR works while DMDA_STENTIL_BOX does not work, for I use the same solver, there shouldn't be any difference.
>>>> 
>>>> If you happen to encounter this kind of bug before and have some idea what's going on, please let me know. Thank you very much! If you want to look at the code I've written, feel free to ask.
>>> 
>> 
>> <test_29_cp.c>
>