On 02/22/2012 05:13 PM, Jed Brown wrote:

On Wed, Feb 22, 2012 at 18:05, Patrick Alken <[email protected] <mailto:[email protected]>> wrote:

Hi all,

I have been trying to track down a problem for a few days with solving a linear system arising from a finite differenced PDE in spherical coordinates. I found that PETSc managed to converge to a nice solution for my matrix at small grid sizes and everything looks pretty good.

But when I try larger more realistic grid sizes, PETSc fails to converge. After trying with another direct solver library, I found that the direct solver found a solution which exactly solves the matrix equation,

This never happens, so what do you mean? You compute the residual and it's similar to what you expect the rounding error to be?

Yes I mean the direct solver residual is around 10e-15. The PETSc residual is 4e00

but when plotting the solution, I see that it oscillates rapidly between the grid points and therefore isn't a satisfactory solution. (At smaller grids the solution is nice and smooth)

What sort of PDE are you solving?

The PDE is: grad(f) . B = g where B is a known vector field, g is a known scalar function, and f is the unknown scalar function to be determined (I am discretizing this equation for f in spherical coords)

I was wondering if this phenomenon is common in PDEs? and if there is any way to correct for it?

I am currently using 2nd order centered differences for interior grid points, and 1st order forward/backward differences for edge points. Would it be worthwhile to try moving to 4th order differences instead? Or would that make the problem worse?

I've even tried smoothing the parameters which go into the matrix entries using moving averages...which doesn't seem to help too much.

Any advice from those who have experience with this phenomenon would be greatly appreciated!

Thanks, Patrick

Along the same lines, also have a look at Essentially Non-Oscillatory (ENO) and Weighted ENO schemes which basically constructed oscillation-free, high-order, 'upwind' approximations to grad(f). Any standard text on FD (such as the one by LeVeque or the one by Strikwerda) have chapters on why central finite differences fail and why you would need an upwind-based method. Mohammad On Wed, Feb 22, 2012 at 4:24 PM, Jed Brown <[email protected]> wrote:

On Wed, Feb 22, 2012 at 18:15, Patrick Alken <[email protected]>wrote:

...

Yes I mean the direct solver residual is around 10e-15. The PETSc residual is 4e00

Did you try -pc_type lu?

...

What sort of PDE are you solving?

The PDE is:

grad(f) . B = g

where B is a known vector field, g is a known scalar function, and f is the unknown scalar function to be determined (I am discretizing this equation for f in spherical coords)

Also known as steady-state advection. This is the most famous test problem in which centered differences fails miserably. If you use an upwind method, it will be stable, but will also be first order. To get higher than first order accuracy, you need a nonlinear spatial discretization (see Godunov's Theorem, "total variation diminishing", or "total variation bounded").