I would like to solve the following PDE: $$ \Delta \mu (x,y,z)=\mu$$ Over the parallelepiped: $x \in [0,100] $, $y \in [0,100] $ and $z \in [-10,0] $
With the following Neumann boundary conditions: $$ \mathbf n \cdot \nabla \mu |_{z=0}= x y, \ \ \ \ \ \ \ \mathbf n \cdot \nabla \mu |_{ rest \ of \ boundary}= 0$$
I use Mathematica 11.0. First, I specify the region:
reg = Parallelepiped[{0, 0, -10}, {{100, 0, 0}, {0, 100, 0}, {0, 0, 10}}]
Then I solve the diffusion equation:
sol = NDSolveValue[Laplacian[mu[x, y, z], {x, y, z}] - mu[x, y, z] == NeumannValue[- x y, z == 0]},
mu, {x, y, z} \[Element] reg]
I get a solution (however, I get the mistake NDSolveValue::femibcnd: No DirichletCondition or Robin-type NeumannValue was specified for {mu}; the result may not be unique.
). Afterwards, I want to check whether the second boundary condition is satisfied (in documentation it is said that if boundary condition is not specified, it is automatically assumed to be Neumann 0). Therefore, I try to calculate the flow of the function at a random point of the boundary at the lower boundary:
D[sol[1, 1, z], z] /. z -> -10
And obtain -0.00125585, and not 0 as I expected. What is wrong with my solution?
Also, I am confused with the result not being smooth. See the following:
Plot[Evaluate[sol[1, 1, z], {z, -10, 0}]]
Another check shows that probably accuracy of the solution is unsufficient even at the boundary, however, I am not sure how to improve it without significant cost to the speed of computation:
Plot[{Evaluate[sol[1, y, 0]], y}, {y, 0, 100}, PlotRange -> All]
I tried using right-hand side in the following form: NeumannValue[- x y, z == 0] + NeumannValue[0, z != 0]
, which doesn't affect the answer.
NDSolveValue::femibcnd: No DirichletCondition or Robin-type NeumannValue was specified for {mu}; the result may not be unique.
? $\endgroup$NDEigensystem
, since you are trying to find the eigenfunctions of the Laplacian for the eigenvalue 1. But that does not explain the result you observe... $\endgroup$NeumannValue[x y, z == 0] + NeumannValue[0, z != 0]
as right hand side of your pde? $\endgroup$NDSolve
employs the finite element element method: The differential equation gets discretized into a finite dimensional linear system and this gets solved afterwards. If the infinite dimensional equations are not well formulated then the finite dimensional system will be badly behaved. But I have correct myself: Since the differential operator is $\Delta -1$, the system does not have a null space at all. Long story short: I do not know why this happens. I would have expected that Ulrich's answer would settle this problem... $\endgroup$