I've been trying to use either solve_ivp or solve_bvp to solve a problem I've been having but I'm not making any progress. I think the code I have here will work but I cannot get the range to be correct. for some reason I cannot understand the range is always going from 0 to x and not from -x to x can someone help me fix this part for solve_ivp?
here is the code reduced to a min
from pylab import *
from scipy.integrate import solve_ivp
from scipy.optimize import brentq
import numpy as np
import itertools
a=1
B=4
L= B+a
Vmax= 50
Vpot = False
N = 1000 # number of points to take
psi = np.zeros([N,2]) # Wave function values and its derivative (psi and psi')
psi0 = array([0,1]) # Wave function initial states
Vo = 50
E = 0.0 # global variable Energy needed for Sch.Eq, changed in function "Wave function"
b = L # point outside of well where we need to check if the function diverges
x = linspace(-B-a, L, N) # x-axis
def V(x):
'''
#Potential function in the finite square well.
'''
if -a <=x <=a:
val = Vo
elif x<=-a-B:
val = Vmax
elif x>=L:
val = Vmax
else:
val = 0
return val
def SE(z, p):
state0 = p[1]
state1 = 1.0*(V(z) - E)*p[0]
return array([state0, state1])
def Wave_function(energy):
global E
E = energy
# odeint(func, y0, t)
# solve_ivp(fun, t_span, y0)
psi = solve_ivp(SE, [-B-a, L], psi0, max_step = ((B+a+L)/(N)))
return psi.y
def main():
# main program
f2 = figure(2)
plot(x, Wave_function(9.8)[0][:1000])
grid()
f2.show()
if __name__ == "__main__":
main()
this is what the code gives me in the end. the right side is okay but the left side is wrong. I depend on both sides working, not for visuals. 
from solve_ivp from -x to +x
No comments:
Post a Comment