# Extract mesh coordinates from a LaplaceSolution Plot

I have made a code for the solution to a particular problem involving the Laplacian, and after plotting in 3D the resultant potential, I have tried to look for a way to extract points from that graph in order to later get the electric field produce by that potential, so I wonder if there is a way to do this, or another way to plot the electric field by getting the following solution to the Laplacian?

## This is the solution

a = 4;
f[x_, y_] := Sin[ArcTan[Abs[y/x]]]
g[x_, y_] := -Sin[ArcTan[Abs[y/x]]]
leqn = Laplacian[u[x, y], {x, y}] == 0;
\[CapitalOmega] =
ImplicitRegion[x^2 + y^2 >= 1, {{x, -a, a}, {y, -a, a}}];
A[x_, y_] := NDSolveValue[{\!$$\*SubsuperscriptBox[\(\[Del]$$, $${x, y}$$, $$2$$]$$u[x, y]$$\) == 0,
DirichletCondition[u[x, y] == f[x, y], y == Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == g[x, y], y == -Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == 0,
x == a || x == -a || y == a || y == -a]},
u, {x, y} \[Element] \[CapitalOmega]]
ContourPlot[A[x, y], {x, y} \[Element] \[CapitalOmega],
PlotRange -> All, FrameLabel -> {x, y}, PlotPoints -> 20,
ColorFunction -> "Rainbow", PlotLegends -> Automatic]
Plot3D[A[x, y], {x, y} \[Element] \[CapitalOmega], PlotRange -> All,
ColorFunction -> "Rainbow", ViewPoint -> {0, 0, \[Infinity]},
PlotLegends -> Automatic]


So I tried using this idea from http://community.wolfram.com/groups/-/m/t/452396, by putting Array[{##, A[##]} &, {10, 10}, {x, y} \[Element] \[CapitalOmega]] Graphics3D[Point /@ %]

But it doesn't seem to work, or I am using it wrong. Does anyone has an idea of how to solve this problem?

• Does Array[{##, A[##]} &, {10, 10}, {{-a, a}, {-a, a}}]; Graphics3D[Point /@ %] give what you need? – kglr Oct 25 '17 at 6:31
• It doesn't seem to work, it gives me two errors: InterpolatingFunction::femdmval: Input value {-0.444444,-0.444444} lies outside the range of data in the interpolating function. And, General::stop: Further output of InterpolatingFunction::femdmval will be suppressed during this calculation. – Mounice Oct 25 '17 at 6:34
• Are you sure about that Sin[ArcTan[Abs[y/x]]? It seems to me you should be using two-argument arctangent. – J. M.'s torpor Oct 25 '17 at 7:12

aa = NDSolveValue[{leqn,
DirichletCondition[u[x, y] == f[x, y], y == Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == g[x, y], y == -Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == 0, x == a || x == -a || y == a || y == -a]},
u, Element[{x, y}, Ω]];

points = Array[If[RegionMember[Ω, {##}], {##, aa[##]} , ##&[]]&,
{10, 10}, {{-a, a}, {-a, a}}];
Graphics3D[{Red, PointSize[Large], Point /@ points}, BoxRatios -> 1]


Show[Plot3D[aa[x, y], Element[{x, y}, Ω],  PlotRange -> All,
ColorFunction -> "Rainbow",  PlotLegends -> Automatic], g3d]


• Using RegionMember in this way can be very expensive. You could generate a rmf=RegionMember[\[CapitalOmega]] and use that in the If statement. But using a discretized version of Omega would be even better. – user21 Oct 28 '17 at 0:00
• @user21, thank you. I will update with your suggestion. – kglr Oct 28 '17 at 13:56
• @kglr is there a new command in MMA 12.2 for this task? – ABCDEMMM Jun 7 at 1:01
• @ABCDEMMM, i don't know if there is a new command that does this. – kglr Jun 7 at 9:21

Here is another alternative:

A = NDSolveValue[{Laplacian[u[x, y], {x, y}] == 0,
DirichletCondition[u[x, y] == f[x, y], y == Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == g[x, y], y == -Sqrt[1 - x^2]],
DirichletCondition[u[x, y] == 0,
x == a || x == -a || y == a || y == -a]},
u, {x, y} \[Element] \[CapitalOmega]]


We extract the values and the mesh (and from that the coordinates) from the interpolating function:

vals = A["ValuesOnGrid"];
mesh = A["ElementMesh"];
coords = mesh["Coordinates"];
data = Join[coords, Partition[vals, 1], 2];
Graphics3D[{Red, Point[data]}]