Korteweg-de Vries 方程

日期2018-10-24（最后修改），2010-09-06（创建）

此页面展示了如何使用 Korteweg-de Vries 方程 在周期域上求解，使用 线方法，空间导数使用伪谱方法计算。在此方法中，导数在频域中计算，首先对数据应用 FFT，然后乘以适当的值，并使用逆 FFT 转换回空间域。这种微分方法由 scipy.fftpack 模块中的 diff 函数实现。

我们离散化空间域，并使用 scipy.fftpack 模块中定义的 diff 函数计算空间导数。在以下代码中，此函数被赋予别名 psdiff，以避免将其与 numpy 函数 diff 混淆。通过仅离散化空间维度，我们得到一个常微分方程组，该方程组在 kdv(u, t, L) 函数中实现。kdv_solution(u0, t, L) 函数使用 scipy.integrate.odeint 来求解此系统。

在 [1]

#!python

import numpy as np
from scipy.integrate import odeint
from scipy.fftpack import diff as psdiff


def kdv_exact(x, c):
    """Profile of the exact solution to the KdV for a single soliton on the real line."""
    u = 0.5*c*np.cosh(0.5*np.sqrt(c)*x)**(-2)
    return u

def kdv(u, t, L):
    """Differential equations for the KdV equation, discretized in x."""
    # Compute the x derivatives using the pseudo-spectral method.
    ux = psdiff(u, period=L)
    uxxx = psdiff(u, period=L, order=3)

    # Compute du/dt.    
    dudt = -6*u*ux - uxxx

    return dudt

def kdv_solution(u0, t, L):
    """Use odeint to solve the KdV equation on a periodic domain.
    
    `u0` is initial condition, `t` is the array of time values at which
    the solution is to be computed, and `L` is the length of the periodic
    domain."""

    sol = odeint(kdv, u0, t, args=(L,), mxstep=5000)
    return sol


if __name__ == "__main__":
    # Set the size of the domain, and create the discretized grid.
    L = 50.0
    N = 64
    dx = L / (N - 1.0)
    x = np.linspace(0, (1-1.0/N)*L, N)

    # Set the initial conditions.
    # Not exact for two solitons on a periodic domain, but close enough...
    u0 = kdv_exact(x-0.33*L, 0.75) + kdv_exact(x-0.65*L, 0.4)

    # Set the time sample grid.
    T = 200
    t = np.linspace(0, T, 501)

    print("Computing the solution.")
    sol = kdv_solution(u0, t, L)


    print("Plotting.")

    import matplotlib.pyplot as plt

    plt.figure(figsize=(6,5))
    plt.imshow(sol[::-1, :], extent=[0,L,0,T])
    plt.colorbar()
    plt.xlabel('x')
    plt.ylabel('t')
    plt.axis('normal')
    plt.title('Korteweg-de Vries on a Periodic Domain')
    plt.show()

Computing the solution.
Plotting.

章节作者：WarrenWeckesser，Warren Weckesser

附件

• kdv.png
• kdv2.png