波函数的 Crank Nicolson 方法 Python

Question

我正在尝试在虚时间（将时间步长乘以单位虚数）中使用曲柄尼科尔森方法传播高斯波包。我为实现此目的而编写的代码如下所示：

import matplotlib.pyplot as plt #this allows you to plot, and changes the name to plt
import numpy as np #this allows you to do math, and changes the name to np
import math
import scipy.linalg as la

def V(x):
#    k = 1
#    v = k*x**4
    v = 0.25*(x-3)**2+0.15*(x-3)**4
    return v

def Psi(x):
    psi = np.exp(-2*(x-3)**2)
    return psi

#Function for computing integral using trapezoid method 
def TrapInt(y, h):
    trap = [(float(y[ii]) + float(y[ii+1])) for ii in range(0, len(y)-1)]
    return float(h)/2*sum(trap)    

N = 1000
L = 3;
h = 0.01
x = np.arange(0,6,h);
t = np.linspace(0,L,300);
t = 1j*t;
dt = t[1] - t[0]
dx = x[1] - x[0]

A = 1j*dt/(2*dx**2)

pot = V(x)

Q = np.zeros([len(x),len(x)],dtype = complex)
P = np.zeros([len(x),len(x)],dtype = complex)

wave = np.zeros([len(x),len(t)],dtype = complex)
wave[:,0] = Psi(x)

B = (1- 2*A - 1j*dt*pot)

for ii in range(0,len(x)-1):
    Q[ii][ii] = -(B[ii])
    P[ii][ii] = (B[ii])
    Q[ii][ii+1] = (2-A)
    P[ii][ii+1] = A
    if ii >= 1:
        Q[ii][ii-1] = -A
        P[ii][ii-1] = A
plt.plot(wave[:,0])
for ii in range(0,len(t)-1):
    one = np.matmul(P,wave[:,ii])
    wave[:,ii+1] = np.matmul(la.inv(Q),one)

我在曲柄尼科尔森方法的实现中似乎找不到任何数学错误；但是，每当我尝试运行时，它都会给出一个错误，指出 Q 是单数（没有逆）。我不确定为什么会这样。任何帮助表示赞赏。谢谢

Answer 1

您从未分配给 Q[-1]。已知零行在某些情况下会产生奇异矩阵。

另外，矩阵不要反复反转。可能根本不反转它，而是存储它的一些分解以允许有效计算 Q^-1x .

波函数的 Crank Nicolson 方法 Python

Crank Nicolson Method on Wave Function Python

python

matrix

numerical-methods