python 미분 방정식 해결 작업(수치 해법)

Python 미분 방정식 풀이(수치 해법)
일부 미분 방정식 에 있어 서 수치 해법 은 원래 의 방정식 을 구하 기 어렵 기 때문에 풀이 에 좋 은 도움 이 된다.
예 를 들 어 방정식:

그러나 우 리 는 그것 의 초기 조건 을 알 게 되 었 다.이것 은 우리 가 세대 간 의 구 해 에 매우 도움 이 되 고 필수 적 인 것 이다.

그러면 지금 우 리 는 Python 으로 이 문 제 를 해결 합 니 다.일반적인 수치 해법 은 오로라 법,암시 적 사다리꼴 법 등 이 있 습 니 다.우 리 는 이런 알고리즘 들 이 중첩 의 정밀도 에 어떤 차이 가 있 는 지 살 펴 보 겠 습 니 다.

```python
```python
import numpy as np
from scipy.integrate import odeint
from matplotlib import pyplot as plt
import os
#  odeint          
#       
class Euler:
    #    ，        ，       
    def __init__(self,h,y0):
        #             
        self.h = h
        self.y0 = y0
        self.y = y0
        self.n = 1/self.h
        self.x = 0
        self.list = [1]
        #    list  , x y    
        self.list2 = [1]
        self.y1 = y0
        #      list2  ， x y1    
        self.list3 = [1]
        self.y2 = y0
        #      list3  ， x y2    
    #      ，    t,x
    def countall(self):
        for i in range(int(self.n)):
            y_dere = -20*self.list[i]
            #      y_dere = -20 * x
            y_dere2 = -20*self.list2[i] + 0.5*400*self.h*self.list2[i]
            #         y_dere2 = -20*x(k) + 0.5*400*delta_t*x(k)
            y_dere3 = (1-10*self.h)*self.list3[i]/(1+10*self.h)
            #        y_dere3 = (1-10*delta_t)*x(k)/(1+10*delta_t)
            self.y += self.h*y_dere
            self.y1 += self.h*y_dere2
            self.y2 =y_dere3
            self.list.append(float("%.10f" %self.y))
            self.list2.append(float("%.10f"%self.y1))
            self.list3.append(float("%.10f"%self.y2))
        return np.linspace(0,1,int(self.n+1)), self.list,self.list2,self.list3
step = input("           :")
step = float(step)
work1 = Euler(step,1)
ax1,ay1,ay2,ay3 = work1.countall()
#    plt
plt.figure(1)
plt.subplot(1,3,1)
plt.plot(ax1,ay1,'s-.',MarkerFaceColor = 'g')
plt.xlabel('   t',fontproperties = 'simHei',fontsize =20)
plt.ylabel('   x',fontproperties = 'simHei',fontsize =20)
plt.title('              '+str(step),fontproperties = 'simHei',fontsize =20)
plt.subplot(1,3,2)
plt.plot(ax1,ay2,'s-.',MarkerFaceColor = 'r')
plt.xlabel('   t',fontproperties = 'simHei',fontsize =20)
plt.ylabel('   x',fontproperties = 'simHei',fontsize =20)
plt.title('                '+str(step),fontproperties = 'simHei',fontsize =20)
plt.subplot(1,3,3)
plt.plot(ax1,ay3,'s-.',MarkerFaceColor = 'b')
plt.xlabel('   t',fontproperties = 'simHei',fontsize =20)
plt.ylabel('   x',fontproperties = 'simHei',fontsize =20)
plt.title('                '+str(step),fontproperties = 'simHei',fontsize =20)
plt.figure(2)
plt.plot(ax1,ay1,ax1,ay2,ax1,ay3,'s-.',MarkerSize = 3)
plt.xlabel('   t',fontproperties = 'simHei',fontsize =20)
plt.ylabel('   x',fontproperties = 'simHei',fontsize =20)
plt.title('        '+str(step),fontproperties = 'simHei',fontsize =20)
ax = plt.gca()
ax.legend(('$Eular$','$fixed Eular$','$trapezoid$'),loc = 'lower right',title = 'legend')
plt.show()
os.system("pause")

