# NDSolveValue for Laplace equation not converging to analytic solution

I'm solving Laplace equation $$\nabla^2 \phi = 0$$ with BC's $$\phi_x(x=\pm 1) = 0,\, \phi_y(y=-h) = 0$$ with a specified BC along the circular arc $$x^2 + (-1 + y)^2 = 4$$, which I call $$\Gamma$$ (so the domain $$\Omega$$ is like a rectangular box with the circular arc $$\Gamma$$ as the lid).

What I've done is taken the analytic solution $$\phi_A = \cos(\pi m (x + 1)/2) \cosh(\pi m (y + h)/2) : m \in \mathbb{N}$$ and evaluated its normal derivative on $$\Gamma$$, and then set this as a Neumann boundary condition for the top lid. Since this does not give a unique solution (no Dirichlet BC) I also give the PDE the BC $$\phi = \phi_A(1,\Gamma(x=1))$$. However, the numerical solution does not converge to the analytic solution. Code is supplied below, which shows the domain $$\Omega$$ for clarity.

h = 2;
r = 2;
yp = 1;
m = 1;
Ω = ImplicitRegion[-1 <= x <= 1 && -h <= y <= 2 && ! (x^2 + (y - yp)^2 <= r^2), {x, y}];
Region[Ω]
ϕA = Cos[(π m)/2 (x + 1)] Cosh[(π m)/2 (y + h)];
g = ϕA /. y -> -Sqrt[r^2 - x^2] + yp /. x -> 1;
dnϕA = Grad[ϕA, {x, y}]. ({x, y}/Sqrt[x^2 + y^2]);
dnϕΓA = dnϕA /. y -> -Sqrt[r^2 - x^2] + yp;
op = Laplacian[ϕD[x, y], {x, y}];
ΓNV = NeumannValue[dnϕA, x^2 + (y - yp)^2 == r^2];
Γ = {DirichletCondition[ϕD[x, y] == g, x == 1 && y == -Sqrt[r^2 - 1^2] + yp]};
Needs["NDSolveFEM"]
mesh = ToElementMesh[Ω, MaxCellMeasure -> 0.01];
mesh["Wireframe"]
ϕ = NDSolveValue[{op == ΓNV, Γ},
ϕD, {x, y} ∈ mesh, InterpolationOrder -> 3];
Plot3D[{ϕ[x, y], ϕA}, {x, y} ∈ Ω,
AxesLabel -> {"x", "y", "ϕ"}]

• I think your Neumann value is wrong: since the circle lid is not centered around (0,0), the unit normal in dn\[Phi]A is not (x,y)/norm(...) but (x,y-1)/norm(...), so dn\[Phi]A = Grad[\[Phi]A, {x, y}].(Normalize[{x, y - 1}]); Jan 14, 2019 at 22:22
• @egwenesedai More precisely dn\[Phi]A = Grad[\[Phi]A, {x, y}].({x, y - yp}/Sqrt[x^2 + (y - yp)^2]); Jan 15, 2019 at 13:50
• Thank you egwene sedai and Alex Trounev. I'm new here and unsure how to upvote your comments, but I agree with you both. (Also, a little unclear if I should be thanking each of you...) Jan 15, 2019 at 21:13
• @JoshMcCraney: something that AFAIK is very welcome on this site in cases where noone answers but people provide useful information in comments is to answer your own question with the information given by others in comments (which also is a good place to mention their help/input)... Jan 15, 2019 at 22:28

Edit: my normal derivative to the circular arc $$\Gamma$$ assumed the center was the origin. After recognizing the center was intentionally above the origin, the normal vector should be changed to
({x, y - yp}/Sqrt[x^2 + (y - yp)^2])