# In a PDE solution let us integrate out one of two coordinates

Here is a trapezoid within which I am solving a simple Laplace equation with the Dirichlet boundary conditions set at its top and bottom:

    Needs["NDSolveFEM"];
coord = {{0, 0}, {1/3, 1}, {2/3, 1}, {1, 0}};
incidents = {{1, 2}, {2, 3}, {3, 4}, {4, 1}};
bm = ToBoundaryMesh["Coordinates" -> coord2,
"BoundaryElements" -> {LineElement[incidents]}];
mesh = ToElementMesh[bm];
mesh["Wireframe"]
sigma = 1;
nds = NDSolveValue[{Inactive[Div][
sigma*Inactive[Grad][u[x, y], {x, y}], {x, y}] == 0,
DirichletCondition[u[x, y] == 0, y == 0],
DirichletCondition[u[x, y] == 1, y == 1]},
u[x, y], {x, y} ∈ mesh] Now I would like to calculate the integral of the square of the result gradient: Evaluate[Grad[nds, {x, y}].Grad[nds, {x, y}]] over x lying within the mesh, and study the dependence of the result on y.

How can I do that?

The result should probably be a nested list of with the sublists structure {yi, resi}, where yiis the point coordinate along y, and resi is the integral value in this point.

What I do understand myself, is that I could divide the initial trapezoid into several parts, and then integrate within each such a part and, thus, make a list. I am curious, if there is a more regular approach?

With two minor changes to the code

bm = ToBoundaryMesh["Coordinates" -> coord,
"BoundaryElements" -> {LineElement[incidents]}];
nds = NDSolveValue[{Inactive[Div][sigma*Inactive[Grad][u[x, y], {x, y}], {x, y}] == 0,
DirichletCondition[u[x, y] == 0, y == 0],
DirichletCondition[u[x, y] == 1, y == 1]}, u, {x, y} ∈ mesh];


the first to correct a typo, and the second for convenience, we obtain

Plot3D[nds[x, y], {x, y} ∈ mesh] and

Plot3D[Evaluate[Grad[nds[x, y], {x, y}].Grad[nds[x, y], {x, y}]], {x, y} ∈
mesh, PlotRange -> All] To obtain the integral described in the question, define

f[y0_] := Quiet@NIntegrate[Evaluate[Grad[nds[x, y], {x, y}].Grad[nds[x, y], {x, y}]]
/. y -> y0, {x, y0/3, 1 - y0/3}]


which gives, for instance,

f[.7]
(* 0.672611 *)

• I see the typo you corrected (coord instead of coord2), but cannot find, what was the second correction. Please kindly let me know. – Alexei Boulbitch Aug 5 '15 at 8:18
• I changed u[x. y], {x, y} ∈ mesh to u, {x, y} ∈ mesh at the end of NDSolveValue. – bbgodfrey Aug 5 '15 at 13:28
• Thank you, I appreciate your above answer very much. – Alexei Boulbitch Aug 5 '15 at 13:37