4
$\begingroup$

For finite element purposes I need a frustum of a cone with a finite wall thickness, i.e. a tapered pipe.

To make a cone is easy

Needs["NDSolve`FEM`"]
c = BoundaryDiscretizeRegion[Cone[{{0, 0, 0}, {-1, 0, 0}}, 0.12]]

Mathematica graphics

mesh = ToElementMesh[c]; mesh["Wireframe"]

Mathematica graphics

So my idea was to use RegionDifference to make a frustum. Here is my module.

ClearAll[coneFrustum];
coneFrustum::usage = 
  "coneFrustum[{pt1, pt2, },r1 ,r2 ,height] generates a frustum of \
cone with axis going from pt1 to pt2 and with radii r1, r2";
coneFrustum[{{x1_, y1_, z1_}, {x2_, y2_, z2_}}, r1_, r2_, L_] := 
 Module[{cL, vec},
  cL = L r1/(r1 - r2);
  vec = ({x2, y2, z2} - {x1, y1, z1})/
   Norm[{x2, y2, z2} - {x1, y1, z1}];
  RegionDifference[
   BoundaryDiscretizeRegion[
    Cone[{{x1, y1, z1}, cL vec}, r1]
    ],
   BoundaryDiscretizeRegion[
    Cone[{{x2, y2, z2}, (cL - L) vec}, r2]
    ]
   ]
  ]

This is what I get for a solid frustum

cf1 = coneFrustum[{{-0.1, 0, 0}, {-0.3, 0, 0}}, 0.12, 0.1, 0.2]

Mathematica graphics

So I have my frustum but I have not managed to take away the part of the cone beyond the frustum. I have a hollow shell of the difference of the cones.

Is there a better way of doing this? I still have to go on to subtract a smaller inner frustum to make my hollow shape.

Edit

As kglr pointed out I made a simple error with my code. The tip has to have the same coordinates. Here is the corrected version.

ClearAll[coneFrustum];
coneFrustum::usage = 
  "coneFrustum[{pt1, pt2}, r1,r2, height] generates a frustum of cone \
with axis going from pt1 to pt2 and with radii r1, r2";
coneFrustum[{{x1_, y1_, z1_}, {x2_, y2_, z2_}}, r1_, r2_, L_] := 
 Module[{cL, vec},
  cL = L r1/(r1 - r2);
  vec = ({x2, y2, z2} - {x1, y1, z1})/
   Norm[{x2, y2, z2} - {x1, y1, z1}];
  RegionDifference[
   BoundaryDiscretizeRegion[
    Cone[{{x1, y1, z1}, cL vec}, r1]
    ],
   BoundaryDiscretizeRegion[
    Cone[{{x2, y2, z2}, cL vec}, r2]
    ]
   ]
  ]

Now we make two frustums and subtract.

cf1 = coneFrustum[{{-0.1, 0, 0}, {-0.3, 0, 0}}, 0.12, 0.1, 0.2];
cf2 = coneFrustum[{{-0.1, 0, 0}, {-0.3, 0, 0}}, 0.1, 0.08, 0.2];
reg = RegionDifference[cf1, cf2]

Mathematica graphics

And to get the mesh

mesh = ToElementMesh[reg]; mesh["Wireframe"]

Mathematica graphics

The method from Henrik Schumacher below using ImplicitRegion is also a good solution.

Are there any good reasons why one could be better than the other?

$\endgroup$

3 Answers 3

4
$\begingroup$
c = BoundaryDiscretizeRegion[Cone[{{0, 0, 0}, {-1, 0, 0}}, 0.12]];
k = .2;
d = BoundaryDiscretizeRegion[Cone[{{-k, 0, 0}, {-1, 0, 0}}, (1 - k) 0.12]];
RegionDifference[c, d]

enter image description here

 

