I'm trying to find where a "3D logarithmic spiral" converges. Background and problem
It's in quotes because I don't know what I'm talking about. I'm more of an artist who likes to dabble in math so bare with me. I know the title is a little confusing so I have made a pretty animation to illustrate what I mean. I assume this is some type of discrete logarithmic spiral. It always starts the origin with no translations, rotations or scale applied. You'll notice in the video, the 2nd cube has 3 axes coming out of it. It's transform is driving the rest of the cubes as seen in this clip. I can translate, rotate, and scale and as long as the translation and rotation are not 0, and the scale is between 0 and 1 it'll make a spiral. Anything else won't produce a spiral as far as I'm aware.
I've searched Google for about 5 hours today and haven't found anything so I either don't know the proper terms, it's never been done before which I doubt, or it's impossible which I don't know but highly doubt.
I'm trying to find the $x, y,$ and $z$ for where this would converge in non-polar/spherical form so I can plug them into my software.
Why I'm here
Some hindrances I'm having are that I've never taken anything on spherical coordinates which I assume would be useful (or even polar coordinates), I don't have a background in linear algebra at all, and I've been out of school about 8 years at this point so I'm very rusty. I've asked this question on Reddit and was given an answer but I don't understand it. They suggested,

If your initial vector is v and your linear transformation is T, then
the point it converges to is
$(1+T+T^2 +T^3...)v=\dfrac v{1-T}.$
For the purpose of computation, you need to write T as a matrix, v as a column vector.

Can come one break that down a bit more? Is $v$ my 2nd cube with the axes poking out? I don't understand what $T$ means at all. Is it even right? I'm ok with not being spoon-fed the answer but I need a little bit more than this. If the first answer is any indication, then is way above my ability.
I'm not even sure which branch of math I need. I know it's some sort of geometric series because it's decreasing and I can see it converging to some 3d point in space. That reminds me of calculus. I've taken Calc 2 at University but that was a long time ago and we never talked about 3d space. Any help would be useful.
Motivation
And in case you're wondering why I'm interested in this. I've been interested in these scaling looping animations lately and I wanted to understand the math behind them. I found some videos by the people who make them (one of the guys mentioned in this video did an animation for Justin Bieber so they're legit) and even they don't understand the math. They're using their artistic skills to get close but that can be time-consuming. I want to know the math behind it. And I refuse to believe the math hasn't been worked out already. The point where I'm trying to find is where if you scale the whole spiral from, it has this unique property of looping perfectly. And if it's scaled at an exponential rate
$\bigg(\dfrac{1}{s}\bigg)^{\frac{n(f - 1)}{a}}$
where $s$ is the scaling factor, $n$ is the number of cycles you want, $f$ is the current frame and $a$ is the length of the animation in frames, then it ends up looking like it's scaling linear. In this example, on the left you get the illusion of it scaling forever but on the right you can see what's really happening. If it doesn't scale at an exponential rate, then it starts slow and gets faster which also breaks the illusion. It's a 1m square that gets shrunk in half each time $(1 + 1/2 + 1/4 + ...)$ so I know the limit is $2$. I'm then scaling the whole spiral from $(2, 0, 0)$ by a factor of $2$ over $50$ frames. It only scales perfectly from $(2, 0, 0)$. This is what happens if I don't use the right point. It slowly drifts away and snaps back each cycle breaking the illusion. So that's why it's important for me for find this location.
Updates
A user on Reddit found a solution that works in 2D and I can confirm that it works. Their solution is as follows:
Let $x'$ and $y'$ be the point where it converges. Then it converges at:
$x' = Ax - By$
$y' = Bx + Ay$
where $A$ is the function $\frac{1 - k cos θ}{1 + k^2 - 2k cos θ}$
$B$ is $\frac{k sin θ}{1 + k^2 - 2k cos θ}$
and where $k$ is the scaling factor and $θ$ is the angle.
Update 2
I think I'm getting closer. I've been told that:

In two dimensions, the transformation T is given by


$T = k R(θ) = k [ cos θ , -sin θ ; sin θ , cos θ ]$,


and therefore


$T^n = k^n [ cos(nθ) , -sin(nθ) ; sin(nθ) , cos(nθ) ]$.


However, this cannot be easily generalized to three dimensions. You would need to explain more precisely what you mean by a rotation in three dimensions. Do you want a single rotation with respect to a given axis? If you want to apply three rotations with respect to the three coordinate axes, in what order do you apply them?

How would I know what order to apply them since I can edit any one of them at any time?
 A: Thanks for asking this, I wouldn't expect someone without a mathematical background to be able to solve such a technical question all by themselves, and admittedly the hints that you have received are rather formal. It is true however, that obtaining explicit answers for all possible parameters of the rotation is more difficult than in 2-D.
