[petsc-users] Dealing with 2nd Derivatives on Boundaries

[Jed Brown]

Second derivative is not a boundary condition for Poisson; that is the
equation satisfied in the interior. Unless you are intentionally attempting
to apply a certain kind of outflow boundary condition (i.e., you're NOT
solving Laplace) then there is a problem with your formulation. I suggest
you revisit the continuum problem and establish that it is well-posed
before concerning yourself with implementation details.


On Thu, Jan 10, 2013 at 9:03 AM, David Scott <d.scott at ed.ac.uk> wrote:

> Hello,
>
> I am solving Poisson's equation (actually Laplace's equation in this
> simple test case) on a 3D structured grid. The boundary condition in the
> first dimension is periodic. In the others there are Von Neumann conditions
> except for one surface where the second derivative is zero. I have
> specified DMDA_BOUNDARY_NONE in these two dimensions and deal with the
> boundary conditions by constructing an appropriate matrix. Here is an
> extract from the Fortran code:
>
>               if (j==0) then
>                 ! Von Neumann boundary conditions on y=0 boundary.
>                 v(1) = 1
>                   col(MatStencil_i, 1) = i
>                   col(MatStencil_j, 1) = j
>                   col(MatStencil_k, 1) = k
>                 v(2) = -1
>                   col(MatStencil_i, 2) = i
>                   col(MatStencil_j, 2) = j+1
>                   col(MatStencil_k, 2) = k
>                 call MatSetValuesStencil(B, 1, row, 2, col, v,
> INSERT_VALUES, ierr)
>               else if (j==maxl) then
>                 ! Boundary condition on y=maxl boundary.
>                 v(1) = 1
>                   col(MatStencil_i, 1) = i
>                   col(MatStencil_j, 1) = j
>                   col(MatStencil_k, 1) = k
>                 v(2) = -2
>                   col(MatStencil_i, 2) = i
>                   col(MatStencil_j, 2) = j-1
>                   col(MatStencil_k, 2) = k
>                 v(3) = 1
>                   col(MatStencil_i, 3) = i
>                   col(MatStencil_j, 3) = j-2
>                   col(MatStencil_k, 3) = k
>                 call MatSetValuesStencil(B, 1, row, 3, col, v,
> INSERT_VALUES, ierr)
>               else if (k==0) then
>
>
> Here the second clause deals with the second derivative on the boundary.
>
> In order for this code to work I have to set the stencil width to 2 even
> though 'j-2' refers to an interior, non-halo
> point in the grid. This leads to larger halo swaps than would be required
> if a stencil width of 1 could be used.
>
> Is there a better way to encode the problem?
>
> David
>
> --
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336.
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130110/fbdd53e5/attachment.html>