How do you specifiy a Neumann (fixed flux normal to face) boundary condition in Fipy?
<p>How do I explicitly set the flux normal to a boundary face in a fipy mesh to be a specific value, without constraining components of flux within the face?</p>

<p>A Neumann boundary condition can be specified as: (1) fixed component of flux normal to a boundary face, or (2) as a complete specification of flux at the face. The default fipy condition is the former (value = 0), but the explicit method (faceGrad.constrain) is the latter. The issue can be understood by adding the following code to the end of the fipy <a href="http://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.mesh20x20.html" rel="nofollow">diffusion.mesh20x20</a> example and noting different face gradient results. </p>

<pre><code>facesNeumann = mesh.exteriorFaces &amp; ~facesTopLeft &amp; ~facesBottomRight
print 'grad(phi) with default Neumann BC...'
print phi.faceGrad.value.T[facesNeumann.value]
phi.faceGrad.constrain(0.,facesNeumann)
DiffusionTerm().solve(var=phi)
print 'and with explicit Neumann BC...'
print phi.faceGrad.value.T[facesNeumann.value]
</code></pre>

<p>In terms of just specifying the flux normal to a boundary for a variable independent of solving an equation, there does not appear to be any way to do that. The syntax might be <code>phi.faceGrad[0].constrain(...)</code> for example, but that currently doesn't work in FiPy. It might also be difficult to implement for arbitrary oriented faces.</p>

<p>However, for practical purposes, the component tangential to a boundary is not used when solving an equation in FiPy, only the normal component has any influence on the solution. For example, take the following code,</p>

<pre><code>import fipy as fp
mesh = fp.Grid2D(nx=2, ny=1)
var = fp.CellVariable(mesh=mesh)
var.constrain(1, mesh.facesLeft)
var.constrain(0, mesh.facesRight)
#var.faceGrad.constrain(0, mesh.facesTop)
fp.DiffusionTerm().solve(var)
print 'face gradient on top plane:',var.faceGrad[0, mesh.facesTop.value]
print 'variable value:',var
</code></pre>

<p>This gives an output of</p>

<pre><code>face gradient on top plane: [-0.5 -0.5]
variable value: [ 0.75  0.25]
</code></pre>

<p>The answer is correct, but the top face gradient is -0.5. However, when the line <code>#var.faceGrad.constrain(0, mesh.facesTop)</code> is uncommented, the result is</p>

<pre><code>face gradient on top plane: [ 0.  0.]
variable value: [ 0.75  0.25]
</code></pre>

<p>The tangential face gradient is now 0, but the answer is the same. The point is that setting the face gradient tangentially (via <code>.constrain</code>) has no influence on the solution.</p>

python mesh pde fipy