Ask Your Question
2

Large System of Quadratic Trigonometric Equations

asked 2023-11-05 16:30:49 +0100

aobara gravatar image

updated 2023-11-11 17:14:40 +0100

Following the answers provided in 2 previous posts https://ask.sagemath.org/question/706... and https://ask.sagemath.org/question/726..., we tried to apply the proposed methods to a yet more complex system of trigonometric equations.

We are now considering a system with one distance variable $l$, and 10 angle variables $\alpha$, $\beta$, $\gamma$, $\delta$, $\rho$, $a$, $b$, $c$, $d$ and $e$. The said trigonometric system of equations is consistuted of 11 equations as described below:

$$ 0 = l \left( {4} \cos{\left({{\alpha}}\right)} + \cos{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} + \cos{\left({{\beta}}\right)} + \cos{\left({{\alpha} - \frac{{2} {\pi}}{{3}}}\right)} \right) - {1} $$ $$ 0 = {4} \sin{\left({{\alpha}}\right)} + \sin{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} + \sin{\left({{\beta}}\right)} + \sin{\left({{\alpha} - \frac{{2} {\pi}}{{3}}}\right)} $$ $$ 0 = {2} \cos{\left({{\alpha}}\right)} - \cos{\left({{\gamma}}\right)} - \cos{\left({{\delta}}\right)} - \cos{\left({{\rho}}\right)} + \cos{\left({{\alpha} + \frac{{2} {\pi}}{{3}}}\right)} $$ $$ 0 = l \left( {2} \sin{\left({{\alpha}}\right)} - \sin{\left({{\gamma}}\right)} - \sin{\left({{\delta}}\right)} - \sin{\left({{\rho}}\right)} + \sin{\left({{\alpha} + \frac{{2} {\pi}}{{3}}}\right)} \right) - {1} $$ $$ 0 = \cos{\left({a}\right)} - \cos{\left({b}\right)} - \cos{\left({c}\right)} + \cos{\left({{\gamma}}\right)} + \cos{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} $$ $$ 0 = \sin{\left({a}\right)} - \sin{\left({b}\right)} - \sin{\left({c}\right)} + \sin{\left({{\gamma}}\right)} + \sin{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} $$ $$ 0 = \cos{\left({c}\right)} + \cos{\left({e}\right)} - \cos{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} - \cos{\left({{\alpha}}\right)} + \cos{\left({{\delta}}\right)} $$ $$ 0 = \sin{\left({c}\right)} + \sin{\left({e}\right)} - \sin{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} - \sin{\left({{\alpha}}\right)} + \sin{\left({{\delta}}\right)} $$ $$ 0 = \cos{\left({b}\right)} + \cos{\left({d}\right)} - \cos{\left({{\rho}}\right)} - \cos{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} - \cos{\left({e}\right)} $$ $$ 0 = \sin{\left({b}\right)} + \sin{\left({d}\right)} - \sin{\left({{\rho}}\right)} - \sin{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} - \sin{\left({e}\right)} $$ $$ 0 = \cos{\left({{\rho}}\right)} + \cos{\left({{\alpha} - \frac{{\pi}}{{3}}}\right)} + \cos{\left({{\alpha}}\right)} - \cos{\left({{\gamma}}\right)} - \cos{\left({{\delta}}\right)} - \cos{\left({d}\right)} - \cos{\left({a}\right)} $$

as usual, we are interested in finding the minimal polynomial of the root $l$ between 0 and 1 for which a numerical appriximation is $l \approx 0.2226926944766917$.

Below is my (unsuccessfull) attempt to use to solve this problem:

#numerical values
l = 0.2226926944766917
alpha = 0.6331728275782392
beta = -1.3273875824789818
gamma = -1.4524093897328
delta = -1.226349124695996
rho = -1.3273875824743633
a = 0.9747473834353503
b = -0.9033592481852397
c = 0.2194297013449713
d = 0.33163441038029523
e = 1.1505721363111867

