[petsc-users] a strange bug on using DMDA_STENCIL_BOX
Sun, Hui
hus003 at ucsd.edu
Thu Jul 10 23:01:11 CDT 2014
Thank you very much for your time, Barry. I really appreciate it!
Best,
Hui
________________________________________
From: Barry Smith [bsmith at mcs.anl.gov]
Sent: Thursday, July 10, 2014 8:41 PM
To: Sun, Hui
Cc: petsc-users at mcs.anl.gov
Subject: Re: [petsc-users] a strange bug on using DMDA_STENCIL_BOX
On Jul 10, 2014, at 10:05 PM, Sun, Hui <hus003 at ucsd.edu> wrote:
> If you run this:
> ./test_29_cp -da_grid_x 41 -da_grid_y 41 -ksp_type bcgs -ksp_monitor_short -ksp_converged_reason
>
> Then DMDA_STENTIL_STAR gives you convergence to the correct solution in 45 steps, while DMDA_STENTIL_BOX does not converge to the correct solution even after 1000 steps.
I saved the matrices and right hand sides for both of these runs (box and star). They are exactly the same. So you are doing everything right.
I also ran both box and star with -pc_type jacobi and -ksp_type gmres and got essentially the exact same convergence (as one would expect in exact arithmetic).
The combination of bcgs and ilu(0) is just messing you up for box stencil. Intuitively one would expect that putting extra zeros into the matrix would result in a “better” ilu(0) preconditioner, but it need not be the case!
I would recommend you avoid ilu as the preconditioner for your problem. Likely you should use PCFIELDSPLIT
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>
>
More information about the petsc-users
mailing list