Find the central point in each Sierpinski triangle

Question

I'm trying to do something like this for higher orders.

What I've done is using a code available on internet to generate the Sierpinski points, and then TriangleCenter to obtain the medium points but in a manual way. For l=4, for example, this become very tedious to do manually

sierpinski[{a_, b_, c_}] := 
  With[{ab = (a + b)/2, bc = (b + c)/2, 
    ca = (a + c)/2}, {{a, ab, ca}, {ab, b, bc}, {ca, bc, c}}];

pts = {{0, 0}, {1, 0}, {1/2, Sqrt[3]/2}} // N;
l = 2;
d = Nest[Join @@ sierpinski /@ # &, {pts}, 
  l]; 

SGg = DeleteDuplicates[Sort[Flatten[d, 1]]];

(*L=1*)
(*MEDIUM={TriangleCenter[{SGg[[1]],SGg[[2]],SGg[[3]]},"Centroid"],\
TriangleCenter[{SGg[[3]],SGg[[5]],SGg[[6]]},"Centroid"],\
TriangleCenter[{SGg[[2]],SGg[[4]],SGg[[5]]},"Centroid"]};*)

(*L=2*)
MEDIUM = {TriangleCenter[{SGg[[1]], SGg[[2]], SGg[[3]]}, "Centroid"], 
   TriangleCenter[{SGg[[3]], SGg[[5]], SGg[[7]]}, "Centroid"], 
   TriangleCenter[{SGg[[7]], SGg[[10]], SGg[[12]]}, "Centroid"], 
   TriangleCenter[{SGg[[14]], SGg[[15]], SGg[[12]]}, "Centroid"], 
   TriangleCenter[{SGg[[2]], SGg[[4]], SGg[[5]]}, "Centroid"], 
   TriangleCenter[{SGg[[10]], SGg[[13]], SGg[[14]]}, "Centroid"], 
   TriangleCenter[{SGg[[4]], SGg[[6]], SGg[[8]]}, "Centroid"], 
   TriangleCenter[{SGg[[8]], SGg[[11]], SGg[[13]]}, "Centroid"], 
   TriangleCenter[{SGg[[6]], SGg[[9]], SGg[[11]]}, "Centroid"]};

SG = Union[SGg, MEDIUM]; 

How could I optimize this?

Answer 1

If I interpret your code correctly, then d contains the coordinates of the vertices of all triangles. To start from there, you may map "TriangleCenter" onto d and then apply "Point" to the result.

Further. To draw the triangles, you may either map "Triangle" on d:

l = 3;
d = Nest[Join @@ sierpinski /@ # &, {pts}, l];
Graphics[{Green, Triangle /@ d, Red, Point[TriangleCenter /@ d]}]

Or you may draw the polygon that makes up the triangle. However, for this you need to add the last stroke that closes the triangle:

Graphics[{Line@Append[#, #[[1]]] & /@ d, Point[TriangleCenter /@ d]}]

Answer 2

One way is to retrieve the MeshCellCentroid property for each SierpinskiMesh cell.

mesh = SierpinskiMesh[3];

centers = AnnotationValue[{mesh, 2}, MeshCellCentroid];

Show[mesh, Graphics[Point[centers]]]

---

Per the comments, we can get exact centers by injecting an undocumented option into MeshRegion. As with anything undocumented, no promises that this will work in future versions.

Internal`InheritedBlock[{MeshRegion},
  Unprotect[MeshRegion];
  MeshRegion[args__] /; FreeQ[{args}, WorkingPrecision -> _, {1}] := 
     MeshRegion[args, WorkingPrecision -> ∞];

  mesh = SierpinskiMesh[3];
]

centers = AnnotationValue[{mesh, 2}, MeshCellCentroid];

centers[[1 ;; 3]]

{{1/16, 1/(16 Sqrt[3])}, {3/16, 1/(16 Sqrt[3])}, {1/8, 1/(4 Sqrt[3])}}