# conversion to non trig variables

c0 = cos(alpha)
s0 = sin(alpha)
c1 = cos(beta)
s1 = sin(beta)
c2 = cos(gamma)
s2 = sin(gamma)
c3 = cos(delta)
s3 = sin(delta)
c4 = cos(rho)
s4 = sin(rho)
c5 = cos(a)
s5 = sin(a)
c6 = cos(b)
s6 = sin(b)
c7 = cos(c)
s7 = sin(c)
c8 = cos(d)
s8 = sin(d)
c9 = cos(e)
s9 = sin(e)

# system of equations
eq1 = 4*c0+(c0+sqrt(3)*s0)/2+c1+(-c0+sqrt(3)*s0)/2-1/l
eq2 = 4*s0+(-sqrt(3)*c0+s0)/2+s1+(-sqrt(3)*c0-s0)/2
eq3 = 2*c0-c2-c3-c4-(c0+sqrt(3)*s0)/2
eq4 = 2*s0-s2-s3-s4+(sqrt(3)*c0-s0)/2-1/l
eq5 = c5-c6-c7+c2+(c0+sqrt(3)*s0)/2
eq6 = s5-s6-s7+s2+(-sqrt(3)*c0+s0)/2
eq7 = c7+c9-(c0+sqrt(3)*s0)/2-c0+c3
eq8 = s7+s9-(-sqrt(3)*c0+s0)/2-s0+s3
eq9  = c6+c8-c4-(c0+sqrt(3)*s0)/2-c9
eq10 = s6+s8-s4-(-sqrt(3)*c0+s0)/2-s9
eq11 = c4+(c0+sqrt(3)*s0)/2+c0-c2-c3-c8-c5




#numerical check
show(eq1.numerical_approx())
show(eq2.numerical_approx())
show(eq3.numerical_approx())
show(eq4.numerical_approx())
show(eq5.numerical_approx())
show(eq6.numerical_approx())
show(eq7.numerical_approx())
show(eq8.numerical_approx())
show(eq9.numerical_approx())
show(eq10.numerical_approx())
show(eq11.numerical_approx())


reset()
from time import time as stime
RR.<l, c0, c1, c2, c3, c4, c5, c6, c7, c8, c9, s0, s1, s2, s3, s4, s5, s6, s7, s8, s9, str3 >=AA[]
Sys = [ # linear equations
    l*(4*c0+(c0+str3*s0)/2+c1+(-c0+str3*s0)/2)-1,
    4*s0+(-str3*c0+s0)/2+s1+(-str3*c0-s0)/2,
    2*c0-c2-c3-c4-(c0+str3*s0)/2,
    l*(2*s0-s2-s3-s4+(str3*c0-s0)/2)-1,
    c5-c6-c7+c2+(c0+str3*s0)/2,
    s5-s6-s7+s2+(-str3*c0+s0)/2,
    c7+c9-(c0+str3*s0)/2-c0+c3,
    s7+s9-(-str3*c0+s0)/2-s0+s3,
    c6+c8-c4-(c0+str3*s0)/2-c9,
    s6+s8-s4-(-str3*c0+s0)/2-s9,
    c4+(c0+str3*s0)/2+c0-c2-c3-c8-c5,
    # Trigonometrics
    str3^2-3,
    s0^2 + c0^2 -1,
    s1^2 + c1^2 -1, 
    s2^2 + c2^2 -1, 
    s3^2 + c3^2 -1, 
    s4^2 + c4^2 -1, 
    s5^2 + c5^2 -1, 
    s6^2 + c6^2 -1, 
    s7^2 + c7^2 -1, 
    s8^2 + c8^2 -1, 
    s9^2 + c9^2 -1]
