Здравствуйте (это мой первый пост в переполнении стека), я пытаюсь вычислить первые 3 энергетических уровня частицы в гармоническом потенциале, используя метод стрелок
Код адаптирован из скриптав «Вычислительной физике» Марка Ньюмана этот скрипт вычислял основное состояние для частицы в коробке.
вот ссылка http://www -personal.umich.edu / ~ mejn / cp / Programs /squarewell.py
вот мой адаптированный код
import numpy as np
from scipy.constants import m_e,hbar,elementary_charge
#define constants
vo = 50
a = 1e-11
#define limits of integration for adaptive runge-kutta method
xi = -10*a
xf = 10*a
N = 1000
dx = (xf-xi)/N
def V(x):#define harmonic potential well
return vo*((x**2)/(a**2))
def f(r,x,E): #schrodinger's equation
psi = r[0]
phi = r[1]
fpsi = psi
fphi = (2*m_e/(hbar**2))*(V(x)-E)*psi
return np.array([fpsi,fphi],float)
def solve(E): #calculates wave function for an energy E
psi = 0.0
phi = 1.0
r = np.array([psi,phi],float)
for x in np.arange(xi,xf,dx): #adaptive runge-kutta method
k1 = dx*f(r,x,E)
k2 = dx*f(r+0.5*k1,x+0.5*dx,E)
k3 = dx*f(r+0.5*k2,x+0.5*dx,E)
k4 = dx*f(r+k3,x+dx,E)
r += (k1+2*k2+2*k3+k4)/6
return r[0]
#finds the energy using secant method
E1 = 0.0
E2 = elementary_charge
psi2 = solve(E1)
target = elementary_charge/1000
while abs(E1-E2)>target:
psi1,psi2 = psi2,solve(E2)
E1,E2 = E2,E2-psi2*(E2-E1)/(psi2-psi1)
print (E2/elementary_charge)
при запуске я получаю эту ошибку
RuntimeWarning: invalid value encountered in double_scalars
E1,E2 = E2,E2-psi2*(E2-E1)/(psi2-psi1)
, что, я думаю, означает, что psi2 и psi1слишком близко друг к другу, но я не совсем уверен, как это исправить