I'm curious to find the shape of a surface bounded between the rungs of a helix, ie the shape of the cloth stretched between the rungs of this child's play tunnel. I'm wondering if we could find it like how we find minimal surfaces: it seems like a tilted but asymmetrical catenoid.

I tried to find this surface using the code posted here but it returns a helicoid. How do I "restrict" the code to find the minimal surface bounded between the rungs? (ie a hollow centre)

enter image description here

Further question

To verify whether the surfaces produced by @Greg's code are minimal surfaces, I compared an analytical catenoid to the surface obtained when passing an identical cylinder through areaGradientDescent. However, I noticed a sharp edge that wasn't minimised: enter image description here

Which could be the reason why the obtained surface is more "minimised" than the catenoid (the outer surface is the analytical solution): enter image description here

How do I fix this? Update: this is resolved by using the original areaGradientDescent

Comparison with analytical solution

As suggested, here are 3 different mesh refinements: From left to right, MaxCellMeasure->{"Length"->0.1, 0.01, 0.005}. The outer one is the analytical solution for a catenoid of the same radius. Comparing the areas, the refined solution returns 5.98 (the increased refinement only changes the area from 5.983 to 5.979) while the analytical is 6.12. While the error may not be very large, about 2.3%, I'd like to understand how the numerical errors arise. Thanks so much!

  • $\begingroup$ I guess you generated the cylinder by a ParametricPlot3D. Then it is rectangular surface piece rolled to a cylinder, but the seems are not glued together. My code will then treat all the edges along that seem as boundary which is why the points there are fixed. $\endgroup$ Commented Nov 14, 2023 at 3:22
  • $\begingroup$ Oh I see! After inspecting the mesh further I noticed the straight edge, using back areaGradientDescent gives the smooth result. However the above observation that it "overminimises" the surface still holds unfortunately... $\endgroup$ Commented Nov 15, 2023 at 13:13
  • $\begingroup$ How much does it "overminimise"? This effect should dissappear with finer and finer meshes. $\endgroup$ Commented Nov 15, 2023 at 14:00
  • $\begingroup$ The effect doesn't seem to be disappearing for the range of refinement that my computer can handle...I'm curious to understand how the errors arise and how we could possibly improve it. Thank you so much! $\endgroup$ Commented Nov 15, 2023 at 14:39
  • $\begingroup$ To address your 'further question'. How did you generate the cylinder? My guess is there's a topological boundary along the length of the tube. FindMeshDefects could possibly hint this to you. $\endgroup$
    – Greg Hurst
    Commented Nov 15, 2023 at 20:27

4 Answers 4


Another approach is to modify areaGradientDescent, defined here, to fix user specified edges in space.

areaGradientDescentFixedEdges[{R_MeshRegion, fixededges_List}, stepsize_: 1., steps_: 10, 
 reassemble_: False] := 
