Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I wanted to draw (lattice-aligned) unit cubes inside a sphere with a given radius centered on the origin, so I wrote this program:

Cube[x_, y_, z_] := Cuboid[{x - 1, y - 1, z - 1}, {x, y, z}];
Coords[r_] := {#1, #2, Floor[Sqrt[r^2 - #1^2 - #2^2]]} & @@ # & /@ 
   With[{t = Sqrt[r^2 - 1]}, 
    Select[Flatten[Table[{x, y}, {x, 1, t}, {y, 1, x}], 1], 
     Norm[#] <= t &]];
Cubes[r_] := 
  Cube @@ # & /@ 
   Union[Flatten[Permute[#, SymmetricGroup[3]] & /@ Coords[r], 1]];
Draw[r_] := 
  Graphics3D[
   Union[Cubes[r], {{Green, Opacity[0.1], Sphere[{0, 0, 0}, r]}}], 
   PlotRange -> {{0, r}, {0, r}, {0, r}}, 
   ViewPoint -> {r, 3 r/4, 3 r/5}];
Draw /@ Sqrt /@ {50, 100, 250, 500, 1000, 2500, 6054}

(Note that it shows only the cubes on the outside -- drawing the inner ones takes too long and usually Mathematica runs out of memory.) It gets the job done, but it's pretty ugly. I'm not a regular Mathematica user, and tips for improving my code?

That is, how can I write idiomatically-better Mathematica code?


Improved code:

Coords[r_] := {#1, #2, Floor[Sqrt[r^2 - #1^2 - #2^2]]} & @@@ 
   With[{t = Sqrt[r^2 - 1]}, 
    Select[Join @@ Table[{x, y}, {x, t}, {y, 1, x}], Norm[#] <= t &]];
Cubes[r_] := (Cuboid /@ (Union @@ (Permute[#, SymmetricGroup[3]] & /@ 
         Coords[r]) - 1));
Draw[r_] := 
  Graphics3D[
   Union[Cubes[r], {{Green, Opacity[0.1], Sphere[{0, 0, 0}, r]}}], 
   PlotRange -> {{0, r}, {0, r}, {0, r}}, 
   ViewPoint -> {r, 3 r/4, 3 r/5}];
Draw /@ Sqrt /@ {50, 100, 250}
share|improve this question
    
Without examining the algorithm this looks like pretty good code to me. In what way do you find it "ugly?" I immediately see a few ways to make it shorter, e.g. {x - 1, y - 1, z - 1} can be written {x, y, z} - 1; are you interested in such things or does that just make it harder to read? –  Mr.Wizard Mar 19 at 15:28
    
@Mr.Wizard: I'll take whatever tips I can get -- coding style, optimizations, etc. I'm new to Mathematica (but not to programming) and so to me it looks strange to see functional code on top of other styles, and I had to break the code into more helper functions than I wanted just to make it readable. It still seems very hard to read (vs. write). –  Charles Mar 19 at 15:40
1  
You can drop Cube since it serves no purpose. E.g. Graphics3D[Cuboid[{2, 3, 4}], Axes -> True] will produce a unit cube aligned with the grid. –  David Carraher Mar 19 at 15:42
    
@DavidCarraher: I did not know that, thanks! –  Charles Mar 19 at 15:43
    
@David I feel foolish for overlooking that! –  Mr.Wizard Mar 19 at 15:44

3 Answers 3

up vote 5 down vote accepted

Focusing only on syntax I think your code looks pretty good, but there are a few things I note that I think could be improved.

  • The method f @@ # & /@ expr can be replaced with f @@@ expr, which is shorthand for Apply[f, expr, {1}]. See Apply.

  • You can use the listability of Plus to write {x, y, z} - 1 instead of {x - 1, y - 1, z - 1}

  • Union[Flatten[expr, 1]] can be written Union @@ expr

Your first three functions with these changes (among others), and David's observation regarding Cuboid:

Cube[x_, y_, z_] := Cuboid[{x, y, z} - 1];

Coords[r_] := {#1, #2, Floor[Sqrt[r^2 - #1^2 - #2^2]]} & @@@ 
   With[{t = Sqrt[r^2 - 1]},
     Select[Join @@ Table[{x, y}, {x, t}, {y, 1, x}], Norm[#] <= t &]];

Cubes[r_] := Cube @@@ Union @@ (Permute[#, SymmetricGroup[3]] & /@ Coords[r]);

Regarding style I think it is fairly personal whether or not something is easy to read, but here is an example of a different style which you may prefer:

Coords[r_] :=
  With[{t = Sqrt[r^2 - 1]},
    Cases[
      Join @@ Table[{x, y}, {x, t}, {y, 1, x}],
      {x_, y_} /; Norm[{x, y}] <= t :>
         {x, y, ⌊ Sqrt[r^2 - x^2 -y^2] ⌋ }
    ]
  ]

Fortunately Mathematica supports many different styles.

share|improve this answer
    
There's such a thing as @@@? Amazing! I thought @@ # & /@ was a ridiculous hack, but I couldn't get it to work otherwise. –  Charles Mar 19 at 15:47
    
@Charles I added an example of a different style, both in semantics and indentation. I'm curious to know if you find this more pleasing. –  Mr.Wizard Mar 19 at 15:58
    
It's very different. I'm going to have to stare at it a while to understand it -- I don't even know all the symbols used, like <= and :>. –  Charles Mar 19 at 16:07
    
@Charles Okay. You can select any of those tokens in the Notebook interface and press F1 to bring up the documentation page. <= is LessEqual, :> is RuleDelayed, and /; is Condition. Once you become familiar with Mathematica patterns I think you will find it is not hard to read. –  Mr.Wizard Mar 19 at 16:10
2  
Without Cube: Coords[r_]:={#1,#2,Floor[Sqrt[r^2-#1^2-#2^2]]}&@@@With[{t=Sqrt[r^2-1]},Select[J‌​oin@@Table[{x,y},{x,t},{y,1,x}],Norm[#]<=t&]]; Cubes[r_]:=(Cuboid/@(Union@@(Permute[#,SymmetricGroup[3]]&/@Coords[r])-1)); –  David Carraher Mar 19 at 16:53

Translate is probably the most efficient way to represent and display such a figure.

With[{r = 30},
 Graphics3D[{
   Translate[Cuboid[],
    Union @@ (Permute[#, SymmetricGroup[3]] & /@ Coords[r])
   ],
   {Green, Opacity[0.1], Sphere[{0, 0, 0}, r]}
  }
 ]
]

Mathematica graphics

share|improve this answer
    
Ah, much better. +1 –  Mr.Wizard Mar 21 at 9:57

New, not the best but can do

r = 50;
da = Pi/(4. r);
sa = da;
ea = Pi/2. - da;

Composition[
  Graphics3D[{#, 
              GeometricTransformation[#, ReflectionTransform /@ {{0, 0, 1}, {0, 1, 1}}]
             }, Axes -> True, PlotRange -> All, ImageSize -> 600] &,
  Cuboid /@ # &,
  {# - Sign[#], #} & /@ # &,
  Select[#, FreeQ[#, 0] &] &,
  DeleteDuplicates,
  Floor,
  Flatten[#, 1] &
  ][
 Table[r {Cos[d] Sin[a], Cos[d] Cos[a], Sin[d]}, {d, sa, ea, da}, {a, 
   sa, ea, da}]
 ]

enter image description here

Old and slow This is not optimal or faster code but it is cleaner:

r = 6;
s = 0;
e = r;
d = .5;

Composition[

  Graphics3D[{#}, Axes -> True, PlotRange -> {{0, r}, {0, r}, {0, r}}] &,
  Flatten,
  Map[If[1 == Times @@ UnitStep[r - Norm /@ #], 
         {Hue@RandomReal[], Cuboid@#[[{1, -1}]]},
         ## &[]
        ] &, 
      #, {3}] &,      
  Map[Flatten[#, 2] &, #, {3}] &,
  Partition[#, {2, 2, 2}, 1] &

  ][Table[{i, j, k}, {i, s, e, d}, {j, s, e, d}, {k, s, e, d}]]

Important part of this code is function Partition[#, {2, 2, 2}, 1] & which does something like {a,b,c} --> {{a,b},{b,c}} but in 3D. It is clean but it is a big choke point. We have to create equally spaced array for this. Including points we know are to far away from origin. You can add Opacity@.7, Sphere[{0, 0, 0}, r] to Graphics3D, look out, opacity slows down rotation. enter image description here

share|improve this answer
    
Thanks especially for the comment on opacity. –  Charles Mar 20 at 18:03

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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