J1 = RR.ideal(Sys)
print("Solution dimension = ", J1.dimension())
t0 = stime()
Sols = J1.variety()
t1 = stime()
print("Number of solutions = ", len(Sols), ", found in %f6.3 seconds."%(t1-t0))
Errs = [abs(s[l]-0.22269269447669) for s in Sols]
Sol = Sols[Errs.index(min(Errs))]
show(Sols)
t2 = stime()
print(Sol[r])
P = Sol[l].minpoly()
t3 = stime()
print("Minimal polynomial (found in %f6.3 seconds) :"%(t3-t2))
print(P)

Permalink for testing purposes:

https://sagecell.sagemath.org/?z=eJx9...

Additionnal Permalink using Elimination_ideal, also fails:

https://sagecell.sagemath.org/?z=eJy1...

In this case sagemathcell websocket kernell crashes at the following line J1 = RR.ideal(Sys). Because no error message is returned we have no idea what is going wrong. It could yet again be a mistake on our side. We have included all numerical checks of the system's equations in order to facilitate further analysis. We are aware that this particular system is now much larger than the ones in previous posts, but do not know if it does not work because of a mistake on our side, or if we are pushing beyond the limits of sagemath?

In advance, thanks for your help.

kr,

edit retag flag offensive close merge delete

Comments

1

J1 = RR.ideal(Sys) works fine at Sagecell. Perhaps, it's time to update your local Sage.

Max Alekseyev gravatar imageMax Alekseyev ( 2023-11-05 17:57:48 +0100 )edit

Oh ok but, actually I do not have a local Sage installation. I am using the online version on https://sagecell.sagemath.org/ ...

aobara gravatar imageaobara ( 2023-11-05 18:10:31 +0100 )edit

I do not see the crash at Sagecell.

Max Alekseyev gravatar imageMax Alekseyev ( 2023-11-05 19:09:00 +0100 )edit

Here is what I see in the browser network console

Kernel: kernel_starting (null)
2kernel.js:107 Kernel: kernel_created (0d2fb748-7165-4a01-9315-4550c765d2fb)
kernel.js:464 Starting WebSockets: 
wss://sagecell.sagemath.orghttps://sagecell.sagemath.org/kernel/0d2fb748-7165-4a01-9315- 
4550c765d2fb
2kernel.js:107 Kernel: kernel_connected (0d2fb748-7165-4a01-9315-4550c765d2fb)
2kernel.js:107 Kernel: kernel_dead (0d2fb748-7165-4a01-9315-4550c765d2fb)

No results returned .

aobara gravatar imageaobara ( 2023-11-05 19:19:10 +0100 )edit

Permalink added in the Question Edit Box under the Sage Code Box

aobara gravatar imageaobara ( 2023-11-05 19:23:03 +0100 )edit

I used you link. Sagecell goes into long computation of dimension, but the ideal is constructed without an issue. No crash.

Max Alekseyev gravatar imageMax Alekseyev ( 2023-11-05 19:27:14 +0100 )edit
1

I also used this link. It starts well, prints the nuerical values of the quantities defining the problem, displaying a rotating image signalling a computation in progress. Then, after about two minutes, the rotating image disappears.

Crash ?

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2023-11-07 08:57:30 +0100 )edit
1

2 minutes is about the time limit for Sagecell. I would not assume a crash but just a timeout. At very least, J1 = RR.ideal(Sys) works fine, which can be seen via adding a print statement after it.

Max Alekseyev gravatar imageMax Alekseyev ( 2023-11-07 14:07:51 +0100 )edit

I don't understand. There is a print statement right after that print("Solution dimension = ", J1.dimension()), however it never prints anything. Why ? if J1 = RR.ideal(Sys) works fine why does the rest of the programm not execute afterwards?

aobara gravatar imageaobara ( 2023-11-07 14:44:16 +0100 )edit
1

Indeed SageMathCell has a time limit for computations. That print statement is not executed because J1.dimension() is a big (Groebner basis) computation that exceeds the timeout. I am running (a slight modification of) your script on a supercomputer with a 7 day time limit. You could also try alternatives like using msolve.

