Thứ tự thứ hai Sympy ode

Question

Tôi có một giải pháp đồng nhất cho một ODE thứ hai đơn giản, khi tôi cố gắng giải quyết các giá trị ban đầu bằng Sympy, trả về cùng một giải pháp. Nó nên thay thế cho y (0) và y '(0) và mang lại một giải pháp không có hằng số, nhưng không. Đây là mã để thiết lập phương trình (nó là một phương trình cân bằng mùa xuân, k = hằng số mùa xuân và m = khối lượng). Một số biểu tượng thừa mà tôi sử dụng ở nơi khác, xin lỗi.

%matplotlib inline 
from sympy import * 
m,k, x,t,s,T, omega,A,B = symbols('m k x t s T omega A B',float=True) 
a = symbols('a', positive=True) 
f,y,g,H,delta,S=symbols('f y g H delta S',cls=Function) 
Eq1 = Eq(m*diff(y(t),t,2)+k*y(t)) 
Eq1 

Kết quả là (chính xác): $ y {\ left (t \ right)} = C_ {1} e^{- t \ sqrt {- \ frac {k} {m}}} + C_ {2} e^{t \ sqrt {- \ frac {k} {m}}} $

y (t) = C1e^(- t√ (−k/m)) + C2e^(t√ (−km)), cũng có y_n = c1.cos (√ (−k/m) t) + c2.sin (√ (−k/m) t).

Khi phương trình này được giải quyết một cách phân tích và chuyển đổi thành giải pháp sử dụng sin và cosin với omega = sqrt (-k/m), rồi c1 = y (0) và c2 = y '(0)/omega. Vì vậy, trong khi giải pháp là một phần liên quan đến số phức, tất nhiên, dsolve chỉ đơn giản là trả về phương trình thuần nhất ban đầu như trên. Mã của tôi để đánh giá ODE tại y (0) và y '(0) là:

Eq1_soln_IVP =dsolve(Eq1,y(t),x0=0, ics={y(0): a, y(t).diff(t).subs(t, 0): a}) 

Tôi đánh giá cao dsolve mà chỉ đơn giản có thể không có khả năng xử lý IVP này phân tích, nhưng tôi sẽ ngạc nhiên nếu điều này là rất dựa trên khả năng khác của nó. Bất kỳ trợ giúp như thế nào vấn đề này và do đó các vấn đề thứ tự phân tích thứ hai khác có thể được giải quyết sẽ được nhiều đánh giá cao. Nub của câu hỏi là xung quanh:

ics={y(0): a, y(t).diff(t).subs(t, 0): a} 

Vì vậy, các giải pháp tôi đã cố gắng, mà Dietrich khẳng định, là:

#Create IVP for y(0) 
expr = Eq(Eq1_soln_IVP.rhs.subs(sqrt(-k/m),I*omega),y(0)) 
#Create IVP for y'(0) 
expr2 = Eq(diff(y(t),t).subs(t,0),expr.lhs.diff(t)) 
#Maps all free variables and solves for each where t = 0. 
solve([expr.subs(t,0),expr2.subs(t,0)]) 

Trong khi đó là "một" giải pháp, điều này dường như là một cách rất phức tạp của tìm y (t) = y (0) cos (omega * t - phi) ... trả lời câu hỏi ngầm về một số hạn chế của người giải quyết này và câu hỏi trực tiếp về cách mà ics arg được giải quyết. Cảm ơn.

Answer 1

Tham số ics trong dsolve() không thực sự hoạt động (Issue 4720), vì vậy bạn phải thực hiện thay thế theo cách thủ công. Bạn có thể thử:

from IPython.display import display 
import sympy as sy 

sy.init_printing() # LaTeX-like pretty printing for IPython 

t = sy.Symbol("t", real=True) 
m, k = sy.symbols('m k', real=True) # gives C_1 Exp() + C_2 Exp() solution 
# m, k = sy.symbols('m k', positive=True) # gives C_1 sin() + C_2 cos() sol. 
a0, b0 = sy.symbols('a0, b0', real=True) 
y = sy.Function('y') 

Eq1 = sy.Eq(m*sy.diff(y(t), t, 2) + k*y(t)) 
print("ODE:") 
display(Eq1) 

print("Generic solution:") 
y_sl0 = sy.dsolve(Eq1, y(t)).rhs # take only right hand side 
display(sy.Eq(y(t), y_sl0)) 

# Initial conditions: 
cnd0 = sy.Eq(y_sl0.subs(t, 0), a0) # y(0) = a0 
cnd1 = sy.Eq(y_sl0.diff(t).subs(t, 0), b0) # y'(0) = b0 

# Solve for C1, C2: 
C1, C2 = sy.symbols("C1, C2") # generic constants 
C1C2_sl = sy.solve([cnd0, cnd1], (C1, C2)) 

# Substitute back into solution: 
y_sl1 = sy.simplify(y_sl0.subs(C1C2_sl)) 
print("Solution with initial conditions:") 
display(sy.Eq(y(t), y_sl1))