# Fast integration of ElementMeshInterpolation[]-functions

In this simplified example I consider a small FEM-mesh

<< NDSolveFEM
poly = Polygon[{{0, 0}, {3/2, 0}, {3/2, 1}, {1 , 1}, {1, 1/2}, {0,1/2}}]
mesh = ToElementMesh[poly, MaxCellMeasure -> 1/5,"MeshElementType" -> "TriangleElement","MeshOrder" -> 1]
mesh["Wireframe"]


For this mesh, I define, using ElementMeshInterpolation[], some basisfunctions φ

p = mesh ["Coordinates"]
φ = Map[ElementMeshInterpolation[mesh, #] &,IdentityMatrix[Length[p ]]];


Trying to integrate products (matrix) of these functions

φφ =Outer[Times, Map[#[x, y] &, φ ], Map[#[x, y] &, φ ]];
M = NIntegrate[φφ, Element[{x, y}, mesh]]; // AbsoluteTiming
(*~6 seconds *)


I notice a large evaluation time.

My question: How is it possible to speed up integration of expressions formed by ElementMeshInterpolation[]-functions?

Thanks!

\$Version

"13.0.0 for Mac OS X ARM (64-bit) (December 3, 2021)"


Many thanks to @Ulrich Neumann for the comments.

Here we provide three alternatives:

1. Options in NIntegrate

2. Changing the InterpolationOrder

3. Chaning the MaxCellMeasure

In each case we provide timings and comparisons of the numerical integration.

## Final version:

• Performing the code of the OP for comparison

We have

<< NDSolveFEM
poly = Polygon[{{0, 0}, {3/2, 0}, {3/2, 1}, {1, 1}, {1, 1/2}, {0,
1/2}}];
mesh = ToElementMesh[poly, MaxCellMeasure -> 1/5,
"MeshElementType" -> "TriangleElement", "MeshOrder" -> 1];
mesh["Wireframe"]


p = mesh["Coordinates"];
φ =
Map[ElementMeshInterpolation[mesh, #] &,
IdentityMatrix[Length[p]]];
φφ =
Outer[Times, Map[#[x, y] &, φ],
Map[#[x, y] &, φ]];


Ιntegrating

MOP = NIntegrate[φφ,
Element[{x, y}, mesh]]; // AbsoluteTiming


{2.14101, Null}

• A suggestion at the level of NIntegrate

With the mesh of the OP untouched

<< NDSolveFEM
poly = Polygon[{{0, 0}, {3/2, 0}, {3/2, 1}, {1, 1}, {1, 1/2}, {0,
1/2}}];
mesh = ToElementMesh[poly, MaxCellMeasure -> 1/5,
"MeshElementType" -> "TriangleElement", "MeshOrder" -> 1];
mesh["Wireframe"]


the following

p = mesh["Coordinates"]
φ =
Map[ElementMeshInterpolation[mesh, #] &,
IdentityMatrix[Length[p]]];
φφ =
Outer[Times, Map[#[x, y] &, φ],
Map[#[x, y] &, φ]];
M1 = NIntegrate[φφ, Element[{x, y}, mesh],
Method -> {"FiniteElement",
"MeshOptions" -> {"MeshElementType" -> TetrahedronElement,
"MaxCellMeasure" -> 1/5}}]; // AbsoluteTiming


returns

{0.542303, Null}

• Comparing the numerical results

A simple comparison

MOP - M1 // Chop[#, 7 10^-3] &


• A suggestion at the level of ElementMeshInterpolation

With the mesh of the OP untouched

<< NDSolveFEM
poly = Polygon[{{0, 0}, {3/2, 0}, {3/2, 1}, {1, 1}, {1, 1/2}, {0,
1/2}}];
mesh = ToElementMesh[poly, MaxCellMeasure -> 1/5,
"MeshElementType" -> "TriangleElement", "MeshOrder" -> 1];
mesh["Wireframe"]


we do

p = mesh["Coordinates"]
φ =
Map[ElementMeshInterpolation[mesh, #, InterpolationOrder -> 1] &,
IdentityMatrix[Length[p]]];
φφ =
Outer[Times, Map[#[x, y] &, φ],
Map[#[x, y] &, \[CurlyPhi]]];
M2 = NIntegrate[φφ,
Element[{x, y}, mesh]]; // AbsoluteTiming


to get

{2.13627, Null}

which is not a huge improvement.

• Comparing the numerical results

As before

MOP - M2


• A suggestion at the level of MaxCellMeasure

We begin with:

<< NDSolveFEM
poly = Polygon[{{0, 0}, {3/2, 0}, {3/2, 1}, {1, 1}, {1, 1/2}, {0,
1/2}}];
mesh = ToElementMesh[poly, MaxCellMeasure -> 1/4,
"MeshElementType" -> "TriangleElement", "MeshOrder" -> 1];
mesh["Wireframe"]


and proceed as before

p = mesh["Coordinates"]
φ =
Map[ElementMeshInterpolation[mesh, #] &, IdentityMatrix[Length[p]]];

φφ =
Outer[Times, Map[#[x, y] &, φ],
Map[#[x, y] &, φ]];
M3 = NIntegrate[φφ,
Element[{x, y}, mesh]]; // AbsoluteTiming


{1.71267, Null}

• Comparing the numerical results

As before

MOP[[All ;; -2, All ;; -2]] - M3 // Chop[#, 10^-3] &


• Thank you for your answer. Evaluation of your code gives step1: 5.84s (9points) and step2: 4.72s (8points). Roughly spoken the speedup depends, not surprising, on the number of points. Still I'm looking for a speedup which is comparable to FEM realization in Mathematica Mar 24 at 7:12
• @UlrichNeumann thanks for taking the time to check. I am a bit confused that the first integration gives roughly 6s, because I copied exactly what you wrote in the OP. And it is evident that I misunderstood your needs, sorry about that.
– bmf
Mar 24 at 7:14
• @UlrichNeumann sorry for the double-comment, but I forgot to say that if there are no follow-ups here, I feel I should delete my response as it is not well-suited. So just let me know if you have additional comments. Thanks :-)
– bmf
Mar 24 at 7:17
• Thank you for your effort. I checked your and my code and corrected my answer to 6s evaluation time. Mar 24 at 8:31
• Thank you for your helpful suggestions. I didn't know NIntegrate-Method-FiniteElement, seems to factor 10 faster! Mar 24 at 22:30