4
$\begingroup$

I have a region r0, obtained from an interpolated function, with the use of RegionPlot. Let's call this function under the conditions to form the region $\phi$. A approximation of $\phi$ can be given as,

ϕ=10 y^(1/2))/x > 1 && (10 y^(1/2))/x < 1.2

I have used the following code to define r0,

r0=DiscretizeGraphics[RegionPlot[ϕ, {x,1,100}, {y,1,100}, PlotPoints -> 100]]

Which looks like this,

RegionPlot[r0]

enter image description here

Now what I'd thought would be the simple part, isn't living up to the expectations. After applying TransformedRegion on the region r0 with the following function f results in a error.

f={#[[2]] , #[[1]]/#[[2]]} &;
TransformedRegion[r0, f]

CoefficientArrays::poly: Region`Private`RFVar$11/Region`Private`RFVar$12 is not a polynomial. >>
TransformedRegion::rnimpl: The function TransformedRegion is not implemented foTransformedRegion::rnimpl: The function TransformedRegion is not implemented fo

While a 'pure' region, with $\phi$, works perfectly fine.

r1= ImplicitRegion[(10 y^(1/2))/x > 1 && (10 y^(1/2))/x < 1.2 && x > 1 &&x < 100 && y > 1 && y < 100, {x, y}];
f = {#[[2]] , #[[1]]/#[[2]]} &;
RegionPlot[DiscretizeRegion[TransformedRegion[r1, f]]]

enter image description here


So my question, is there a way of transforming an arbitrarily obtained region, other than using TransformedRegion since it does not seem to agree with the DiscretizeGraphics/RegionPlot combination? In my particular problem I am using a map $f:(x,y)\rightarrow(\alpha y,x/y)$ with $\alpha\in\mathbb{R}$, which is only surjective (or onto) and not injective (or one-to-one).

Furthermore it seems that TransformedRegion is broken in combination with RegionPlot according to this link, but I am not sure if this is the case with this example.

$\endgroup$

1 Answer 1

10
$\begingroup$

The OP's updated example

The OP's example exhibits some numerical problems about which the fastidious ToElementMesh and even some System functions complain. Since the OP is dealing with the System` Region* functions to produce graphics, I'll assume the warnings can be ignored as long as the functions do not fail. There are two things that lead to problems in the OP's example. (1) The RegionPlot is so fine that somehow it generates parts that lead to error messages. (2) The transformation reverses the orientation and in some cases the order-2 "triangles" (which are really hexagons from a Graphics point of view) self-intersect due to the curvature of the images of the lines under the transformation.

A simpler approach here is to work with the boundary of the region and manually transform the coordinates. We can also use a linear ElementMesh since the interest seems to be graphics and not PDEs.

ϕ = (10 y^(1/2))/x > 1 && (10 y^(1/2))/x < 1.2;
reg = DiscretizeGraphics@ RegionPlot[ϕ, {x, 1, 100}, {y, 1, 100}, PlotPoints -> 100];
breg = BoundaryMesh[reg];

BoundaryMesh::brepl: There are components in [reg] having dimension lower than the embedding dimension 2 that will not be included in the boundary representation. >>

Needs["NDSolve`FEM`"];
bmesh = ToBoundaryMesh[
   "Coordinates" -> Function[{x, y}, {y, x/y}] @@@ MeshCoordinates[breg],
   "BoundaryElements" -> {LineElement @@ Thread[MeshCells[breg, 1], Line]}
   ];
emesh = ToElementMesh[bmesh, "MeshOrder" -> 1, 
   MaxCellMeasure -> {"Area" -> Infinity}];

tmesh = MeshRegion[emesh];
Show[tmesh, AspectRatio -> 1/2, Frame -> True]

Mathematica graphics

My original example

Without a supplied example, let's use this:

reg = DiscretizeGraphics@ RegionPlot[(x^2 + y^2)^2 - 12 x y <= 1, {x, -2, 2}, {y, -2, 2}]

Mathematica graphics

I get a different error, suggesting that TransformedRegion is not implemented for mesh regions (well, sort of suggesting something like it).

TransformedRegion[reg, Function[{x, y}, {x - y^2/2, y - x^2/2}]];

TransformedRegion::rnimpl: The function TransformedRegion is not implemented for (reg) . >>

Such functionality is built into the FEM utilities, though.

Needs["NDSolve`FEM`"]

emesh = ToElementMesh[reg, "MeshOrder" -> 1];  (* order 1 added in update *)

xfn = Function[{x, y}, x - y^2/2];
yfn = Function[{x, y}, y - x^2/2];
xifn = ElementMeshInterpolation[{emesh}, xfn @@@ emesh["Coordinates"]];
yifn = ElementMeshInterpolation[{emesh}, yfn @@@ emesh["Coordinates"]];
dmesh = ElementMeshDeformation[emesh, {xifn, yifn}];

treg = MeshRegion@dmesh

Mathematica graphics

Notes: (1) Since we're transforming a mesh region, the polygonal boundary will be transformed into a polygonal boundary. ToElementMesh may subdivide an original boundary segment and the subdividing points will be mapped onto the image of the segment; they will not be mapped on to the boundary of the original region defined by the equation. (2) A transformation that is not 1-1 or reverses orientation will generate an error.

Side note on OP's example and ElementMeshDeformation

It turns out that if we apply the original approach to the boundary mesh breg, ToElementMesh complains but produces an ElementMesh that can be translated to a MeshRegion without difficulty. (The problem is that the ordering of the vertices of the triangles is wrong.)

emesh = ToElementMesh[breg, "MeshOrder" -> 1];

xfn = Function[{x, y}, y - x];
yfn = Function[{x, y}, x/y - y];
xifn = ElementMeshInterpolation[{emesh}, xfn @@@ emesh["Coordinates"]];
yifn = ElementMeshInterpolation[{emesh}, yfn @@@ emesh["Coordinates"]];
dmesh = ElementMeshDeformation[emesh, {xifn, yifn}]

ElementMesh::femimq: The element mesh has insufficient quality of -0.545049. A quality estimate below 0. may be caused by a wrong ordering of element incidents or self-intersecting elements. >>

(*  ElementMesh[{{1., 100.}, {0.83334, 9.99999}}, {TriangleElement["<" 12707 ">"]}]  *)

Despite the error, we can still use the ElementMesh to get the desired MeshRegion.

tmesh = MeshRegion@dmesh;
Show[tmesh, Frame -> True, AspectRatio -> 1/2]
(* output similar to above *)
$\endgroup$
6
  • $\begingroup$ First off, thank you for your clear instructions, and it should have worked, but sadly the transformation I am working with is everything but injective. Which you predicted would result in an error, as I can confirm, and an incorrect region (saying "the element mesh has insufficient quality"). Is there anyway to get this working with non-injective maps, such as $(x,y)\rightarrow)(y,x/y)$? $\endgroup$
    – user19218
    Commented Jun 15, 2015 at 8:48
  • $\begingroup$ @user19218 Your actual example has some perniciously thin triangles. See if my updated method works for you. $\endgroup$
    – Michael E2
    Commented Jun 15, 2015 at 15:03
  • 1
    $\begingroup$ That certainly did the trick! Using the boundary of the region is an rather elegant solution I must add, since they do map one-to-one, the original answer could be used, I tip my hat off to you sir. And indeed I neglected to mention my interest was solely in the shape, nothing PDE related. $\endgroup$
    – user19218
    Commented Jun 15, 2015 at 15:39
  • 1
    $\begingroup$ May I use a modified version of the above boundary mesh method in a coding project licensed under MIT-0 (much more permissive than CC BY-SA 3.0 for this post)? I will of course link this answer in the code comments. $\endgroup$
    – yawnoc
    Commented Apr 7, 2021 at 5:27
  • 2
    $\begingroup$ @yawnoc It's fine with me — it would be an honor. Thanks. And thanks for asking. I'm not a licensing expert, and I'm not sure I'm the one who can grant permission. If CC BY-SA 3.0, allows you to do this, then fine; if not, then I don't what to say. If my OK makes any difference, then you have it. $\endgroup$
    – Michael E2
    Commented Apr 7, 2021 at 5:38

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.