$\endgroup$
5
  • 1
    $\begingroup$ And RegionDifference[RegionDifference[c, d], RegionDifference[a, b]] with a = BoundaryDiscretizeRegion[Cone[{{0, 0, 0}, {-1, 0, 0}}, 0.1]]; b = BoundaryDiscretizeRegion[ Cone[{{-k, 0, 0}, {-1, 0, 0}}, (1 - k) 0.1]]; cuts out the interior nicely. Well done (+1). $\endgroup$ Oct 16, 2018 at 18:52
  • 1
    $\begingroup$ Why does yours work and not mine? I think it is the same idea... $\endgroup$
    – Hugh
    Oct 16, 2018 at 18:54
  • $\begingroup$ @Hugh, the two cones should have the same tip, so if you change (cL - L) vec to cL vec your code gives the same result. $\endgroup$
    – kglr
    Oct 16, 2018 at 19:06
  • $\begingroup$ Thank you @Henrik. $\endgroup$
    – kglr
    Oct 16, 2018 at 19:07
  • $\begingroup$ Rats! Simple error. Back to school for me. Thanks. $\endgroup$
    – Hugh
    Oct 16, 2018 at 19:09
3
$\begingroup$

You mean something like this?

n = 36;
rout1 = 1.;
rout2 = 0.8;
rin1 = .9;
rin2 = 0.7;
h1 = 1.;
h2 = 0.;
R = MeshRegion[
 Join[
  Join[CirclePoints[rout1, n], ConstantArray[h1, {n, 1}], 2], 
  Join[CirclePoints[rout2, n], ConstantArray[h2, {n, 1}], 2], 
  Join[CirclePoints[rin1, n], ConstantArray[h1, {n, 1}], 2], 
  Join[CirclePoints[rin2, n], ConstantArray[h2, {n, 1}], 2]
  ],
 With[{
   edges = Partition[#, 2, 1, 1] & /@ Partition[Range[4 n], n]
   },
  Polygon@Join[
    Join[edges[[1]], Reverse /@ edges[[2]], 2],
    Join[edges[[3]], Reverse /@ edges[[4]], 2],
    Join[edges[[1]], Reverse /@ edges[[3]], 2],
    Join[edges[[2]], Reverse /@ edges[[4]], 2]
    ]
  ]
 ]

enter image description here

Needs["NDSolve`FEM`"]
mesh = ToElementMesh[R]; mesh["Wireframe"]

enter image description here

This is quick and robust but it is not so easy to create adaptive meshes with this method.

Edit:

Here is another way utilizing ImplicitRegion.

BoundaryDiscretizeRegion[
 ImplicitRegion[{0.5 <= z <= 1, 0.1 z^2 <= x^2 + y^2 <= 0.12 z^2}, {{x, -1, 1}, {y, -1, 1}, {z, 0.5, 1}}],
 MaxCellMeasure -> (1 -> 0.02)
 ]

BoundaryDiscretizeRegion has a bit of trouble to detect all parts of the mantle for larger values of MaxCellMeasure, so this is not super robust. But is allows the use of MeshRefinementFunction in order to create adaptive meshes.

enter image description here

$\endgroup$
2
  • $\begingroup$ Gosh, that was so fast I would guess you had done this before. Thanks. You seem to be recommending building from polygons. What happens if I want a fine mesh at one location that is part of a curved surface rather than built-up from polygons? I think with your method I have to put in lots of polygons rather than use a mesh refinement function. $\endgroup$
    – Hugh
    Oct 16, 2018 at 18:48
  • $\begingroup$ You're welcome. And yeah, you're right. A true adaptive discretization is not possible this way. =/ At least not that easily. $\endgroup$ Oct 16, 2018 at 18:50
1
$\begingroup$

In recent version you can directly use Boolean operations to get quite good results:

Needs["NDSolve`FEM`"]
outer = Cone[{{0, 0, 0}, {-1, 0, 0}}, 0.12];
inner = Cone[{{0, 0, 0}, {-0.8, 0, 0}}, 0.10];
box = Cuboid[{-1/4, -1, -1}, {0, 1, 1}];
reg = RegionDifference[outer, inner]; reg = 
 RegionIntersection[reg, box];
MeshRegion[ToBoundaryMesh[reg]]

enter image description here

$\endgroup$

Your Answer

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

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