EDIT 2: I found a better way to express the formulas below and make them much easier to parse. The new way of looking at the rotation doesn't involve computing the complex eigenvectors $\mathbf{e}_{\pm}$. Consider a rotation of angle $\Theta$, around an axis $\mathbf{n}$. Then it is well known that
$$R_{\mathbf{n}, \Theta}\mathbf{a}=\mathbf{a}_{\perp}\cos\Theta+ (\mathbf{n}\times\mathbf{a})\sin\Theta+\mathbf{a}_{\|}$$
so applying this rule to our situation we find that the spiral converges to
$$\mathbf{r}_{\infty}=\mathbf{a}_{\perp}\frac{1-\lambda \cos\Theta}{1-2\lambda \cos\Theta+\cos^2\Theta}+ (\mathbf{n}\times\mathbf{a})\frac{\lambda \sin\Theta}{1-2\lambda \cos\Theta+\cos^2\Theta}+\frac{\mathbf{a}_{\|}}{1-\lambda}$$
where now everything can be computed only with knowledge of the formulas presented below for $\mathbf{n}, \cos\Theta$:
$$1+2\cos\Theta=c_1c_2+c_2c_3+c_1c_3-s_1s_2s_3$$
$$\mathbf{n}\propto (s_2 - c_1c_3 s_2 + s_1s_3~,~ c_3s_1-c_2s_1+c_1s_2s_3~,~ 1+s_1s_2c_3+c_1c_2-c_2c_3-c_1c_3)$$
$$R=R_x(\theta_1)R_y(\theta_2)R_z(\theta_3)$$
$$s_i=\sin\theta_i, c_i=\cos\theta_i$$
$$\mathbf{a}_{\|}=(\mathbf{n}\cdot\mathbf{a})\mathbf{n}~~,~~ \mathbf{a}_{\perp}=\mathbf{a}-\mathbf{a}_{\|}$$
Click on the spoiler below to reveal an extensive discussion about the original solution, which is more cumbersome, and derivations for the expressions above.

ORIGINAL: We can construct the sequence of vectors describing the center of the n-th cube by noting that the difference of two consecutive vectors is


 $$\mathbf{r}_n-\mathbf{r}_{n-1}=(\lambda R)^{n-1}\mathbf{a} ~~,~n>0~~,~\mathbf{r}_0=0$$

where $R$ is the rotation applied to the body, $\lambda$ the scale factor, and $\mathbf{a}$ the vector that you originally chose to translate with respect to. This recursion can be easily solved for any $n$ by adding all the equations up to the $n$-th one together and summing the geometric series formally:

$$\mathbf{r}_n=(1-\lambda R)^{-1}(1-(\lambda R)^n)\mathbf{a}$$

To evaluate this formal series we need access to the eigenvalues and eigenvectors of the rotation matrix. It is a well known fact that the eigenvalues are $\{1,e^{i\Theta}, e^{-i\Theta}\}$ with corresponding normalized eigenvectors denoted by $\{\mathbf{n},\mathbf{e}_+, \mathbf{e}_-\}$. The first eigenvector is the axis of the rotation, the unique vector with real components that remains invariant under applications of the rotation matrix. The other two eigenvectors are complex. We can decompose the vector in terms of those 3 vectors:

$$\mathbf{a}=x_n\mathbf{n}+ x_+\mathbf{e}_+ +x_- \mathbf{e}_-$$

This however solves the problem because it is true that under certain conditions that are fulfilled for this problem

$$f(R)\mathbf{n}=f(1)\mathbf{n}~~,~~ f(R)\mathbf{e}_{\pm}=f(e^{\pm i\Theta})\mathbf{e}_{\pm}$$

and thus we can write the solution to the problem formally as

$$\mathbf{r}_{n}=\frac{1-\lambda^n}{1-\lambda}x_n\mathbf{n}+ \frac{1-\lambda^n e^{in\Theta}}{1-\lambda e^{i\Theta}}x_+\mathbf{e}_+ +\frac{1-\lambda^n e^{-in\Theta}}{1-\lambda e^{-i\Theta}}x_-\mathbf{e}_- $$

which converges to


$$\mathbf{r}_{\infty}=\frac{1}{1-\lambda}x_n\mathbf{n}+ \frac{1}{1-\lambda e^{i\Theta}}x_+\mathbf{e}_+ +\frac{1}{1-\lambda e^{-i\Theta}}x_-\mathbf{e}_-$$