Module[{method, faces, bndedges, bndvertices, pts, intvertices, pat,
  flist, A, S, solver}, 
 method = If[reassemble, "Pardiso", "Multifrontal"];
 pts = MeshCoordinates[R];
 faces = MeshCells[R, 2, "Multicells" -> True][[1, 1]];
 bndedges = Union[
   Join @@ MeshCells[R, {1, fixededges}, "Multicells" -> True][[All, 1]], 
   Developer`ToPackedArray[Region`InternalBoundaryEdges[R][[All, 1]]]
 bndvertices = Union @@ bndedges;
 intvertices = Complement[Range[Length[pts]], bndvertices];
 pat = Flatten[getLaplacianCombinatorics[faces], 1];
 flist = Flatten[faces];
 Do[A = LaplaceBeltrami[pts, flist, pat];
  If[reassemble || i == 1, 
   solver = LinearSolve[A[[intvertices, intvertices]], Method -> method]];
  pts[[intvertices]] -= stepsize solver[(A.pts)[[intvertices]]];, {i, 1, steps}];
 S = MeshRegion[pts, MeshCells[R, 2], PlotTheme -> "LargeMesh"];

And then take a similar approach to my other answer here:

helixPoints[l_, p_, r_, z0_:0] := Table[{r*Cos[2π/p * θ], r*Sin[2π/p * θ], z0 + θ}, {θ, 0.0, Round[l/p, 0.025p], 0.025p}]

minimalHelixTube[l_, p_, r_] :=
  Block[{helix1, helix2, n, m, mr, ldist, fixed, mrminimal},
    helix1 = helixPoints[l+4p, p, r, -2p];
    n = Length[helix1];
    m = Round[1/0.025];
    mr = MeshRegion[
        Polygon[Join @@ ({{#1, #1-m, #1-m+1}, {#1, #1-m+1, #1+1}}& /@ Range[m+1, n-1])]
    mr = DiscretizeRegion[mr, MaxCellMeasure -> {"Length" -> 0.1}];
    ldist = RegionDistance[MeshRegion[helix1, Line[Range[n]]]];
    fixed = Pick[
      Range[MeshCellCount[mr, 1]], 
      Threshold[ldist[AnnotationValue[{mr, 1}, MeshCellCentroid]]],
    mrminimal = areaGradientDescentFixedEdges[{mr, fixed}, 0.25, 20, False];
    DiscretizeRegion[mrminimal, {{-r, r}, {-r, r}, {0, l}}]

And test:

minimalHelixTube[6, 1, 1]

enter image description here

  • $\begingroup$ Thanks so much Greg! I think the link to "here" got lost... could you attach it? $\endgroup$ Commented Nov 13, 2023 at 8:39
  • $\begingroup$ mathematica.stackexchange.com/a/158356/4346 $\endgroup$
    – Greg Hurst
    Commented Nov 13, 2023 at 13:02
  • $\begingroup$ I was comparing your code to the analytical solution of a catenoid as a sanity check, by minimising the surface of an equivalent cylinder, but encountered an oddity. Could you take a look at the Further Question in the updates? Thanks! $\endgroup$ Commented Nov 14, 2023 at 3:10
  • $\begingroup$ Awesome! Good job! $\endgroup$ Commented Nov 14, 2023 at 3:23

We can use smoothMeshRegion from the 'Code Dump' section in this answer to make a region that resembles the desired object:

helixPoints[l_, p_, r_, z0_:0] := Table[{r*Cos[2π/p * θ], r*Sin[2π/p * θ], z0 + θ}, {θ, 0.0, Round[l/p, 0.025p], 0.025p}]

helixRidgedTube[l_, p_, r_, r2_] :=
  Block[{helix1, helix2, n, m, mr},
    helix1 = helixPoints[l+4p, p, r, -2p];
    helix2 = helixPoints[l+4p, p, r2, -1.5p];
    n = Length[helix1];
    m = Round[1/0.025];
    mr = MeshRegion[
        Join[helix1, helix2],
            Polygon[Join @@ ({{#1, #1+n, #1+n+1}, {#1, #1+n+1, #1+1}}& /@ Range[n-1])],
            Polygon[Join @@ ({{#1, #1+n-m, #1+n-m+1}, {#1, #1+n-m+1, #1+1}}& /@ Range[m+1, n-1])]
    mr = DiscretizeRegion[mr, MaxCellMeasure -> {"Length" -> 0.1}];
      smoothMeshRegion[mr, "VertexPenalty" -> 0.0125],
      {{-r, r}, {-r, r}, {0, l}}


helixRidgedTube[10, 1, 2, 1]

enter image description here


Not an answer, just too big for a comment.

I cannot reproduce OP's problems with the "overminimization". With my orignial code from here I get the following results:

radius = 1.;
height = 2.;
meshsize = 0.01;

cylinder0 = DiscretizeRegion[
    x^2 + y^2 == radius^2,
    {{x, - 1.1 radius, 1.1 radius}, {y, -1.1 radius, 1.1 radius}, {z, -height/2, height/2}}
   MaxCellMeasure -> {1 -> meshsize}

cells = MeshCells[cylinder0, 2, "Multicells" -> True];
(*The true minimal surface, at least if h/r is small.*)
catenoid = MeshRegion[({x, y, z} |-> {Cosh[z] x, Cosh[z] y, z}) @@@ MeshCoordinates[cylinder0], cells];

(*A cylinder with the same boundary conditions as the catenoid.*)
cylinder = MeshRegion[({x, y, z} |-> {Cosh[height/2] x, Cosh[height/2] y, z}) @@@ MeshCoordinates[cylinder0], cells];

discreteMinimizer = areaGradientDescent[cylinder, 1., 20., False];

Print["Catenoid area = ", Area[catenoid]]
 MeshRegion[catenoid, MeshCellStyle -> ({2, All} -> Orange)],
 PlotRange -> {All, {0, All}, All}

Initial area = 19.3909

Final area = 17.6773

Catenoid area = 17.6773

enter image description here


Bend the solenoid to a torus. This is a vast field

wikipedia helicoids

Minimal Art

Crazy bubbles


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.