Suppose I have a set of points $\{(x_1,x_2) : g(x_1,x_2) > f(x_1,x_2)\}$, how would I turn this into a mesh that acts as a region for PDE solvers, or even to use as a plotting region?

I have tried using Implicit Region, but my function $f$ is an interpolating function and so is defined over a Disk region, which is not a square domain.

Hence if I run

omega = ImplicitRegion[g[x1,x2] > f[x1,x2], {{x1, -1, 1}, {x2, -1, 1}}]

and try to plot my region


It doesn't give me anything, and trying to use this omega for use as a plotting domain / PDE domain doesn't work either. However, I can plot my region by running

set = f[x1,x2] > g[x1,x2]
RegionPlot[set, {x1, -1, 1}, {x2, -1, 1}]

I suppose this means RegionPlot is more robust in the sense that it doesn't care if my function $f$ is indeterminate at certain points.

I have also tried to change the values of $x_1,x_2$ in the implicit region definition to where my function $f$ is defined i.e. commands like

omega = ImplicitRegion[g[x1,x2] > f[x1,x2], {{x1,x2} \in Disk[0,0,r]}]

dont seem to work for ImplicitRegion.

Here is a minimal working example:

sol = NDSolve[{-Laplacian[u[x, y], {x, y}] == 1, 
    DirichletCondition[u[x, y] == 0, True]}, 
   u, {x, y} \[Element] Disk[]];
g[x1_, x2_] := x1^2 + x2^2;
f[x1_, x2_] := u[x1, x2] /. sol[[1]];
setfn[x1_, x2_] := g[x1, x2] > f[x1, x2];
omega = ImplicitRegion[setfn[x1, x2], {x1, x2}];
Plot3D[{f[x1, x2], g[x1, x2]}, {x1, x2} \[Element] Disk[]]

My problem is that the Region[omega] part of the code gives me some weird Region[Embedding Dimension 2] output instead of actually plotting the region, whereas

RegionPlot[setfn[x1, x2], {x1, -1, 1}, {x2, -1, 1}]

plots the region, even though $f$ is undefined at points outside the disk.

I want to use such an omega as the domain to solve my PDE over, for example

 NDSolve[{-Laplacian[u[x, y], {x, y}] == 0, 
  DirichletCondition[u[x, y] == Sin[x y], True]}, u, {x, 
   y} \[Element] omega]

Also tried to extend my interpolating function with values that are irrelevant to my problem

f[x1_, x2_] := 
 Piecewise[{{u[x1, x2] /. sol[[1]], {x1, x2} \[Element] Disk[]}, {0, 

so it is defined everywhere. This doesn't seem to change it. Perhaps ImplicitRegion does not like interpolating functions at all?

  • 1
    $\begingroup$ Have you tried with DiscretizeRegion[omega]? $\endgroup$
    – Domen
    Sep 5, 2021 at 12:54
  • 1
    $\begingroup$ It's harder to help if you don't include an MWE. $\endgroup$
    – Michael E2
    Sep 5, 2021 at 13:00
  • $\begingroup$ Updated with a MWE $\endgroup$
    – jcm
    Sep 5, 2021 at 21:59
  • 1
    $\begingroup$ RegionPlot[DiscretizeRegion[omega, MaxCellMeasure -> 0.001], AspectRatio -> Automatic] $\endgroup$
    – cvgmt
    Sep 6, 2021 at 0:59
  • $\begingroup$ cvgmt's solution worked - maybe the maxcellmeasure is an essential option for interpolated functions. Cheers mate!! $\endgroup$
    – jcm
    Sep 6, 2021 at 4:38

1 Answer 1


For a subsequent FEM analysis you can use

mesh = ToElementMesh[omega];

enter image description here

For other region related computations use:


enter image description here


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.