Solve a PDE over a region defined by a Bézier patch
You might create a NearestFunction
to help pick the particular boundary you want. You can use it to mark the boundary elements of an ElementMesh
(FEM).
plot = ParametricPlot[
bezierfunc[ξ, η], {ξ, 0, 1}, {η, 0, 1}];
edges = Map[
First@Cases[
Normal@ParametricPlot[#, {t, 10^-5, 1 - 10^-5},
PlotPoints -> 100], Line[p_] :> p, Infinity] &,
{bezierfunc[t, 0],
bezierfunc[1, t],
bezierfunc[1 - t, 1],
bezierfunc[0, 1 - t]
}
];
markerRules = Join @@ Thread /@ Thread[edges -> Range[4]];
markernf = Nearest[markerRules];
plotregion = DiscretizeGraphics[plot];
mesh = ToElementMesh[
plotregion,
"BoundaryMarkerFunction" :> (First@*markernf@*Mean /@ # &)];
Show[
mesh["Wireframe"],
mesh["Wireframe"["MeshElement" -> "BoundaryElements",
"MeshElementMarkerStyle" -> Red]],
Frame -> True, FrameTicks -> All
]
It's a little hard to see, but the η == 0
boundary elements have ElementMarker
equal to 1
.
sol = NDSolveValue[{D[u[x, y], x, x] + D[u[x, y], y, y] == -3,
DirichletCondition[u[x, y] == 0, ElementMarker == 1]},
u, {x, y} ∈ mesh];
Plot3D[sol[x, y], {x, y} ∈ Ω]
Here is another way, which is more straightforward than my other answer. At first, I got stumped by couple of things, including, it turns out, a Bug in ArcLength?, and I didn't have time to explore the issues.
Instead of using a "BoundaryMarkerFunction"
we can list the markers directly in LineElement[elements, markers]
. We can make a fairly general function that create a boundary mesh from a list of edges, marking each edge in turn by its index in the list. The edges are specified by a list of points representing a polygonal path that is part of the boundary of the region.
edgesToBoundaryMesh[edges_, opts : OptionsPattern[ToBoundaryMesh]] :=
With[{lengths = Length /@ edges, bcoords = Join @@ Most /@ edges},
ToBoundaryMesh[
"Coordinates" -> bcoords,
"BoundaryElements" -> {LineElement[
Partition[Range@Length[bcoords], 2, 1, 1],
Flatten[MapIndexed[ConstantArray[#2, #1] &, lengths - 1]]]},
opts
]]
Given a list of parametrizations defining the edges, such as
{bezierfunc[#, 0] &, bezierfunc[1, #] &, bezierfunc[1 - #, 1] &, bezierfunc[0, 1 - #] &}
it would be nice to have a function that creates a discretization of the edges. This can be done with ParametricPlot
as the OP shows in the question.
edges = Map[First@Cases[Normal@ParametricPlot[#[t], {t, 0, 1}], Line[p_] :> p, Infinity] &,
{bezierfunc[#, 0] &, bezierfunc[1, #] &, bezierfunc[1 - #, 1] &, bezierfunc[0, 1 - #] &}];
mesh = ToElementMesh@edgesToBoundaryMesh[edges];
ElementMeshWireframe@mesh
The result is okay, if a little nonuniform due to the recursive subdivision of ParametricPlot
where the boundary has greater curvature. (One could also simply use Table
instead of ParametricPlot
to generate the points directly.) The size of the cells in the middle (and overall) can be controlled by the ToElementMesh
option "MaxCellMeasure"
. The initial boundary can be controlled by ParametricPlot
options PlotPoints
and MaxRecursion
.
The ElementMarker
for the boundary for the DirichletCondition
is 1
. The solution to the PDE is obtained with
sol = NDSolveValue[{D[u[x, y], x, x] + D[u[x, y], y, y] == -3,
DirichletCondition[u[x, y] == 0, ElementMarker == 1]},
u, {x, y} ∈ mesh]
Addendum
For some reason, I was initially curious about what seemed to be somewhat excessive refinement by ParametricPlot
and also with controlling the mesh size on each segment of the boundary. When the edges have greatly different lengths, this will be desirable. If we assume that the parametrizations do not change their velocity much, then it seemed one might obtain a more uniform mesh with Table
or Range
without too much trouble. The subdivision of the interval {t, 0, 1}
of the parameterizations should depend on the arc length of the edge and maximum cell measure.
parametricEdgesToBoundaryMesh[edgefns_,
opts : OptionsPattern[ToBoundaryMesh]] :=
With[{lengths = Map[ArcLength@DiscretizeRegion@Line@Table[#[t], {t, 0., 1., 1/32}] &, edgefns]},
With[{ninteverals = Map[Ceiling[(1/If[NumericQ[#], #, Max[lengths]/50.] &@
OptionValue["MaxBoundaryCellMeasure"]) #] &, lengths]},
edgesToBoundaryMesh[
MapThread[Function[t, #[t]] /@ Rescale[Range[0, #2], {0., #2}] &, {edgefns, ninteverals}],
opts]
]];
ClearAll[parametricPatch];
SetAttributes[parametricPatch, HoldAll];
parametricPatch[f_, domain : {{u1_, u2_}, {v1_, v2_}} : {{0, 1}, {0, 1}}] :=
Activate[
Function /@ Inactive[f] @@@ Transpose@MapThread[Rescale,
{Transpose[{{#1, v1}, {u2, #1}, {u2 - #1, v2}, {u1, v2 - #1}}],
{{u1, u2}, {v1, v2}}}]
];
The utility parametricPatch
takes a parametrization f
over a rectangular domain
and returns a list of parametrizations of the edges over the unit interval.
bmesh = parametricEdgesToBoundaryMesh[
parametricPatch[bezierfunc], "MaxBoundaryCellMeasure" -> 2.];
mesh = ToElementMesh[bmesh];
mesh["Wireframe"]
This has 254 triangles compared to over 3000 with the default ParametricPlot
. This is not a great advantage in this particular case. I think the important consideration here is to divide the edge segments of the boundary into similarly sized line elements.