5
$\begingroup$

This question already has an answer here:

I'm attempting to create a certain shape. It's essentially a rectangular prism with an array of cylinders cut out of it, and then some half shells of cylinders along the sides. My attempt is basically this:

  1. Create a list of these cylinders, each defined by ImplicitRegion, then combine them into a single region with RegionUnion.
  2. Define the rectangular prism, also with ImplicitRegion
  3. Define the total region with ImplicitRegion by including the points that are a part of the prism, but not the cylinder array
  4. Plot this resulting region with RegionPlot3D

Here's the code for that, just plotting the cylinder array to begin:

 hCub = 1.5;
W = 5;
L = 8;
CYLrad = W/20;
shellthickness = CYLrad/6;
deltah = hCub*0;
pitch = 1;
xnum = 5;
ynum = 6;
xinnercylnum = xnum - 2;
xpitch = W/(xnum - 1);
ypitch = xpitch*Sqrt[3.]/2;
Clear[x, y, z];
firstcylpos = {xpitch, ypitch, 0};
cyltable = 
  RegionUnion@Flatten[#, 1] &@
   Table[Table[
     ImplicitRegion[
      0 <= (x - (j - 1)*xpitch)^2 + (y - i*ypitch)^2 <= CYLrad^2 && 
       0 <= (z - 0) <= hCub + deltah(*&&
      x\[LessEqual]cylcenter[[1]]*), {x, y, z}], {i, 0, ynum, 1}], {j,
      xnum}];
innercyltable = 
  RegionUnion@Flatten[#, 1] &@
   Table[Table[
     ImplicitRegion[
      0 <= (x - j*xpitch)^2 + (y - i*ypitch)^2 <= CYLrad^2 && 
       0 <= (z - 0) <= hCub + deltah(*&&
      x\[LessEqual]cylcenter[[1]]*), {x, y, z}], {i, ynum}], {j, 
     xinnercylnum}];
yzface1cylshelltable = 
  RegionUnion@
   Table[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= (x - W)^2 + (y - i*ypitch)^2 <= 
       CYLrad^2 && 0 <= (z - 0) <= hCub + deltah && x <= W, {x, y, 
      z}], {i, ynum}];
yzface2cylshelltable = 
  RegionUnion@
   Table[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= (x - 0)^2 + (y - i*ypitch)^2 <= 
       CYLrad^2 && 0 <= (z - 0) <= hCub + deltah && x >= 0, {x, y, 
      z}], {i, ynum}];
xzface2cylshelltable = 
  RegionUnion@
   Table[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= (x - i*xpitch)^2 + (y - 0)^2 <= 
       CYLrad^2 && 0 <= (z - 0) <= hCub + deltah && y >= 0, {x, y, 
      z}], {i, xinnercylnum}];
allcylshells = 
  RegionUnion[yzface1cylshelltable, yzface2cylshelltable, 
   xzface2cylshelltable];
cub = ImplicitRegion[
   0 <= x <= W && 0 <= y <= L && 0 <= z <= hCub, {x, y, z}];
cubnocyls = RegionDifference[cub, cyltable];
Timing@RegionPlot3D[{Evaluate@cubnocyls, Evaluate@innercyltable, 
   Evaluate@allcylshells}, 
  PlotStyle -> {{Black, Lighting -> "Neutral", 
     Opacity@1}, {Lighting -> "Neutral", 
     White}, {Lighting -> "Neutral", White}}, PlotPoints -> 20, 
  Boxed -> False, ImageSize -> 800]

Doing this takes about 12 second and produces:

enter image description here

This isn't an especially complicated structure (relatively), and it's at a much lower resolution than I'll eventually need to do (PlotPoints->20). I cranked it up to PlotPoints->100, it took ~11 minutes, and it still looked pretty bad.

I'm wondering if I'm just doing this in a very naive way, and there's a more efficient way to do it. Is there?

$\endgroup$

marked as duplicate by e.doroskevic, user9660, MarcoB, Edmund, gpap Jul 19 '16 at 8:57

This question has been asked before and already has an answer. If those answers do not fully address your question, please ask a new question.

  • 5
    $\begingroup$ I think you've asked a similar question already :s I think this was the reference for more information: mathematica.stackexchange.com/questions/48486/… $\endgroup$ – e.doroskevic Jul 14 '16 at 22:03
  • $\begingroup$ @E.Doroskevic I saw that other thread, and it indeed answers my previous thread, but I don't think it answers this one. For example, I tried using ContourPlot3D on my array of cylinders (using their logical expressions, not their ImplicitRegion's), and it still took ~90s, for a pretty bad rendering of them. $\endgroup$ – YungHummmma Jul 15 '16 at 17:21
