The boundary conditions, i.e. Eqs. 2, are defined in
subroutine bcond(x,y,ipiece,c,g,itype). The point
is passed in through the real variables x and
y, and the boundary segment number (see section
5.3) is passed in through
the integer variable ipiece. The user defines the
functions 
and 
from Eq. 2 and returns
them through the real variables c and g, and
sets the integer variable itype to be 1, 2 or 3 to
indicate whether the boundary condition is Dirichlet, natural
(Neumann), or mixed (Robin), respectively.
With 
, the second form of the boundary condition
is the `natural' boundary condition for the differential
operator in Eq. 1. Often, the desired boundary condition
is the Neumann condition 
.
In many cases it is possible to represent the Neumann condition
in the form of the natural condition.
If the outward normal to the boundary of the domain makes an
angle 
with the x-axis, then
and
.
The operator for Neumann boundary conditions is
Contrast this to the operator for the natural boundary conditions
(Laplace, Poisson or Helmholtz
equation), then the natural boundary condition is Neumann.
but they are not identically 1,
multiply the boundary condition by p to get the correct form.
For example, with the mixed condition
the correct form for the mixed natural boundary condition is
So the coefficient of u is 
and the right hand
side is 
and in subroutine bcond you set
c = 
and g = 
.
and the sides of the domain are parallel to
the x and y axes, then from the definition of
in terms of 
, the mixed boundary
condition can be expressed as
and the sides are not parallel to the axes,
it is not, in general, possible to express Neumann conditions.
However, in some special cases it might be possible.