How to calculate the maximum value of numerical solutions u_h

General
1

Hi, everybody！

I want to know how to calculate the maximum value of a numerical solution of dicontinuous Galerkin method？

Thanks,
Yongbin

2

Hi Yongbin,

In general computing the exact maximum of a (piecewise) polynomial is quite hard, however you can approximate the maximum as follows: you sample your function on a (sufficiently large) set of points, and then you take the maximum of the sampled values.

For example, in two dimensions, I can sample the solution on the vertices and on midpoints of edges. The following function generate a list of points that can be used for evaluation:

def computeQuadPoints(mesh, fes):
# create integration points for fast evaluation
eps = 1e-4
ir = IntegrationRule(points=[(0+eps,0+eps),(0+eps,1-2eps), (1-2eps,0+eps), (0.5,0+eps), (0.5-eps, 0.5- eps), (0+eps,0.5)], weights=[1/6,1/6,1/6,1/6,1/6,1/6])
xpts =
ypts =
zpts =
for e in fes.Elements():
trafo = e.GetTrafo()
for p in ir :
xpts.append(trafo(p).point[0])
ypts.append(trafo(p).point[1])
zpts.append(0)
mips = mesh(xpts, ypts, zpts)
npts = 6
return mips, npts

It is inspired from a code by Gousheng Fu. Then, you can evaluate gfu on mips and take the maximum:

vals = gfu(mips)
maxu = np.max(vals)

I hope this solves your problem.

Enrico

3

there is a function to get all the mapped integration points from a integration rule on a mesh, so after the intrule line you can also do this one liner:

ir = ...
return mesh.MapToAllElements(ir, VOL)

instead of VOL you can also give a Region to only find maxima on certain regions.

4

Thank you very much for your advice.
Yongbin Han

5

Thank you very much!
Yongbin Han