rburing gravatar imagerburing ( 2023-11-07 15:13:04 +0100 )edit

1 Answer

Sort by ยป oldest newest most voted
2

answered 2024-01-05 14:47:13 +0100

aobara gravatar image

updated 2024-01-05 16:31:39 +0100

I finally got to a result using a completely different approach. First I noticed $\rho$ and $\beta$ are one and the same variable as numerical approximative values suggest. That brings down the number of varialbles from 11 to 10. My sincere appologies for not picking that up earlier. Then, rewriting the set of equations as a linear system over the complex field yields:

image description

Then, inverting the matrix in the left hand side yields:

image description

This way, we can eliminate the five variables $\beta$, $\gamma$, $c$, $d$ and $e$, by computing the norms of the right hand sides and the left hand sides. The norms on the left hand sides are of course equal to one, which gets rid of 5 variables out of 10. We are now left with 5 equations of the 5 variables $l, \alpha, \delta, a$ and $b$:

image description

At this point, I tried to use Groebner elimination hoping that having less variables would make it possible to solve before timeout. However, this did not succeed so I went for a different approach.

With the first 3 equations we can linearize any sines and cosines powers of $\alpha$, $\delta$ and $a$. With the last 2 equations we can solve for $\cos(b)$ and $\sin(b)$ and elimiminate variable $b$ writing $\cos(b)^2+\sin(b)^2=1$. This yields a polynomial expression in $l$, as well as powers of sines and cosines of the variables $\alpha$, $\delta$ and $a$.

We then proceed to linearize this expression using the first 3 equations in the system until we are only left with an equation that is polynomial in $l$ and linear in $cos(\alpha)$, $cos(a)$ and $cos(\delta)$.

Since the expression is linear in $cos(\delta)$, we solve for $cos(\delta)$ and substitute into equation (2). This leaves us with a rational equation in $cos(a)$ and $cos(\alpha)$, which we linearize just as earlier in order to obtain a linear equation in $cos(a)$ and $cos(\alpha)$.

We repeat that process in order to eliminate $cos(a)$ and $cos(\alpha)$.

As a result, we obtain a 136 degree polynomial in $l$, for which our numerical l value is a root.

The result is accessible in the following Permalink:

Final "l" Polynomial

We still have a final question though... Can we be absolutely certain that this is the minimum polynomial for $l$? Or could we still factor out a smaller degree minimum polynomial?

Thanks for your help!!

edit flag offensive delete link more

Comments

1

@aobara :

Then, rewriting the set of equations as a linear system over the complex field yields [ ... ]

How do you pull this from your original set of equations ?

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2024-01-06 12:22:04 +0100 )edit

As an example take equation (1) and (2) from the original set of equations. Add (1) + l i * (2) yields:

$$ 0 = l (4 e^{i \alpha} + e^{i \left( \alpha-\frac{\pi}{3}\right)} +e^{i \beta} + e^{i \left( \alpha-\frac{2\pi}{3}\right)}) - 1 $$

which is equivalent to the first equation of the system in matrix form.

aobara gravatar imageaobara ( 2024-01-06 16:17:55 +0100 )edit
1

In your code, abs(r[0]-0.3) < 0.1 is satisfied by 4 different roots, and you are taking just the first of them. You need a better estimate to make sure it's the one.

Max Alekseyev gravatar imageMax Alekseyev ( 2024-01-08 12:40:42 +0100 )edit

Yes indeed! Thanks for pointing that out. Luckily it did return the correct root but by chance I guess. Did not notice these 4 roots were so close.

aobara gravatar imageaobara ( 2024-01-08 18:18:01 +0100 )edit
1

@aobara : Indeed..Nice ! I didn't view that...

Thank you very much !

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2024-01-08 18:53:31 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2023-11-05 16:30:49 +0100

Seen: 265 times

Last updated: Jan 05