17
$\begingroup$

Searching the web for information about the affine transformation, I found the one page, which called my attention for the tree that show and is this

tree

but unfortunately do not give information about the algorithm to create it, I would like to ask help to make one the same or very similar, maybe someone knows where to get more information about it. Thanks in advance, here is the link to the page mentioned above

Improve this question
$\endgroup$
6
  • $\begingroup$ This and this. $\endgroup$ Commented Aug 21, 2017 at 23:54
  • $\begingroup$ @ corey979 Thank you very much for the information, read it with great attention to understand the methods for the generation of this type of graphics. $\endgroup$ Commented Aug 22, 2017 at 4:49
  • 2
    $\begingroup$ In case you don't know, it's common courtesy to link to the corresponding Community question if you ask something here and at Community. (Likewise, if you post there, link to your question here.) $\endgroup$ Commented Aug 22, 2017 at 10:37
  • $\begingroup$ I want to apologize to everyone for not putting the link that leads to the community as you suggest, I really did not know, but I will not forget, greetings. $\endgroup$ Commented Aug 27, 2017 at 0:39
  • $\begingroup$ You might be interested in so-called "algorithmic botany" or "computational botany". $\endgroup$ Commented Aug 6, 2019 at 8:16

2 Answers 2

Reset to default
39
$\begingroup$

First, an idomatic, but slow version.

s1 = 1/GoldenRatio // N; s2 = 1/GoldenRatio // N; stem = {0., 0., 1.}; thickness = 0.15; branches = Table[RotationMatrix[2. k Pi/3., {0, 0, 1}].{Cos[Pi/4.], 0., Sin[Pi/4.]}, {k, 0, 2}]; data0 = {Join[{{0., 0., 0.}}, {stem}, branches, {{thickness, 1., 0.}}]}; iteration[data_] := Block[{U}, Flatten[Table[ U = data[[j]]; Table[ Join[{U[[1]] + U[[2]]}, {U[[i]]}, s1 U[[3 ;; 5]].RotationMatrix[{U[[i]], U[[2]]}], {s2 U[[6]]}], {i, 3, 5}], {j, 1, Length[data]} ], 1 ] ]

This generates the tree structure.

result = NestList[iteration, data0, 6]; // AbsoluteTiming (* {0.211536, Null} *)

This generates the tree plot.

t = 0.5; colfun[x_] := ColorData["Rainbow"][t + (1 - t) x]; plot[U_] := {colfun[U[[6, 2]]],Table[Sphere[U[[1]] + t U[[2]], U[[6, 1]] (1 - t) + t s2 U[[6, 1]]], {t, 0.0, 0.9, 0.1}]}; Graphics3D[ Flatten[plot /@ Flatten[result, 1]], Lighting -> "Neutral", Background -> Black, Boxed -> False, SphericalRegion -> True ]

enter image description here

A faster version is generated with Compile and some handcraft:

