सिम्पी दूसरा ऑर्डर ode

Question

मेरे पास एक साधारण दूसरे ऑर्डर ओडीई के लिए एक सजातीय समाधान है, जो जब मैं सिम्पी का उपयोग करके प्रारंभिक मानों के लिए हल करने का प्रयास करता हूं, तो वही समाधान देता है। इसे y (0) और y '(0) के लिए प्रतिस्थापित करना चाहिए और स्थिरांक के बिना समाधान उत्पन्न करना चाहिए, लेकिन नहीं। समीकरण स्थापित करने के लिए कोड यहां दिया गया है (यह एक वसंत संतुलन समीकरण है, के = वसंत स्थिरांक और एम = द्रव्यमान)। कुछ अनावश्यक प्रतीक जो मैं कहीं और उपयोग करता हूं, क्षमा करें।

%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 

परिणाम (सही) है: $ y {\ छोड़ दिया (टी \ right)} = C_ {1} ई^{- टी \ sqrt {- \ frac {कश्मीर} {मीटर}}} + सी_ {2} ई^{टी \ sqrt {- \ frac {k} {m}}} $

वाई (टी) = सी 1e^(- टी√ (-के/एम)) + सी 2e^(t√ (-km)), जिसमें y_n = c1.cos (√ (-k/m) t) + c2.sin (√ (-k/m) t) है।

जब यह समीकरण विश्लेषणात्मक रूप से हल किया जाता है, और omega = sqrt (-k/m) के साथ साइन और कोसाइन का उपयोग करके समाधान में परिवर्तित होता है, तो c1 = y (0) और c2 = y '(0)/omega। इसलिए, जबकि समाधान आंशिक रूप से जटिल संख्याओं को शामिल कर रहा है, निश्चित रूप से, डिस्लोव ऊपर के रूप में मूल समरूप समीकरण देता है। y पर स्तोत्र का मूल्यांकन करने के मेरे कोड (0) और वाई '(0) है:

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

मुझे खुशी है कि dsolve बस विश्लेषणात्मक इस आईवीपी को संभालने में सक्षम नहीं हो सकता है, लेकिन अगर यह इतना है कि मैं आश्चर्य होगा इसकी दूसरी क्षमता के आधार पर। इस समस्या के बारे में कोई मदद और इसलिए अन्य विश्लेषणात्मक दूसरी ऑर्डर समस्याओं को हल किया जा सकता है, इसकी सराहना की जाएगी। सवाल का थक्का के आसपास है:

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

तो समाधान मैंने कोशिश की है, जो इस बात की पुष्टि Dietrich, था:

#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)]) 

हालांकि यह है "एक" समाधान है, इसका एक बहुत ही जटिल तरीका लगता है y (टी) = y (0) cos (omega * t - phi) ढूंढ रहा है ... जो इस सॉल्वर की कुछ सीमाओं के बारे में निहित प्रश्न का उत्तर देता है और इस बारे में सीधा सवाल है कि कैसे ics तर्क हल किया गया है। धन्यवाद।

Answer 1

पैरामीटर dsolve() में ics वास्तव में काम नहीं करता है (Issue 4720) है, तो आप मैन्युअल रूप से प्रतिस्थापन क्या करना है। आप कोशिश कर सकते हैं:

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))