Matplotlib 不会为函数设置动画

Matplotlib won't animate a function

我正在尝试制作以下波函数的动画:

出于某种原因,我的代码适用于 n=0,但不适用于任何其他 n。我检查了各种 n、x 和 t 的函数值,它似乎工作正常但由于某种原因 matplotlib 没有动画。这是代码:

%matplotlib qt

import numpy as np
from sympy import hermite
from scipy.special import gamma
import matplotlib.pyplot as plt
from matplotlib import animation


def psi(n, x, t):
    A = (np.pi ** -0.25) * (1 / np.sqrt((2 ** n) * gamma(n + 1)))
    E = n + 0.5
    f = A * hermite(n, x) * np.exp(-(x ** 2) / 2) * np.cos(E * t)
    return f


def animar(f, x0=0, xf=1, dx=0.01, t0=0, tf=1, dt=0.01, ym=-2, yM=2):
    nf = int((xf - x0) / dx + 1)
    nt = int((tf - t0) / dt + 1)
    x = np.linspace(x0, xf, nf)
    t = np.linspace(t0, tf, nt)

    fig, ax = plt.subplots()

    ax.set_xlim((x0, xf))
    ax.set_ylim((ym, yM))

    line, = ax.plot([], [], lw=2)

    def init():
        line.set_data([], [])
        return line,

    def animate(i):
        y = f(x, i)
        line.set_data(x, y)
        return line,

    anim = animation.FuncAnimation(fig, animate, init_func=init,
                                   frames=5 * t, interval=20, blit=True)
    plt.show()
    return anim


F = lambda x,t: psi(1, x, t)
anim = animar(F, x0=-3, xf=3, tf=3, ym=-1, yM=1)

这是一个修复程序。将 psi 函数重新定义为:

def psi(n, x, t):
    A = (np.pi ** -0.25) * (1 / np.sqrt((2 ** n) * gamma(n + 1)))
    E = n + 0.5
    herms = np.array([hermite(n, xelem) for xelem in x])
    f = A * herms * np.exp(-(x ** 2) / 2) * np.cos(E * t)
    return f

基本上hermite好像有点weird/buggy。如果 x 是一个 numpy 数组,那么 hermite(n=0, x) 就可以正常工作。但是,如果 n>0 那么 hermite 会在 x 不是单个值时痛苦地抱怨。因此,我们通过一次输入一个从 xhermite 的值来解决这个问题,然后根据这些结果制作 herms 数组。

量子力学是很多年前的事了,所以我不能告诉你 hermite 相对于数组输入的不一致是否有正当理由。我的猜测是,这只是 sympy 中的一个稍微有问题的实现。无论如何,一旦你修复了 psi 中的这两行,动画似乎就可以正常工作了。

我不认为你想象征性地计算埃尔米特多项式。相反,为什么不使用 numpy.polynomial.hermite.hermval? Also maybe use scipy.special.factorial 来计算阶乘。

import numpy as np
from numpy.polynomial.hermite import hermval
from scipy.special import factorial
import matplotlib.pyplot as plt
from matplotlib import animation


def psi(n, x, t):
    A = (np.pi ** -0.25) * (1 / np.sqrt((2 ** n) * factorial(n)))
    E = n + 0.5
    nar = np.zeros(n+1)
    nar[-1] = 1
    f = A * hermval(x,nar) * np.exp(-(x ** 2) / 2) * np.cos(E * t)
    return f


def animar(f, x0=0, xf=1, dx=0.01, t0=0, tf=1, dt=0.01, ym=-2, yM=2):
    nf = int((xf - x0) / dx + 1)
    nt = int((tf - t0) / dt + 1)
    x = np.linspace(x0, xf, nf)
    t = np.linspace(t0, tf, nt)

    fig, ax = plt.subplots()

    ax.set_xlim((x0, xf))
    ax.set_ylim((ym, yM))

    line, = ax.plot([], [], lw=2)

    def init():
        line.set_data([], [])
        return line,

    def animate(i):
        y = f(x, i)
        line.set_data(x, y)
        return line,

    anim = animation.FuncAnimation(fig, animate, init_func=init,
                                   frames=5 * t, interval=20, blit=True)
    plt.show()
    return anim


F = lambda x,t: psi(3, x, t)
anim = animar(F, x0=-3, xf=3, tf=3, ym=-1, yM=1)