Решение для нулей задачи Штурма Лиувилля в Симпи

Question

Я пытаюсь найти нули задачи Штурма Лиувилля в sympy

, вот мой код

from sympy import *                                                            
init_printing(use_latex='mathjax')                                             

# Function Definitions                                                         
u = Function('u')                                                              
X = Function('X')                                                              

# Variables and Dummies                                                        
x  = Dummy('x')                                                                                                                                               
mu  = Symbol('mu', nonzero=True, positive=True)                                
C1 = Symbol('C1', nonzero=True)                                                
C2 = Symbol('C2', nonzero=True)                                                

print("STURM LIOUVILLE PROBLEM IN X")                                        
print("GENERAL SOLUTION")                                                      
# --- --- ---#                                                                 
print("")                                                                      
gen_sol = dsolve(Eq(Derivative(X(x), x, 2) + (mu**2)*X(x), 0), X(x))           
print(pretty(gen_sol))                                                         
print("")                                                                      
print("SOLVE FOR BOUNDARY")                                                    
bc_2 = gen_sol.rhs.diff(x).subs(x,5).subs(Symbol('C1')*mu, Symbol('C2'))       
print("DETERMINE THAT C2 != 0 (NON TRIVIAL)")                                  
print(pretty(bc_2))                                                            
print("")                                                                      
print("SOLVE FOR REMAINING SOLUTION")                                          
bc_2 = bc_2.subs(Symbol('C2'), 1)                                              
print(pretty(bc_2))                                                            

print("")                                                                      
linear_solutions = solve([bc_2], [mu])                                         
print(pretty(linear_solutions))                               

К сожалению, я получаю следующую NotImplemented Error.

Есть ли другой способ решения нулей в sympy, который может справиться с этой операцией?

STURM LIOUVILLE PROBLEM IN X
GENERAL SOLUTION

X(x) = C₁⋅sin(x⋅μ) + C₂⋅cos(x⋅μ)

SOLVE FOR BOUNDARY
DETERMINE THAT C2 != 0 (NON TRIVIAL)
-C₂⋅μ⋅sin(5⋅μ) + C₂⋅cos(5⋅μ)

SOLVE FOR REMAINING SOLUTION
-μ⋅sin(5⋅μ) + cos(5⋅μ)

Traceback (most recent call last):
  File "stack_overflow.py", line 31, in <module>
    linear_solutions = solve([bc_2], [mu])
  File "/home/cgould/.local/lib/python3.8/site-packages/sympy/solvers/solvers.py", line 1176, in solve
    solution = _solve_system(f, symbols, **flags)
  File "/home/cgould/.local/lib/python3.8/site-packages/sympy/solvers/solvers.py", line 1943, in _solve_system
    raise NotImplementedError('could not solve %s' % eq2)

(РЕДАКТИРОВАТЬ ДЛЯ ЯСНОСТИ) Идея в том, что я пытаюсь найти нули из ниже

-μ⋅sin(5⋅μ) + cos(5⋅μ) = 0

, что примерно означает:

cot(5μ) = μ 

Answer 1

solve будет решать только уравнение, для которого было запрограммировано решение в замкнутой форме, что дает точное решение для символов c. nsolve найдет числовое приближение:

>>> μ = mu
>>> nsolve(-μ*sin(5*μ) + cos(5*μ), 0)
0.262767543298580

Поскольку ваше уравнение имеет несколько корней, вам придется использовать другое начальное предположение, чтобы получить другое решение (или использовать метод деления пополам и регион в который вы хотите найти root). Хороший способ сделать это - найти первые два корня, а затем использовать разницу между ними в качестве предположения о расстоянии до следующего root: метод продолжения первого порядка.

>>> r = []
>>> r.append(nsolve(-μ*sin(5*μ) + cos(5*μ),0))
>>> r.append(nsolve(-μ*sin(5*μ) + cos(5*μ),.5))
>>> for i in range(5):
...   r.append(nsolve(-μ*sin(5*μ) + cos(5*μ),2*r[-1] - r[-2]))
...
>>> [i.round(2) for i in r]
[0.26, 0.81, 1.38, 1.98, 2.59, 3.20, 3.82]