오로라 법 에 대한 중첩 방법 은:

오로라 법의 중첩 방법 개선:

암시 적 사다리꼴 법:

서로 다른 보폭 에 대해 그 구 해 의 정밀도 도 크게 다 를 것 이다.나 는 먼저 몇 장의 결과 도 를 놓 았 다.

보충:python 기반 미분 방정식 수치 해법 회로 모델
환경 패키지 설치
numpy 설치(range 조절 에 사용)와 matplotlib 설치(그림 그리 기 에 사용)
명령 줄 에 입력

pip install numpy 
pip install matplotlib

회로 모형 과 미분 방정식
모델 1
무 손해,커 패 시 터 전압 5V,커 패 시 터 0.01F,인덕턴스 0.01H 의 병렬 공진 회로
회로 모형 1

미분 방정식
모형 2
저항 손실 이 있 는 커 패 시 터 전압 은 5V 이 고 커 패 시 터 는 0.01F 이 며 인덕턴스 는 0.01H 의 병렬 공진 이다.
회로 모형 2

미분 방정식
python 코드
모델 1

import numpy as np
import matplotlib.pyplot as plt
 
L = 0.01  #     F
C = 0.01  #     L
u_0 = 5   #       
u_dot_0 = 0
 
def equition(u,u_dot):#    
    u_double_dot = -u/(L*C)
    return u_double_dot
 
def draw_plot(time_step,time_scale):#       
    u = u_0
    u_dot = u_dot_0  #            
    time_list = [0] #  lis
    Votage = [u] #  list
    plt.figure()
    for time in np.arange(0,time_scale,time_step):#              
        u_double_dot = equition(u,u_dot) #    
        u_dot = u_dot + u_double_dot*time_step #    
        u = u + u_dot*time_step #  
        time_list.append(time) #    
        Votage.append(u) #    
        print(u)
    plt.plot(time_list,Votage,"b--",linewidth=1) #  
    plt.show()
    plt.savefig("easyplot.png")
 
if __name__ == '__main__':
    draw_plot(0.0001,1)

모형 2

import numpy as np
import matplotlib.pyplot as plt
 
L = 0.01  #     F
C = 0.01  #     L
R = 0.1   #   
u_0 = 5   #       
u_dot_0 = 0
 
def equition(u,u_dot):#    
    u_double_dot =(-R*C*u_dot -u)/(L*C)
    return u_double_dot
 
def draw_plot(time_step,time_scale):#       
    u = u_0
    u_dot = u_dot_0  #            
    time_list = [0] #  lis
    Votage = [u] #  list
    plt.figure()
    for time in np.arange(0,time_scale,time_step):#              
        u_double_dot = equition(u,u_dot) #    
        u_dot = u_dot + u_double_dot*time_step #    
        u = u + u_dot*time_step #  
        time_list.append(time) #    
        Votage.append(u) #    
        print(u)
    plt.plot(time_list,Votage,"b-",linewidth=1) #  
    plt.show()
    plt.savefig("result.png")
 
if __name__ == '__main__':
    draw_plot(0.0001,1)

수치 분해 결과
모델 1

세로 축 은 커 패 시 터 양 끝 전압 이 고 가로 축 은 시간 과 공식 계산 이 일치 합 니 다.
모델 2 결과

세로 축
커 패 시 터 양 끝 전압,가로 축 은 시간 제목 입 니 다.
마지막 으로 우 리 는 조절 저항 에 따라 서로 다른 상태 에 도달 할 수 있다.

R=0.01,저항 부족

R=1.7,임계 저항

R=100,과 저항
이상 은 개인 적 인 경험 이 므 로 여러분 에 게 참고 가 되 기 를 바 랍 니 다.여러분 들 도 저 희 를 많이 응원 해 주시 기 바 랍 니 다.