Given $\color{blue}{t = 2\pi/p}$ for the appropriate prime $p=4m+1$.
I. Sine
$$\begin{align} \frac{5+\sqrt{5}}8 &=\sin^2(t)\\ \frac{13+\sqrt{13}}8 &=\sin^2(t)+\sin^2(3t)+\sin^2(4t)\\ \frac{17+\sqrt{17}}8 &=\sin^2(3t)+\sin^2(5t)+\sin^2(6t)+\sin^2(7t)\\ \frac{29+\sqrt{29}}8 &=\sum_{k=1}^7\sin^2(a_k\, t) \end{align}$$
with the seven $a_k = 1,4,5,6,7,9,13.$ And so on for other prime $p=4m+1.$ For the opposite sign, one uses the remaining integers $b_k \leq \frac{p-1}2.$ For example,
$$\frac{17-\sqrt{17}}8 =\sin^2(t)+\sin^2(2t)+\sin^2(4t)+\sin^2(8t)$$
where the $a_k$ are simply $2^n$. (The next Fermat prime $p=257$ isn't so nice since it has $256/4 = 64$ sine terms.)
II. Cosine
This uses the same set of multipliers $a_k$.
$$\begin{align} \frac{3-\sqrt{5}}8 &=\cos^2(t)\\ \frac{11-\sqrt{13}}8 &=\cos^2(t)+\cos^2(3t)+\cos^2(4t)\\ \frac{15-\sqrt{17}}8 &=\cos^2(3t)+\cos^2(5t)+\cos^2(6t)+\cos^2(7t)\\ \frac{27-\sqrt{29}}8 &=\sum_{k=1}^7\cos^2(a_k\, t) \end{align}$$
with the same seven $a_k = 1,4,5,6,7,9,13.$ For the opposite sign,
$$\frac{15+\sqrt{17}}8 =\cos^2(t)+\cos^2(2t)+\cos^2(4t)+\cos^2(8t)$$
III. Conclusion
Given prime $p=4m+1$ and $t = 2\pi/p.$ The pattern clearly is,
$$\frac{p\pm\sqrt{p}}8 = \sum_{k=1}^m \sin^2(a_k\, t)$$
$$\frac{(p-2)\mp\sqrt{p}}8 = \sum_{k=1}^m \cos^2(a_k\, t)$$
Question: I used Mathematica's integer relations to find the above examples. But, for any prime $p=4m+1$, what is a clever and faster algorithm to derive the correct set of $a_k$?