How can I draw a pyramid of spheres like this?

Question

I want to draw this figure such that the height of the pyramid of spheres is h (this picture with h = 5. I only tried

r = 1;
s1 = Graphics3D[{Sphere[{0, 0, 0}, r], Sphere[{2 r, 0, 0}, r], 
   Sphere[{r, r*Sqrt[3], 0}, r]}, Boxed -> False]

How can I draw it?

Answer 1

Manipulate[
Graphics3D[{Specularity[White,50],
Table[{Hue[i/n],Sphere[{k-j/2,i/√3-√3j/2,-√6i/3},r]},{i,0,n-1},{j,0,i},{k,0,j}]
},Lighting->{{"Directional", White, ImageScaled[{2,2,2}]}},
Method->{"SpherePoints"->150},Background->Black],{{n,8},3,15,1},{{r,0.5},0.1,1}]

n = 6;
mat = {{0,1,-1},{0,-1,0},{1,0,0}} . {{1,0,0}, {1/2,√3/2,0}, {1/2,√3/6,√6/3}};
pts = 2 Reverse@GroupTheory`Tools`Multisets[Range[n - 1, 0, -1], 3] . mat;
Graphics3D[{
  {Opacity[0.5, Blue], Specularity[1, 20], Sphere[pts]},
  MapIndexed[Inset[Style[Tr@#2, 20, White, Bold], #] &, pts]
  }, ImageSize -> 600, Boxed -> False, Background -> Black]

Answer 2

e1 = {1, 0, 0};
e2 = RotationTransform[π/3, {0, 0, 1}]@e1;
e3 = SolveValues[{Norm@{x, y, z} == 1, 
     VectorAngle[{x, y, z}, e2] == π/3, 
     VectorAngle[{x, y, z}, e1] == π/3, z > 0}, {x, y, z}][[1]];
(* {e1,e2,e3}={{1,0,0},{1/2,Sqrt[3]/2,0},{1/2,1/(2 \
Sqrt[3]),Sqrt[2/3]}} *)
n = 5;
sol = SolveValues[{0 <= {x, y, z} <= n, 0 <= x + y + z <= n}, {x, y, 
    z}, Integers];
Graphics3D[{Blue, Specularity[White, 10], 
  Table[Ball[pt, 1/2], {pt, sol . {e1, e2, e3}}]}, Boxed -> False]

• We can also test another inequalities.

sol = SolveValues[{-n <= {x, y, z} <= n, -n <= x + y + z <= n}, {x, y,
     z}, Integers];

Answer 3

Minimally modifying my answer here:

With[{n = 5}, 
     Graphics3D[Sphere[Flatten[Table[{{2, 1, -1}, {0, Sqrt[3], -1/Sqrt[3]},
                                      {0, 0, -2 Sqrt[2/3]}} . {i, j, k},
                                     {k, 0, n}, {j, 0, k}, {i, 0, k - j}], 2]], 
                Boxed -> False]]

Answer 4

Here's an approach...

Start with a tetrahedron. The vectors from the peak to each of the other vertices make up a sort of generator set. So, given a starting point, adding each of these vectors to it gives us a new layer of points. Adding our generator vectors to each of these gives us another new layer (there will be duplicates this time, but we'll worry about that later). We can do this until we have all the layers we want.

We can easily get the data we need to start this off:

tetraVertices = PolyhedronData["Tetrahedron", "VertexCoordinates"]

It doesn't really matter which one we consider the peak, but it turns out that the first vertex in this list is the top one, and that appeals to me, so let's use it to derive the three vectors we need:

tetraEdgeVectors = (# - tetraVertices[[1]] &) /@ Drop[tetraVertices, 1]

This whole process of add vectors to the previous layer of points sounds like Nest, and actually we want all layers, so we'll use NestList:

PyramidCenters[n_Integer?Positive] :=
  Catenate[
    NestList[
      DeleteDuplicates[Plus @@@ Tuples[{#, tetraEdgeVectors}]] &, 
      {{0, 0, 0}}, 
      n - 1]]

I bundled a bunch of steps together there. I haven't worried about performance for just a few layers; I just went for clarity (well, such as it is :) ).

Now, we need to know the radius of each sphere. We can get the original edge length and divide by two:

PolyhedronData["Tetrahedron", "EdgeLengths"]
(* gives {1}, so the radius we want is .5  *)

Put this all together as a graphic:

Graphics3D[Sphere[PyramidCenters[5], .5]]