Cube root of two $\sqrt[3]2$ continued fraction I know there is a nice way of getting the continued fraction expansion of quadratic irrationals mainly because they recur after a point, and if they recur after a point they are quadratic irrationals. When constructing the expansion you can multiply by conjugates (kind of), e.g. 
$\sqrt 3 =1+\sqrt 3 -1 = 1+\frac {1}{\frac {\sqrt 3 +1}{2}} $
Where you use $(\sqrt 3 - 1)(\sqrt 3 +1)=2$.
Are there identities that would help with the construction for $ \sqrt[3]{2} $?
One I thought was useful in the first step to get [1; 3,...] was 
$ (\sqrt[3]{2}-1)( \sqrt[3]{4} + \sqrt[3]{2}+1 )=1$,
So you get:
$ \sqrt[3]{2}=1+( \sqrt[3]{2}-1 )=1+\frac {1}{ \sqrt[3]{4} + \sqrt[3]{2}+1 }= 1+\frac {1}{3+ (\sqrt[3]{4} + \sqrt[3]{2}-2)} $
Thanks for the help.
 A: Starting from the column vector $(1,0,0,-2)$, consider the following steps:
Step a) Repeat multiplication by the matrix $A$
$$A=\begin{bmatrix}
1&0&0&0\\
3&1&0&0\\
3&2&1&0\\
1&1&1&1
\end{bmatrix}$$
while the coefficients of the resulting vector have different signs.
Step b) Reverse the coefficients of the vector, or equivalently multiply by
$$B=\begin{bmatrix}
0&0&0&1\\
0&0&1&0\\
0&1&0&0\\
1&0&0&0
\end{bmatrix}$$
Then the number of times you multiply by $A$ in step a gives the partial quotients of continued fraction of $\sqrt[3]{2}$.

For, starting from $(1,0,0,-2)$, successive multiplication by $A$ gives:
\begin{align}
(1,0,0,-2)
&\xrightarrow A\color{red}{(1,3,3,-1)}\\
&\xrightarrow A(1,6,12,6)
\end{align}
hence in step a we multiply by $A$ one time only, because $(1,6,12,6)$ have positive coefficients only, hence the first partial quotient is $1$:
$$\sqrt[3]{2}=1+\cdots$$
Apply step b to $(1,3,3,-1)$ we get $(-1,3,3,1)$.
Then applying step a to $(-1,3,3,1)$, successive multiplication by $A$ gives:
\begin{align}
(-1,3,3,1)
&\xrightarrow A(-1,0,6,6)\\
&\xrightarrow A(-1,-3,3,11)\\
&\xrightarrow A\color{red}{(-1,-6,-6,10)}\\
&\xrightarrow A(-1,-9,-21,-3)\\
\end{align}
hence the second partial quotient is $3$:
$$\sqrt[3]{2}=1+\frac 1{3+}\cdots$$
and so on...

This algorithm holds for every algebraic number of third degree which is the only positive root of it minimal polynomial.
For higher degree the matrix $A$ is enlarged as in the Tartaglia-Pascal triangle; for example for fourth degree:
$$A=\begin{bmatrix}
1&0&0&0&0\\
4&1&0&0&0\\
6&3&1&0&0\\
4&3&2&1&0\\
1&1&1&1&1
\end{bmatrix}$$

For the intuition behind this algorithm.
Then vector $(1,0,0,-2)$ corresponds to the polynomial $p(x)=x^3-2$.
Multyplication by $A$ corresponds to $p(x)\mapsto p(x+1)$, while revesing in step b corresponds to $p(x)\mapsto x^3p(1/x)$.
Finally Descartes's signs rule provide the stopping criterion in step a.
A: The accepted answer looks like based on Vincent's continued fractions method (1836). Downside is it's inefficiency. Say, the root is at $0.000001$ so $a_0=0$. In order to calculate the next term $a_1$ you have to invert the polynomial and the root of the inverted polynomial appears at $1000000$. This means at the next stage you have to perform $p(x)↦p(x+1)$ translation 1000000 times.
Instead one should find the lowest bound $b_0$ of the positive roots and perform a $p_0(x)↦p_0(x+b_0+1) = p_1(x)$ translation. We should also know that the lowest bound is just a bound and most probably doesn't yield an exact figure that pinpoints the smallest positive root. This means, we may have to perform subsequent lowest bound attempts on the translated polynomials up until there is no sign change remains in the coefficients of the last obtained polynomial, a.k.a. $var(p_n(x)) = 0$ which means there is no positive root according to Descartes Rule of Sign Change. Then the $i^{th}$continued fraction coefficient $a_i$ would be the sum of all horizontal translation amounts less $1$ such as $a_i = b_{0,i}+b_{1,i}+...+b_{n,i}+n$. Obviously if the subsequent translations yield a no sign change at the polynomial coefficients and at the same time $p(0) = 0$ case occurs then the continued fraction ends here. If not then we need to perform $p(\frac{1}{x-1})$ translation (one step to the right and invert) which is two translations one after the other like $p(x-1)$ and $p(\frac{1}{x})$. Keep in mind that inverting a polynomial like $p(\frac{1}{x})$ is as simple as reversing the coefficients of the polynomial. Such as if $p(x) = x^4+7x^2+2x-5$ then $p(\frac{1}{x}) = -5x^4+2x^3+7x+1$ but we are only interested in the roots and should make the initial coefficient 1. Accordingly dividing all coefficients by $-5$ we can safely say that the $p(\frac{1}{x})$ of our interest would actually be $p(\frac{1}{x}) = x^4-\frac{2}{5}x^3-\frac{7}{5}x-\frac{1}{5}$.
Now there are two crucial parts to this problem.

*

*Efficiently and most precisely finding the lowest bound.

*Translating the polynomial by an arbitrary amount ($b_{lower}+1$) in it's extended form.

1. Lowest Positive Bound
There are several linear and quadratic algorithms to find the lower positive bound. Linear ones resolve faster but yielding a coarse lower bound while quadratic ones yield much finer bounds. Keep in mind that a sharper lower positive bound would minimize the translation count dramatically. According to my calculations on overall performance, I prefer a quadratic one. If you need to know more about them please drop a comment and i will extend the answer.
2. Translation of a Polynomial by an Arbitrary Amount
This is as simple as performing subsequent synthetic divisions, the degree of the polynomial many times. It's known as the Ruffini-Horner Method. A description of how it works can be found in one of my SO answers.
