Set of equations describing general triangle only numerically solvable using sympy

Question

I have two triangles sharing side c. Knowns are angles alpha,alpha' and side a and a'. I am looking for b and b'. From the set of 3 cosine rules I can find e.g. b' as shown in the code below. Unfortunately there are 2 solutions after every solve, from which I pick the correct one manually and in the end I have to solve numerically for the solution, which is correct.

Can I somehow deduce which solution after solving is the correct one or put some constraints on solve? Is there a closed solution to my problem?

from mpmath import radians
from sympy import N, symbols, cos, init_printing,Eq,solve,simplify,nsolve

init_printing()

alphaPrime, alpha, a, aPrime, b, bPrime,c   = symbols('αprime α a aprime b bprime,c')

eq1=Eq( a**2          , b**2+c**2-2*b*c*cos(alpha))
eq2=Eq( aPrime**2     , bPrime**2+c**2-2*bPrime*c*cos(alphaPrime))
eq3=Eq( (a+aPrime)**2 , bPrime**2+b**2-2*b*bPrime*cos(alpha+alphaPrime))

eq1c=solve(eq1,c,simplify=False)
eq2c=solve(eq2,c,simplify=False)
Eq1MinusEq2=eq1c[1]-eq2c[0]
eq1_eq2_b=solve(Eq1MinusEq2,b,simplify=False)

eq3_b=solve(eq3,b,simplify='False')
eq_1n2_Minus_eq3=eq1_eq2_b[1]-eq3_b[1]
display(simplify(eq_1n2_Minus_eq3))
expr=eq_1n2_Minus_eq3.subs([(aPrime,1.6547),(a,1.1),(alphaPrime,radians(27.031)), 
(alpha,radians(25.268))])
display(nsolve(expr,bPrime,1)) ##works

gives me: 3.45613973389693+1.86640158050948⋅10−8𝑖, which's real part is correct.

Answer 1

These are your equations:

In [156]: eq1
Out[156]: 
 2    2                   2
a  = b  - 2⋅b⋅c⋅cos(α) + c 

In [157]: eq2
Out[157]: 
  2     2                     2
a′  = b′  - 2⋅b′⋅c⋅cos(α′) + c 

In [158]: eq3
Out[158]: 
        2    2                          2
(a + a′)  = b  - 2⋅b⋅b′⋅cos(α + α′) + b′ 

We want to eliminate c and then solve for a and a' so let's eliminate c first:

In [159]: eq12 = resultant(eq1, eq2, c)

In [160]: eq12
Out[160]: 
 4      2   2      2  2      2                          2   2    2          2   2     4       2  2  ↪
a  - 2⋅a ⋅a′  - 2⋅a ⋅b  + 4⋅a ⋅b⋅b′⋅cos(α)⋅cos(α′) - 4⋅a ⋅b′ ⋅cos (α′) + 2⋅a ⋅b′  + a′  - 4⋅a′ ⋅b ⋅ ↪

↪    2          2  2       2                           2   2    4      3                        2   ↪
↪ cos (α) + 2⋅a′ ⋅b  + 4⋅a′ ⋅b⋅b′⋅cos(α)⋅cos(α′) - 2⋅a′ ⋅b′  + b  - 4⋅b ⋅b′⋅cos(α)⋅cos(α′) + 4⋅b ⋅b ↪

↪  2    2         2   2    2          2   2         3                    4
↪ ′ ⋅cos (α) + 4⋅b ⋅b′ ⋅cos (α′) - 2⋅b ⋅b′  - 4⋅b⋅b′ ⋅cos(α)⋅cos(α′) + b′

Now eq12 and eq3 are two equations for the two unknowns a and a'.

We can eliminate say a' to find an equation only for a. This equation is a quartic with two double roots that we can simplify down to nice expressions:

In [167]: eq_a = resultant(eq12, eq3, aPrime)

In [178]: r1, r2 = roots(eq_a, a)

In [179]: r2.expand(trig=True).simplify()
Out[179]: 
      ____________________________________________
     ╱  2 ⎛ 2                          2⎞    2    
    ╱  b ⋅⎝b  - 2⋅b⋅b′⋅cos(α + α′) + b′ ⎠⋅sin (α) 
   ╱   ────────────────────────────────────────── 
  ╱                                    2          