Now what we need is expressions for $\Theta, n, x_-, x_+$ in terms of the way you apply your rotation. From what I gathered from your simulation "Euler XYZ" angle mode, the system is applying first a rotation around the z-axis, then one around the y-axis, and finally one around the $x$-axis, with angles you can control. Then the representation of the rotation matrix in terms of the 3 rotation angles in the simulation is

$$R=R_x(\theta_1)R_y(\theta_2)R_z(\theta_3)={\begin{bmatrix}c_{2}c_{3}&-c_{2}s_{3}&s_{2}\\c_{1}s_{3}+c_{3}s_{1}s_{2}&c_{1}c_{3}-s_{1}s_{2}s_{3}&-c_{2}s_{1}\\s_{1}s_{3}-c_{1}c_{3}s_{2}&c_{3}s_{1}+c_{1}s_{2}s_{3}&c_{1}c_{2}\end{bmatrix}}$$

where $c_i=\cos\theta_i~,~ s_i=\sin\theta_i$. We can find the angle $\Theta$ of rotation around the axis $\mathbf{n}$ from the trace of the rotation matrix

$$1+2\cos\Theta=c_1c_2+c_2c_3+c_1c_3-s_1s_2s_3$$

and the axis of rotation is given by

$$\mathbf{n}\propto (s_2 - c_1c_3 s_2 + s_1s_3~,~ c_3s_1-c_2s_1+c_1s_2s_3~,~ 1+s_1s_2c_3+c_1c_2-c_2c_3-c_1c_3)$$

The other two vectors are given by

$$\mathbf{e}_+=(\mathbf{e}_-)^*\propto((e^{i\Theta}->!c_1c_3)s_2+s_1s_3~,~e^{i\Theta}c_2s_1-s_1c_3-c_1s_2s_3~,~c_1(c_2->!e^{i\Theta}c_3)+e^{i\Theta}(e^{i\Theta}-c_2c_3+s_1s_2s_3))$$

You can also compute the components of $\mathbf{a}$ in the basis of >!eigenvectors by inverting the Gram matrix.

Comment 1: It is not necessary to derive these relations explicitly, since all one has to do is numerically compute all the eigenvectors and eigenvalues for each set of angles desired and then the convergence point can be found on a case by case basis.

Comment 2: The spiral is indeed a logarithmic spiral in some sense. This can be made explicit for example by considering angles that can be written as $\Theta=\frac{2\pi \ell}{q}$, with $\ell, q$ integers with $gcd(\ell,q)=1$. Then note that

$$\frac{|\mathbf{r}_{n+2q}-\mathbf{r}_{n+q}|}{|\mathbf{r}_{n+q}-\mathbf{r}_{n}|}=\lambda^q$$

which means that the distances between adjacent windings of the spiral decays exponentially. However it is not exactly the same as a 2-D logarithmic spiral.

Comment 3: The effect that you are noticing for certain values of the angles is basically the 3-D spiral degenerating into an exact 2-D logarithmic spiral. This can be achieved by choosing your translation vector $\mathbf{a}$ such that any one of it's components $(x_n, x_+, x_-)$ is zero.  If two of it's components are zero, the spiral should degenerate into a straight line! If you have the rotation matrix, it is easy to engineer the translation vector so as to arrange the 2-D spiral to appear. The inverse problem is much more difficult to solve because of the non-linearity of the equations involved.

EDIT 1: To find the components $(x_n,x_+, x_-)$ one needs to define the complex inner product $\mathbf{a}^*\cdot \mathbf{b}\equiv a^H b=a^*_1b_1+a^*_2b_2+a^*_3b_3$. Then taking inner products of the equation defining the components of $\mathbf{a}$ shows that in order to find the components we need to solve the following system of equations

$$\begin{pmatrix}1&\mathbf{n}^*\cdot\mathbf{e_+}&\mathbf{n}^*\cdot\mathbf{e_-}\\\mathbf{e_+}^*\cdot\mathbf{n}&1&\mathbf{e_+}^*\cdot\mathbf{e_-}\\\mathbf{e_-}^*\cdot\mathbf{n}&\mathbf{e_-}^*\cdot\mathbf{e_+}&1\end{pmatrix}\begin{pmatrix}x_n\\x_+\\x_-\end{pmatrix}=\begin{pmatrix}\mathbf{n}^*\cdot\mathbf{a}\\\mathbf{e_+}^*\cdot\mathbf{a}\\\mathbf{e_-}^*\cdot\mathbf{a}\end{pmatrix}$$

This can be done numerically easily, by inverting the Gram matrix on the LHS.

Note: This formula works when the eigenvectors are normalized, the ones stated above have not been, so they have to be divided by their norm $||\mathbf{n}||=\sqrt{\mathbf{n^*}\cdot \mathbf{n}}$.

