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

Sanjay Govindjee s_g at berkeley.edu
Thu Jan 10 13:18:09 CST 2013


Perhaps more simply, consider the 1D version of your problem:

u''=0 on [0,1]   with "Boundary Conditions"  u'(0)=a and u'(1) =b

The solution is u(x) = a x + c  (for any! value of c).
This is not going to solve happily :(


On 1/10/13 8:27 AM, Jed Brown wrote:
> ... and you describe a "boundary condition" as "second derivative is 
> zero", which is not a boundary condition, making your problem 
> ill-posed. (Indeed, consider the family of problems in which you 
> extend the domain in that patch and apply _any_ boundary conditions in 
> the extended domain. All of those solutions are also solutions of your 
> problem with "second derivative is zero" on your "boundary".)
>
>
> On Thu, Jan 10, 2013 at 9:44 AM, David Scott <d.scott at ed.ac.uk 
> <mailto:d.scott at ed.ac.uk>> wrote:
>
>     All right, I'll say it differently. I wish to solve
>       div.grad phi = 0
>     with the boundary conditions that I have described.
>
>
>     On 10/01/2013 15:28, Jed Brown wrote:
>
>         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
>         <mailto:d.scott at ed.ac.uk>
>         <mailto:d.scott at ed.ac.uk <mailto: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.
>
>
>
>
>     -- 
>     Dr. D. M. Scott
>     Applications Consultant
>     Edinburgh Parallel Computing Centre
>     Tel. 0131 650 5921
>
>
>     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/7d8f204c/attachment.html>


More information about the petsc-users mailing list