# NDSolve ignores my NeumannValue boundary conditions

I am trying to solve a simple linear differential equation for $$f(x,y)$$ on a square with area $$L\times L =1$$.

I consider $$(\partial_x^2 + \partial_y^2)f + \partial_x \partial_y f = 0$$ with the following boundary conditions: $$\partial_x f(0,y) = 1$$, $$\partial_y f(x,0) = 0$$, $$\partial_y f(x,1) = 0$$, $$f(1,y) = 0$$.

To solve this differential equation, I have used the NDSolveValue method as

sysf = { D[f[x, y], {x, 2}] +  D[f[x, y], {y, 2}] +
D[f[x, y], x, y] ==
NeumannValue[1, x == 0] + NeumannValue[0, y == 0] +
NeumannValue[0, y == 1],
DirichletCondition[f[x, y] == 0, x == 1]};

solf = NDSolveValue[sysf,
f , {x, y} \[Element] Rectangle[{0, 0}, {1, 1}]];

fplot = ContourPlot[
Abs[solf[x, y]], {x, y} \[Element] Rectangle[{0, 0}, {1, 1}],
AxesLabel -> {x, y,  Abs[f]},  PlotLabel -> "f",
ColorFunction -> "BlueGreenYellow", PlotRange -> All,
PlotLegends -> Automatic];
GraphicsRow[{fplot}]


I get no warnings or errors, but the result is clearly incorrect. The gradient normal on the boundaries $$y=0$$ and $$y=1$$ are finite, and in conflict with my enforced boundary conditions.

Apparently, NDSolveValue has ignored my NeumannValue boundary conditions, what could be the reason for this?

• The solution looks correct to me. Maybe you have miss understand what a NeumannValue does; did you have a look at the 'Details' section on the ref page? Also note that the Neumann 0 values at y==0 or y==1 can be left out, as a Neumann 0 is the implicit default value. Dec 6, 2023 at 12:26
• Did you check which NeumannValue belong to your pde-problem? Dec 6, 2023 at 13:11

Your pde problem possesses the following NeumannValue's

{state} =NDSolveProcessEquations[{
D[f[x, y], {x, 2}]  + D[f[x, y], {y, 2}] +D[f[x, y], x, y]  ==NeumannValue[1, x == 0]  ,
DirichletCondition[f[x, y] == 0, x == 1] },
f, {x, y} \[Element] Rectangle[]];

NDSolveFEMGetInactivePDE[state]


$$\left\{\nabla _{\{x,y\}}\cdot \left(-\left( \begin{array}{cc} -1 & -\frac{1}{2} \\ -\frac{1}{2} & -1 \\ \end{array} \right).\nabla _{\{x,y\}}f(x,y)\right)\right\}$$

Last formula (normal brakets) shows the flux. Thus the NeumannValue at x==0 is  Derivative[1, 0][solf][0, y] + 1/2 Derivative[0, 1][solf][0, y]

Plot[  Derivative[1, 0][solf][0, y] + 1/2 Derivative[0, 1][solf][0, y], {y, 0, 1}, PlotRange -> {0.9, 1.1} ]


Plot shows NeumannValue[1,x==0]` is quite good fullfilled.