Here is how i would proceed. There will be some human calculation routine too, but the "technique" is a combined human+machine activity. The structure will come in a second, but so far, let us experiment something.
First of all, i asked sage about the Galois group of the number field generated by a root of the polynomial
$$
P = x^8 + x^7 - 7x^6 - 6x^5 + 15x^4 + 10x^3 - 10x^2 - 4x + 1
$$
to be factorized. It is a cyclic group of order eight. I also asked for the subfields, and sage mentions a subfield of degree four over $\Bbb Q$, related to the polynomial
$$
Q = x^4 + x^3 - 6 x^2 - x + 1
\ .
$$
(Code and results are postponed to not disturb here the flow. A structural way to see it follows below.)
Let us compute its roots. This is a "twisted reciprocal" polynomial, one can substitute $X=ix$ first to get a reciprocal polynomial. The substitution
$$
t = x-\frac 1x
$$
allows to relate $Q/x^2$ with
$$
T = t^2 + t -4\ .
$$
(This polynomial is also present in the sage list of subfields...)
The roots of $T$ are $t_{\pm}=(-1\pm \sqrt{17})/2$. The polynomials $x^2-t_\pm x-1$ have roots
$$
\begin{aligned}
x_{\boxed\pm}
&=\frac 12\left(t_\pm\boxed\pm\sqrt{t_\pm^2+4}\right)
\\
&=\frac 12\left(t_\pm\boxed\pm\sqrt{8-t_\pm}\right)
\\
&=
\frac 14\left(-1\pm\sqrt{17}\right)
\boxed\pm
\frac 12\sqrt{\frac 12\left(17\mp\sqrt{17}\right)}
\ .
\\
4x_{\boxed\pm}
&=
\left(-1\pm\sqrt{17}\right)
\boxed\pm
\sqrt{2\left(17\mp\sqrt{17}\right)}
\ .
\end{aligned}
$$
Time for the structure. I learned it from an older book of Hans Rademacher,
typed it for myself some years ago. The idea is to use Gaussian periods.
Let $K=\Bbb Q(w)$ be the cyclotomic field generated by the primitive root of order $17$, here fixed and denoted by $w$ instead of the longer $\zeta_{17}$, which is hard to type. The Galois group of $K$ over $\Bbb Q$ has elements sending $w$ to some other primitive root $w^k$, so $k\in 1,2,\dots 16$.
My choice is $3$ as a generator of the cyclic multiplicatie group $(\Bbb Z/17)^\times$, so let us use the Frobenius morphism
$$ Fw = w^3\ .$$
So $F^2w = FFw=Fw^3=(w^3)^3=w^9$, and $F^3w = FF^2w=Fw^9=(w^3)^9 = w^{27}$, and in general $F^jw = w^{3^j}$.
The idea of the Gaussian periods is simple. We associate the following sums:
$$
\begin{aligned}
S &= F^0w + Fw + F^2w + F^3w+\dots+F^{14}w+F^{15}w\\
&=w+w^2+w^3+\dots+w^{16}\text{ (in other order)}\\
&=-1\ ,\\[2mm]
&\qquad\text{ and we split it in two subsums,}\\[2mm]
S_0 &= F^0w + F^2w + F^4w+\dots+F^{14}w\\
S_1 &= Fw + F^3w + F^5w+\dots+F^{15}w\\[2mm]
&\qquad\text{ and we split each in two subsums,}\\[2mm]
S_{00} &= F^0w + F^4w + F^8w+F^{12}w\\
S_{10} &= F^2w + F^6w + F^{10}w+F^{14}w\\[2mm]
S_{01} &= Fw + F^5w + F^9w+\dots+F^{13}w\\
S_{11} &= F^3w + F^7w + F^{11}w+\dots+F^{15}w\\[2mm]
&\qquad\text{ and we split each in two subsums,}\\[2mm]
S_{000}&= F^0w + F^8w\\
S_{100}&= F^4w + F^{12}w\\
&\qquad\text{and so on...}
\end{aligned}
$$
The indices are my choice, they suggest binary written numbers. So $S_{01}$ is a sum of terms $F^jw$ with $j$ congruent to (the binary phone number) $01_{[2]}=1$ modulo $4$, $S_{101}=F^5w+F^{13}w$ and so on.
There are obvious splitting relations, for instance $S_0+S_1=S=-1$, $S_{00}+S_{10}=S_0$, sum relations. But there are also product relations.
Instead of writing them here explicitly, a hard jax job, i will use sage to "type" them. Sage is a mathematically structurally thinking language, so i hope the dry exposition gives the information in a human way.
sage: K.<w> = CyclotomicField(17)
sage: K
Cyclotomic Field of order 17 and degree 16
sage: w^17
1
sage: S = sum( [ w^(3^j) for j in [0..15] ] )
sage: S0 = sum( [ w^(3^j) for j in [0..15] if j%2 == 0 ] )
sage: S1 = sum( [ w^(3^j) for j in [0..15] if j%2 == 1 ] )
sage: S00 = sum( [ w^(3^j) for j in [0..15] if j%4 == 0 ] )
sage: S01 = sum( [ w^(3^j) for j in [0..15] if j%4 == 1 ] )
sage: S10 = sum( [ w^(3^j) for j in [0..15] if j%4 == 2 ] )
sage: S11 = sum( [ w^(3^j) for j in [0..15] if j%4 == 3 ] )
sage: S000 = sum( [ w^(3^j) for j in [0..15] if j%8 == 0 ] )
sage: S001 = sum( [ w^(3^j) for j in [0..15] if j%8 == 1 ] )
sage: S010 = sum( [ w^(3^j) for j in [0..15] if j%8 == 2 ] )
sage: S011 = sum( [ w^(3^j) for j in [0..15] if j%8 == 3 ] )
sage: S100 = sum( [ w^(3^j) for j in [0..15] if j%8 == 4 ] )
sage: S101 = sum( [ w^(3^j) for j in [0..15] if j%8 == 5 ] )
sage: S110 = sum( [ w^(3^j) for j in [0..15] if j%8 == 6 ] )
sage: S111 = sum( [ w^(3^j) for j in [0..15] if j%8 == 7 ] )
sage: S0+S1, S0*S1
(-1, -4)
sage: S0.minpoly()
x^2 + x - 4
sage: S00*S10, S01*S11
(-1, -1)
sage: S00.minpoly()
x^4 + x^3 - 6*x^2 - x + 1
sage: S10.minpoly()
x^4 + x^3 - 6*x^2 - x + 1
sage: S01.minpoly()
x^4 + x^3 - 6*x^2 - x + 1
sage: S11.minpoly()
x^4 + x^3 - 6*x^2 - x + 1
sage: S000*S100
w^14 + w^12 + w^5 + w^3
sage: S00, S10, S01, S11
(-w^15 - w^14 - w^12 - w^11 - w^10 - w^9 - w^8 - w^7 - w^6 - w^5 - w^3 - w^2 - 1,
w^15 + w^9 + w^8 + w^2,
w^14 + w^12 + w^5 + w^3,
w^11 + w^10 + w^7 + w^6)
sage: S000*S100 == S01
True
(The first long expression staring with -w^15 - w^14 ... is only the computer version to hide the "missing terms" of the minimal polynomial of $w$, the four that not appear in the list, since it thinks it is a good idea to not use $w^{16}$...)
Note that in the "naive approach" from the beginning was also dealing with polynomials that appear as minimal polynomials of the one or the other $S_?$.
The OP wants to understand the splitting of the polynomial $P$, which is the minimal polynomial (over rationals) of
$$S_{000}=w+w^{16}=w+\frac 1w\ .$$
And now we can give four factors of $P$. With the computer aid of sage:
sage: R.<x> = PolynomialRing(K)
sage: R
Univariate Polynomial Ring in x over Cyclotomic Field of order 17 and degree 16
sage: (x-S000)*(x-S100)
x^2 + (w^15 + w^14 + w^12 + w^11 + w^10 + w^9 + w^8 + w^7 + w^6 + w^5 + w^3 + w^2 + 1)*x + w^14 + w^12 + w^5 + w^3
sage: (x-S000)*(x-S100) == x^2 - S00*x + S01
True
sage: (x-S010)*(x-S110)
x^2 + (-w^15 - w^9 - w^8 - w^2)*x + w^11 + w^10 + w^7 + w^6
sage: (x-S010)*(x-S110) == x^2 - S10*x + S11
True
sage: (x-S001)*(x-S101)
x^2 + (-w^14 - w^12 - w^5 - w^3)*x + w^15 + w^9 + w^8 + w^2
sage: (x-S001)*(x-S101) == x^2 - S01*x + S10
True
sage: (x-S011)*(x-S111)
x^2 + (-w^11 - w^10 - w^7 - w^6)*x - w^15 - w^14 - w^12 - w^11 - w^10 - w^9 - w^8 - w^7 - w^6 - w^5 - w^3 - w^2 - 1
sage: (x-S011)*(x-S111) == x^2 - S11*x + S00
True
sage: (x-S00)*(x-S10)*(x-S01)*(x-S11)
x^4 + x^3 - 6*x^2 - x + 1
And here explicitly, because it is structural (the "traces", coefficients in $x$ on the RHSs are easy, we splitted them as the idea to proceed, for the computation of the free coefficients there is also a Galois pattern...):
$$
\begin{aligned}
(x-S_{000})(x-S_{100}) &= x^2 - S_{00}x+ S_{01}\ ,\\
(x-S_{010})(x-S_{110}) &= x^2 - S_{10}x+ S_{11}\ ,\\
(x-S_{001})(x-S_{101}) &= x^2 - S_{01}x+ S_{10}\ ,\\
(x-S_{011})(x-S_{111}) &= x^2 - S_{11}x+ S_{00}\ ,\\[2mm]
&\text{where}\\
(x-S_{00})
(x-S_{10})
(x-S_{01})
(x-S_{11})
&=x^4 + x^3 - 6x^2 - x + 1\ ,
\end{aligned}
$$
and its roots were computed in advance.
I have to submit, but there remains a lot to say, maybe just mention that this structure made the 19 years old Gauss to choose mathematics, not physics, one main reason for still missing a unification of matter in this field...
Postponed:
The sage code used to obtain the polynomial $Q$, the Galois theoretical information, and the decomposition of $P$.
sage: R.<x> = PolynomialRing( QQ )
sage: P = x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
sage: K.<v> = NumberField( P )
sage: v.minpoly()
x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
sage: G = K.galois_group()
sage: G
Galois group of Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
sage: G.is_cyclic()
True
sage: K.subfields()
[
(Number Field in v0 with defining polynomial x + 1, Ring morphism:
From: Number Field in v0 with defining polynomial x + 1
To: Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
Defn: -1 |--> -1, None),
(Number Field in v1 with defining polynomial x^2 + x - 4, Ring morphism:
From: Number Field in v1 with defining polynomial x^2 + x - 4
To: Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
Defn: v1 |--> -v^7 - v^6 + 6*v^5 + 6*v^4 - 10*v^3 - 9*v^2 + 5*v + 1, None),
(Number Field in v2 with defining polynomial x^4 + x^3 - 6*x^2 - x + 1, Ring morphism:
From: Number Field in v2 with defining polynomial x^4 + x^3 - 6*x^2 - x + 1
To: Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
Defn: v2 |--> v^4 - 4*v^2 + v + 2, None),
(Number Field in v3 with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1, Ring morphism:
From: Number Field in v3 with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
To: Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
Defn: v3 |--> v, Ring morphism:
From: Number Field in v with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
To: Number Field in v3 with defining polynomial x^8 + x^7 - 7*x^6 - 6*x^5 + 15*x^4 + 10*x^3 - 10*x^2 - 4*x + 1
Defn: v |--> v3)
]
sage: K2, morphism, _ = K.subfields()[2]
sage: K2
Number Field in v2 with defining polynomial x^4 + x^3 - 6*x^2 - x + 1
sage: K2.gen()
v2
sage: for f, mul in P.base_extend(K2).factor():
....: print f
....:
x^2 - v2*x - 1/2*v2^3 + 3*v2 - 3/2
x^2 + (-v2^3 - v2^2 + 6*v2 + 1)*x - 1/2*v2^3 - v2^2 + 2*v2 + 3/2
x^2 + (1/2*v2^3 - 3*v2 + 3/2)*x + v2^3 + v2^2 - 6*v2 - 1
x^2 + (1/2*v2^3 + v2^2 - 2*v2 - 3/2)*x + v2
Best Answer
There are several ways to do it. One straightforward way is to factor the octic into two quartics over a square root extension, then use the quartic formula. Quartics have a cubic resolvent, but since the $17$th root of unity is construcible, then we know this cubic will also factor over a square root.
I. Factoring the octic
Given,
$$x^8+x^7-7x^6-6x^5+15x^4+10x^3-10x^2-4x+1= 0$$
Factor this into two quartics over $\sqrt{17}$ choosing,
$$x^4 + \frac{1 - \sqrt{17}}2x^3 - \frac{3 + \sqrt{17}}2x^2 + (2 + \sqrt{17})x - 1 = 0$$
such that $x = 2\cos \left(\frac{2\pi}{17} \right)$. Depress the quartic with the transformation,
$$x = y+\frac{-1 + \sqrt{17}}8$$
to get,
$$y^4 + \frac{-51 - 5\sqrt{17}}{16}y^2 + \frac{17 + 7\sqrt{17}}{16}y + \frac{221 + 43 \sqrt{17}}{512} = 0$$
II. Depressed quartic
Given,
$$y^4+py^2+qy+r = 0$$
then its cubic resolvent is simply,
$$ u^3 + 2p u^2 + (p^2 - 4r)u - q^2 =0$$
With any non-zero root $u$, then the depressed quartic's formula for the four roots are,
$$y_k = \frac12\left(\sqrt{u}\pm\sqrt{-(2p+u)-\frac{2q}{\sqrt{u}}}\right)$$
for $\pm\sqrt{u}$ and you are done.
III. Cubic resolvent
Given $(p,q,r)$ we find the resolvent as,
$$u^3 + \frac{-51 - 5 \sqrt{17}}8u^2 + \frac{323 + 53 \sqrt{17}}{32}u - \left(\frac{17 + 7\sqrt{17}}{16}\right)^2 = 0$$
this has the root,
$$u = \frac{17-\sqrt{17}}{8}$$
Using the formula and reversing the transformation, we find,
$$x = \frac{-1+\sqrt{17}}8+\frac12\sqrt{u}+\frac12\sqrt{\sqrt{17}+2u-\sqrt{\frac{17}{u}}-\sqrt{u}}$$
After some minor algebraic manipulation, we end up with the original form,
$$\small{x_1 = 2\cos \left(\frac{2\pi}{17} \right) = \frac{1}{8} \left( -1 + \sqrt{17} + \sqrt{ 2 \left(17 - \sqrt{17} \right)} + 2 \sqrt{ 17 + 3 \sqrt{17} - \sqrt{2 \left(17- \sqrt{17} \right)} - 2 \sqrt{2 \left(17+ \sqrt{17} \right)} } \right)}$$
IV. Conjugate
Using one of its conjugates,
$$\small{x_2 = 2\cos \left(\frac{8\pi}{17} \right) = \frac{1}{8} \left( -1 + \sqrt{17} + \sqrt{ 2 \left(17 - \sqrt{17} \right)} \color{red}{-} 2 \sqrt{ 17 + 3 \sqrt{17} - \sqrt{2 \left(17- \sqrt{17} \right)} - 2 \sqrt{2 \left(17+ \sqrt{17} \right)} } \right)}$$
then,
$$t=x_1+x_2 = 2\cos \left(\frac{2\pi}{17} \right)+2\cos \left(\frac{8\pi}{17} \right) = \frac{1}{4} \left( -1 + \sqrt{17} + \sqrt{ 2 \left(17 - \sqrt{17} \right)}\right)$$
a root of,
$$t^4 + t^3 - 6t^2 - t + 1 = 0$$
which was the second question of the OP.