我如何使用matplotlib为logistic方程K = 1生成矢量场图？

Question

我正在研究Strogatz的《非线性动力学和混沌》，我在第二章练习2.8.1中遇到了麻烦。 （教育工作者的旗帜：我毕业了，所以这不是一堂课，我只是想重新回到微分方程的数值求解中），这是一个非常简单的微分方程，在不同的初始条件下，我可以绘制各个解曲线但我正在尝试使用颤动或流图将单个解决方案叠加在向量字段之上。

我的问题在于了解如何将发现的dy / dx形式的矢量场图转换为dy / dx形式，并将其转换为Strogatz的书主要解决的dx / dt形式。

考虑到逻辑函数中定义的x向量仅是一维的，我很难解说出u和v在颤动或流图中的表达方式，因为问题似乎只有u流动。它可能超级简单，并且正在考虑中，但是任何指导或帮助将不胜感激！

到目前为止，我有以下内容：

here

# 2.8.1
# Plot the vector field and some trajectories for xdot = x(1-x) given 
# some different initial conditions for the logistic equation with carrying 
# capacity K = 1

# dx/dt = x(1-x)  

# Imports:
from __future__ import division
from scipy import *
import numpy as np
import pylab
import matplotlib as mp
from matplotlib import pyplot as plt  
import sys
import math as mt

def logistic(x,t):
    return np.array([x[0]*(1-x[0])])

def RK4(t0 = 0, x0 = np.array([1]), t1 = 5 , dt = 0.01, ng = None):  
    tsp = np.arange(t0, t1, dt)
    Nsize = np.size(tsp)
    X = np.empty((Nsize, np.size(x0)))
    X[0] = x0

    for i in range(1, Nsize):
        k1 = ng(X[i-1],tsp[i-1])
        k2 = ng(X[i-1] + dt/2*k1, tsp[i-1] + dt/2)
        k3 = ng(X[i-1] + dt/2*k2, tsp[i-1] + dt/2)
        k4 = ng(X[i-1] + dt*k3, tsp[i-1] + dt)
        X[i] = X[i-1] + dt/6*(k1 + 2*k2 + 2*k3 + k4)
    return X

def tplot():
    t0 = 0
    t1 = 10
    dt = 0.02
    tsp = np.arange(t0,t1,dt)
    X = RK4(x0 = np.array([2]), t1 = 10,dt = 0.02, ng = logistic)
    Y = RK4(x0 = np.array([0.01]), t1 = 10,dt = 0.02, ng = logistic)
    Z = RK4(x0 = np.array([0.5]), t1 = 10,dt = 0.02, ng = logistic)
    P = RK4(x0 = np.array([3]), t1 = 10,dt = 0.02, ng = logistic)
    Q = RK4(x0 = np.array([0.1]), t1 = 10,dt = 0.02, ng = logistic)
    R = RK4(x0 = np.array([1.5]), t1 = 10,dt = 0.02, ng = logistic)
    O = RK4(x0 = np.array([1]), t1 = 10,dt = 0.02, ng = logistic)
    pylab.figure()
    pylab.plot(tsp,X)
    pylab.plot(tsp,Y)
    pylab.plot(tsp,Z)
    pylab.plot(tsp,P)
    pylab.plot(tsp,Q)
    pylab.plot(tsp,R)
    pylab.plot(tsp,O)
    pylab.title('Logistic Equation - K=1')
    pylab.xlabel('Time')
    pylab.ylabel('Xdot')
    pylab.show()

print tplot()

Answer 1

要从导数绘制斜率（例如dx / dt），您可以先找到dx / dt，然后使用固定的[[dt计算dx 。然后，在每个感兴趣的（t，x）处，绘制从（t，x）到（t + dt，x + dx）的小线段。这里是您的方程式的一个例子

dx / dt = x（1-x）

。 （Strogatz图片中没有箭头，因此我也将其删除了。）image here

Answer 2

wonkybadonk：由于绘制轨迹的斜率与绘制矢量场的斜率不同，似乎是由于您的矢量场不够陡峭。确保dx = dxdt * dt; （逐点相乘，而不是点积）并且您添加了“ angles ='xy'”作为颤振参数。 （请参阅tom10帖子）。