12
$\begingroup$

We can work in 2D and use RegionProduct to extrude into 3D. I've modified your code slightly:

hCub = 1.5;
W = 5;
L = 8;
CYLrad = W/20;
shellthickness = CYLrad/6;
deltah = hCub*0;
pitch = 1;
xnum = 5;
ynum = 6;
xinnercylnum = xnum - 2;
xpitch = W/(xnum - 1);
ypitch = xpitch*Sqrt[3.]/2;
firstcylpos = {xpitch, ypitch, 0};

mr = DiscretizeRegion[ImplicitRegion[
  0 <= x^2 + y^2 <= CYLrad^2, {x, y}]];
cyltable = 
  RegionUnion@Flatten[#, 1] &@
   Table[Table[TransformedRegion[mr, 
   TranslationTransform[{(j - 1)*xpitch, i*ypitch}]], {i, 0, ynum, 1}], {j,xnum}];
innercyltable = 
  RegionUnion@Flatten[#, 1] &@
   Table[Table[
     TransformedRegion[mr, 
   TranslationTransform[{j*xpitch, i*ypitch}]], {i, ynum}], {j, 
     xinnercylnum}];

mr = DiscretizeRegion[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= (x - W)^2 + y^2 <= 
       CYLrad^2 && x <= W, {x, y}]];
    yzface1cylshelltable = 
  RegionUnion@
   Table[TransformedRegion[mr, 
   TranslationTransform[{0, i*ypitch}]], {i, ynum}];

mr = DiscretizeRegion[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= x^2 + y^2 <= 
       CYLrad^2 && x >= 0, {x, y}]];
yzface2cylshelltable = 
  RegionUnion@
   Table[TransformedRegion[mr, 
   TranslationTransform[{0, i*ypitch}]], {i, ynum}];

mr = DiscretizeRegion[ImplicitRegion[(CYLrad - 
          shellthickness)^2 <= x^2 + y^2 <= 
       CYLrad^2 && y >= 0, {x, y}]];
xzface2cylshelltable = 
  RegionUnion@
   Table[TransformedRegion[mr, 
   TranslationTransform[{i*xpitch, 0}]], {i, xinnercylnum}];

allcylshells = 
  RegionUnion[yzface1cylshelltable, yzface2cylshelltable, 
   xzface2cylshelltable];

cub = DiscretizeGraphics[Rectangle[{0, 0}, {W, L}]];
cubnocyls = RegionDifference[cub, cyltable];

extruded = RegionProduct[#, Line[{{0}, {hCub + deltah}}]]& /@ 
  {cubnocyls, innercyltable, allcylshells};

styles = {{Black, Lighting -> "Neutral", 
     Opacity@1}, {Lighting -> "Neutral", 
     White}, {Lighting -> "Neutral", White}};

Show[MapThread[MeshRegion[#1, BaseStyle -> #2] &, {extruded, styles}]]

enter image description here

$\endgroup$
  • $\begingroup$ This looks pretty great and ran fairly quickly for me. However, when I run exactly this code, it's displaying with all these lines all over the place. I looked for a place to increase the PlotPoints or something but couldn't find anything. Here is a pic of what I'm getting, any idea what I could change? imgur.com/ccE9oY3 $\endgroup$ – YungHummmma Jul 15 '16 at 17:35
  • $\begingroup$ @YungHummmma Add EdgeForm[] inside all 3 lists inside styles. Versions 10.4 and higher automatically hide the lines when there's too many, so I imagine you're using an older version? $\endgroup$ – Chip Hurst Jul 15 '16 at 17:39
  • $\begingroup$ Hmm, I'm actually using version 10.4.1. I'll try that in a minute. $\endgroup$ – YungHummmma Jul 15 '16 at 18:00
  • $\begingroup$ Hi from long ago, but I actually needed to do this again, and adding EdgeForm[] solves the problem. I think this works great, thank you! $\endgroup$ – YungHummmma Jan 4 '17 at 21:41
  • $\begingroup$ That was a long minute! Anyway, I guess it's version 11 that hides the edges automatically, not 10.4 $\endgroup$ – Chip Hurst Jan 4 '17 at 22:09

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