╲╱              (b⋅sin(α) + b′⋅sin(α′))           

In [180]: r1.expand(trig=True).simplify()
Out[180]: 
       ____________________________________________
      ╱  2 ⎛ 2                          2⎞    2    
     ╱  b ⋅⎝b  - 2⋅b⋅b′⋅cos(α + α′) + b′ ⎠⋅sin (α) 
-   ╱   ────────────────────────────────────────── 
   ╱                                    2          
 ╲╱              (b⋅sin(α) + b′⋅sin(α′))  

These are the two solutions for a. Corresponding solutions for a' can be found similarly and are the same apart from interchanging primes which only affects the final sin(α) part of the numerator.

The multiple solutions are to do with the equations being satisfied by negative angles and lengths. If we want to require everything to be positive then the formulae above can simplify because most factors under the square root are squares of positive quantities.

The key technique used above is polynomial resultants:

https://en.wikipedia.org/wiki/Resultant https://docs.sympy.org/latest/modules/polys/reference.html#sympy.polys.polytools.resultant

EDIT: I see that b and b' are the unknowns. Okay so adapted:

In [59]: eq_b = resultant(eq12, eq3, bPrime)

In [60]: from sympy.simplify.fu import TR5

In [61]: eq_b2 = TR5((eq_b).expand(trig=True)).expand()

In [64]: r1, r2 = roots(eq_b2, b)

In [74]: r1
Out[74]: 
       ___________________________________________________________________________________________________________________________________________________________________________________________________________________________________
      ╱                                                                        2         2 ⎛ 2    2                                             2    2   ⎞    2                                                                           
     ╱                                                                        a ⋅(a + a′) ⋅⎝a ⋅sin (α′) - 2⋅a⋅a′⋅sin(α)⋅sin(α′)⋅cos(α + α′) + a′ ⋅sin (α)⎠⋅sin (α′)                                                                       
-   ╱   ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── 
   ╱     4    4          3              3                      2   2    4       2          2   2    3       3                      2   2    2       4          2   2    2       2             3    3                            4    4    
 ╲╱     a ⋅sin (α′) - 4⋅a ⋅a′⋅sin(α)⋅sin (α′)⋅cos(α + α′) - 4⋅a ⋅a′ ⋅sin (α)⋅sin (α′) - 8⋅a ⋅a′ ⋅sin (α)⋅sin (α′)⋅cos(α + α′) - 4⋅a ⋅a′ ⋅sin (α)⋅sin (α′) + 6⋅a ⋅a′ ⋅sin (α)⋅sin (α′) - 4⋅a⋅a′ ⋅sin (α)⋅sin(α′)⋅cos(α + α′) + a′ ⋅sin (α) 

In [75]: r2
Out[75]: 
      ___________________________________________________________________________________________________________________________________________________________________________________________________________________________________
     ╱                                                                        2         2 ⎛ 2    2                                             2    2   ⎞    2                                                                           
    ╱                                                                        a ⋅(a + a′) ⋅⎝a ⋅sin (α′) - 2⋅a⋅a′⋅sin(α)⋅sin(α′)⋅cos(α + α′) + a′ ⋅sin (α)⎠⋅sin (α′)                                                                       
   ╱   ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── 
  ╱     4    4          3              3                      2   2    4       2          2   2    3       3                      2   2    2       4          2   2    2       2             3    3                            4    4    
╲╱     a ⋅sin (α′) - 4⋅a ⋅a′⋅sin(α)⋅sin (α′)⋅cos(α + α′) - 4⋅a ⋅a′ ⋅sin (α)⋅sin (α′) - 8⋅a ⋅a′ ⋅sin (α)⋅sin (α′)⋅cos(α + α′) - 4⋅a ⋅a′ ⋅sin (α)⋅sin (α′) + 6⋅a ⋅a′ ⋅sin (α)⋅sin (α′) - 4⋅a⋅a′ ⋅sin (α)⋅sin(α′)⋅cos(α + α′) + a′ ⋅sin (α) 

There might be ways to simplify those.