I am trying to find the exact radicals of the sine and cosine of (m/n)*pi. I have been using sympy to do this and made this pretty good script:
from sympy import *
import sympy
from IPython.display import display
from sympy.abc import k, m, n, x
import math
from sympy import init_printing
init_printing()
print("You are solving for sin and cos of (m/n)pi")
m = int(input("Insert m: "))
n = int(input("Insert n: "))
gcd = sympy.gcd(m,n) #Getting the GCD and dividing by it will make the algorithm have less steps
m = m/gcd
n = n/gcd
prime = primefactors(n) #this outputs an array of all of the factors of n
cosp = sympify(-1) #this is necessary, we set the value of cos(pi) = -1, and format it for sympy
r = int(1)
for p in prime:#the nested loops will go through every factor. The for loop each of the factors,
while n%p == 0:#the while loop for the repetition of the factors
r = int(r*p)#the r variable serves as a comparison for the values of the cosine
n = int(n/p)#dividing by the prime we're checking for will complete the loop
if p == 2:
cosp = sqrt((1+cosp)/2)#simple bisection of the cosine formula
else:
cosp = (Sum((-1)**k*factorial(p)/(factorial(p-2*k)*factorial(2*k))*x**(p-2*k)*(1-x**2)**k , (k, 0, int((p-1)/2)))).doit()-cosp #check out this explanation: https://brilliant.org/wiki/expansions-of-certain-trigonometric-functions/. This is a version of the cosine formula, where x symbolizes cosine pi, and we bring the cosine value to the other side to get "0 = F(x)"
expandeq = expand(cosp)
res = solve(expandeq) #we solve for the roots of the equation
for root in res:
if N(sympy.cos(sympy.pi/r), 15) == re(N(root, 15)): #to remove all extra solutions, we check the approximation of each root and verify it is equal to cos(pi/r)
if str(root)[0] != "C":
cosp = root
sinp = expand(sqrt(1-(cosp**2)))
cosm = expand(Sum((-1)**k*factorial(m)/(factorial(m-(2*k))*factorial(2*k))*cosp**(m-2*k)*sinp**(2*k) , (k, 0, int(m/2))).doit(), trig = False) #finally, we use the formula again to multiply the cosine of pi/n by m
sinm = expand(sympy.sqrt(1-(cosm**2)), trig=False)
print("Approximation of sine of pi*",m, "/",r,": " , round(math.sin(math.pi*m/r), 15))
display(sinm)
print("Approximation of of cosine pi*",m, "/",r,": " , round(math.cos(math.pi*m/r), 15))
display(cosm)
and then opened this issue, highlighting the limits of the script, such as: calculating for the sin of m/25, or simply struggling with factors of numbers greater than 3, as it produces quintic (or higher) solvable equation. https://stackoverflow.com/questions/77970382/using-sympy-to-calculate-the-exact-values-of-sin-and-cos-of-m-npi
Here I want to ask about the mathematics side of this: every polynomial will be cyclotomic (https://en.m.wikipedia.org/wiki/Root_of_unity#Cyclotomic_polynomials). Is there a standard way to solve these? Any way to implement it into the code? Or any other information that may be useful (or interesting)? Thank you!