citeration2 = With[{scale1 = s1, scale2 = s2, part = Compile`GetElement}, Compile[{{U, _Real, 2}}, Block[{A, u, v, w}, v = {part[U, 2, 1], part[U, 2, 2], part[U, 2, 3]}/Sqrt[part[U, 2, 1]^2 + part[U, 2, 2]^2 + part[U, 2, 3]^2]; Table[ u = {part[U, i, 1], part[U, i, 2], part[U, i, 3]}/Sqrt[part[U, i, 1]^2 + part[U, i, 2]^2 + part[U, i, 3]^2]; w = { -(part[u, 3] part[v, 2]) + part[u, 2] part[v, 3], part[u, 3] part[v, 1] - part[u, 1] part[v, 3], -(part[u, 2] part[v, 1]) + part[u, 1] part[v, 2] }; w = {part[w, 1], part[w, 2], part[w, 3]}/Sqrt[part[w, 1]^2 + part[w, 2]^2 + part[w, 3]^2]; A = { { part[u, 1] part[v, 1] + part[w, 1]^2 + (part[u, 3] part[w, 2] - part[u, 2] part[w, 3]) (part[v, 3] part[w, 2] - part[v, 2] part[w, 3]), part[u, 2] part[v, 1] + part[w, 1] part[w, 2] + (-(part[u, 3] part[w, 1]) + part[u, 1] part[w, 3]) (part[v, 3] part[w, 2] - part[v, 2] part[w, 3]), part[u, 3] part[v, 1] + part[w, 1] part[w, 3] + (part[u, 2] part[w, 1] - part[u, 1] part[w, 2]) (part[v, 3] part[w, 2] - part[v, 2] part[w, 3]) }, { part[u, 1] part[v, 2] + part[w, 1] part[w, 2] + (part[u, 3] part[w, 2] - part[u, 2] part[w, 3]) (-(part[v, 3] part[w, 1]) + part[v, 1] part[w, 3]), part[u, 2] part[v, 2] + part[w, 2]^2 + (-(part[u, 3] part[w, 1]) + part[u, 1] part[w, 3]) (-(part[v, 3] part[w, 1]) + part[v, 1] part[w, 3]), part[u, 3] part[v, 2] + part[w, 2] part[w, 3] + (part[u, 2] part[w, 1] - part[u, 1] part[w, 2]) (-(part[v, 3] part[w, 1]) + part[v, 1] part[w, 3]) }, { part[u, 1] part[v, 3] + part[w, 1] part[w, 3] + (part[v, 2] part[w, 1] - part[v, 1] part[w, 2]) (part[u, 3] part[w, 2] - part[u, 2] part[w, 3]), part[u, 2] part[v, 3] + part[w, 2] part[w, 3] + (part[v, 2] part[w, 1] - part[v, 1] part[w, 2]) (-(part[u, 3] part[w, 1]) + part[u, 1] part[w, 3]), part[u, 3] part[v, 3] + (part[u, 2] part[w, 1] - part[u, 1] part[w, 2]) (part[v, 2] part[w, 1] - part[v, 1] part[w, 2]) + part[w, 3]^2 } }; Join[{part[U, 1] + part[U, 2]}, {part[U, i]}, scale1 U[[3 ;; 5]].A, {scale2 part[U, 6]}], {i, 3, 5}]], CompilationTarget -> "C", RuntimeAttributes -> {Listable}, Parallelization -> True, RuntimeOptions -> "Speed" ] ]; iteration2[data_] := Flatten[citeration2[data], 1]; result2 = NestList[iteration2, data0, 6]; // AbsoluteTiming Max[Abs[result2 - result]] (* {0.001042, Null} *) (* 1.33227*10^-15 *) result3 = NestList[iteration2, data0, 9]; // AbsoluteTiming Graphics3D[ Flatten[plot /@ Flatten[result3, 1]], Lighting -> "Neutral", Background -> Black, Boxed -> False, SphericalRegion -> True ] (* {0.018179, Null} *)

enter image description here

The slow part is the rendering by Mathematica, though...

Improve this answer
$\endgroup$
6
  • $\begingroup$ What is getp? $\endgroup$ Commented Aug 22, 2017 at 4:03
  • $\begingroup$ Ah. Wait a sec. I am still building on this... $\endgroup$ Commented Aug 22, 2017 at 4:04
  • $\begingroup$ @J.M. Hm. I guess, it's done. Have a look. $\endgroup$ Commented Aug 22, 2017 at 4:31
  • 3
    $\begingroup$ @HenrikSchumacher Astonishing, the tree is the same as the one on the website, mentioned separately are the improvements you propose. I am learning more and more with these kinds of answers. Thanks for your answer $\endgroup$ Commented Aug 22, 2017 at 4:57
  • $\begingroup$ You're welcome! Actually, it was quite some fun! $\endgroup$ Commented Aug 22, 2017 at 5:54
11
$\begingroup$

Slow version

Clear["`*"]; s=1./GoldenRatio; thickness = 0.15; next[{a_,b_}]:=Table[{a,b}//TranslationTransform[b-a]// RotationTransform[Pi/4,{Cos[2k Pi/3],Sin[2k Pi/3],0},b]// ScalingTransform[{1,1,1}s,b],{k,3}]; n=5; pts=NestList[Join@@next/@#&,N@{{{0,0,0},{0,0,1}}},n];//AbsoluteTiming Graphics3D[{Tube[Join@@pts,0.02]}] Graphics3D[MapIndexed[With[{id=#2[[1]]},{ColorData["Rainbow",1-id/10], MapIndexed[Sphere[#,t=#2[[1]]/10; k=thickness s^id;k(1-t) + t k s]&, Subdivide[#[[1]],#[[2]],9]]}]&,pts,{2}]]

Faster version

Clear["`*"]; s=1./GoldenRatio; thickness = 0.15; next=Table[{a,b}//TranslationTransform[b-a]// RotationTransform[Pi/4,{Cos[2k Pi/3],Sin[2k Pi/3],0},b]// ScalingTransform[{1,1,1}s,b],{k,3}]/. Thread[{a,b}->Table[Indexed[A,{x,y}],{x,2},{y,3}]]/. expr_:>Compile[{{A,_Real,2}}, expr,RuntimeAttributes->{Listable} ]; n=9; pts=MapIndexed[Flatten[#,#2[[1]]-1]&,NestList[next,N@{{{0,0,0},{0,0,1}}},n]];//AbsoluteTiming (* {0.0085782, Null} *)

enter image description here enter image description here enter image description here

Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.