避免绘制 ODE 的不同解 ODEint

Question

我正在尝试绘制相平面，我希望它看起来不错。然而，方程组的一些解由于初始条件而发散。有什么方法可以在解决方案发散时按顺序制作 try/except 链，但不会绘制它。这是我的代码：

import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate import odeint
import pylab as pl 

def aux_func(x):  
  y = x[0]-x[1]
  if (np.abs(y) <= 1):
    f = y**3 + 0.5*y 
  else: 
    f = 2*y - np.sign(y)
  return f 

def function(x,t): 
  x1_dot = x[1]  
  x2_dot = -x[1] - aux_func(x)  
  return [x1_dot,x2_dot]

ts = np.linspace(0, 20, 300)  
ic_1 = np.linspace(-1,1,10)
ic_2 = np.linspace(-1,1,10) 

for r1 in ic_1:
  for r2 in ic_2:
    x0 = (r1,r2)
    try:
      xs = odeint(function, x0, ts)
      plt.plot(xs[:,0], xs[:,1],"r-",linewidth=.8) 
    except: 
      pass

# Nombre de los ejes, limites,
plt.xlabel("$x_1$", fontsize=12)
plt.ylabel("$x_2$", fontsize=12)
# plt.tick_params(labelsize=10)
# plt.xticks(np.linspace(0,1,11))
# plt.yticks(np.linspace(0,1,11))
plt.xlim(-1, 1)
plt.ylim(-1, 1)
# Grafica el campo vectorial
X1, X2 = np.mgrid[-1:1:20j,-1:1:20j]
u=X2
d= X1-X2

t = np.zeros(np.shape(d))
for i in range(len(d)):
    for j in range(len(d[0])):
        if np.abs(d[i][j]) > 1:
            t[i][j]= 2*d[i][j]-0.5*np.sign(d[i][j])
        else:
            t[i][j] =d[i][j]**3 + 0.5*d[i][j]
v=-X2-t

pl.quiver(X1, X2, u, v, color = 'b',width = .002)
plt.grid()
plt.title('Plano de Fase Punto 1')
#plt.savefig('FasePunto4.png')
plt.show()

代码绘制如下：

感谢帮助。

Answer 1

完全避免错误的发散就可以解决，不需要异常处理

这是一个不连续的 ODE，可能会导致不寻常的效果，例如滑动模式。快速解决这个问题的一种方法是通过实施一个混合区来缓和跳跃，在该混合区中矢量场快速但连续地从一个阶段变化到另一个阶段（参见了解其他通用 work-arounds）。对此的更改可以实现为

def aux_func(x):  
  def softsign(u): return np.tanh(1e4*u)
  y = x[0]-x[1]
  h = 0.5*(1+softsign(y**2-1)
  # h is about zero for |y|<1 and about 1 for |y|>1
  f1 = y**3 + 0.5*y          # for |y|<1
  f2 = 2*y - softsign(y)     # for |y|>1, note the second mollification
  return (1-h)*f1+h*f2

无需进一步更改代码即可得出图

请注意 pylab 已过时，它的所有功能也可以通过 plt=matplotlib.pyplot 访问。

避免绘制 ODE 的不同解 ODEint

Avoiding plotting ODEs divergent solutions ODEint

python

numpy

visualization

ode