Wrong solution to a quadratic equation-CodePudding

I'm facing a mathematical/programming problem which I don't really understand. In the following α,β,θc are given parameters, ϕ is the variable. I have to solve the equation Cp=0, where:

and:

The solutions that I have computed are:

Coding the above solutions with Python/Numpy, I find that:

If β=0, then both solutions are correct: Cp(ϕ1)=Cp(ϕ2)=0.
However, if β!=0 then Cp(ϕ1)!=0 (which is wrong) and Cp(ϕ2)=0. Why? I checked again and again the signs and terms on my coded expressions, they seem fine to me...

Here is the code:

import numpy as np
theta_c = np.deg2rad(10)
alpha = np.deg2rad([1, 2, 4, 6, 8, 10, 12, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 85])
beta = np.deg2rad(15)

lamb = np.cos(alpha) * np.cos(beta)
nu = -np.sin(beta)
omega = np.sin(alpha) * np.cos(beta)

common_1 = lamb * nu * np.tan(theta_c)
common_2 = omega * np.sqrt(omega**2   nu**2 - lamb**2 * np.tan(theta_c)**2)
common_3 = omega**2   nu**2

phi_1 = np.arcsin((-common_1   common_2) / common_3)
phi_2 = np.arcsin((-common_1 - common_2) / common_3)

Cp = lambda phi: lamb * np.sin(theta_c)   nu * np.cos(theta_c) * np.sin(phi) - omega * np.cos(phi) * np.cos(theta_c)

print(Cp(phi_1))
# OUTPUT: this is wrong!
# [-0.02357418 -0.04432021 -0.0779439  -0.10278244 -0.12111675 -0.13494261
#  -0.1457658  -0.15858164 -0.17562354 -0.19112922 -0.20703468 -0.22411243
#  -0.24265822 -0.26274848 -0.28434621 -0.3073487  -0.33161092 -0.35695785
#  -0.38319138 -0.41009444 -0.43743318 -0.4649581 ]

print(Cp(phi_2))
# OUTPUT: all values are very close to zero. Good solution.
# [-1.04083409e-17 -4.16333634e-17  0.00000000e 00 -1.38777878e-17
#  -2.77555756e-17 -2.77555756e-17  0.00000000e 00  0.00000000e 00
#   5.55111512e-17 -1.11022302e-16  1.11022302e-16  5.55111512e-17
#  -1.11022302e-16  0.00000000e 00 -3.88578059e-16 -1.11022302e-16
#  -2.77555756e-16  1.66533454e-16  5.55111512e-17  2.22044605e-16
#   5.55111512e-17 -2.22044605e-16]

Alternatively, I could use fsolve from Scipy optimize, which I have already verified it computes the correct solutions. However, this procedure will be executed many times by another iterative algorithm, so I'd rather use the analytical solutions.

CodePudding user response：

Seems like your analytical solution has a mistake.

I rewrote the target function in form

Cp(ϕ) = 2 * (A   B sin ϕ - C cos ϕ)^2

Solutions of Cp = 0 come from A B sin ϕ - C cos ϕ = 0 equation, which I solved using

I played with value of β, seems like the code above works good for both zero and non-zero degrees. Also, I didn't care about function domains, you should reproduce the derivation to ensure that all parameters and solution are well-defined.