# python Simpson integration [closed]

Currently, I am trying to apply a Crank-Nicolson method on a function that I want to evolve. However, I am only facing one drawback and it is the step when I normalize the initial function using the scipy integrate function.

It tells me that it runs into a

runtime error: invalid value encountered in dive

This only happens when I change the sigma parameter in the function (it works when sigma=0.05 but I require to be sigma=5e-7). If anyone could give me some insight about this I will appreciate it. Here is the relevant part of the code:

import scipy as sp
import numpy as np
from scipy import integrate, sparse, linalg,special
import scipy.sparse.linalg
import pylab as pl
from multiprocessing import Pool
import multiprocessing

def crank_nicol(tf):
ax=-10. #punto inicial
bx=10. #punto final
nx= 500#numero de puntos en el x-grid
Lx=bx-ax #longitud total
delx= Lx/nx #x-step size=Lx/nx
dt= 1.e-3#T/mt # time step delta t
mt=int((2*np.pi*tf/dt)/400)
nux= 1j*dt/(4.0*(delx**2)) #nu q son iguales a mis c_i=b_i

gridx = sp.zeros(nx) #space grid
igridx = sp.array(range(nx)) # time grid
psi = sp.zeros(nx)
pot = sp.zeros(nx)
gridx = delx*(igridx - nx/2)
x=gridx
##############
sigma=0.05
alpha=2.66
R_ave=2.75
r_e=.969
xmin=(1/2.66)*np.arccosh(0.5*np.exp(2.66*1.012/2.0))
####################################################
psi = sp.pi**(-0.25)*sp.exp(-(0.5/sigma**2)*(gridx+xmin)**2) #-1j*3*(gridx-10))

psi /= sp.sqrt(sp.integrate.simps(psi*np.conjugate(psi), dx=delx))# THIS IS THE PART WITH THE ERROR

print(nux)

return x, psi*np.conjugate(psi),psi

times=np.linspace(10,200,20)
if __name__ == '__main__':
n_times = 10
pool = multiprocessing.Pool(processes=20)
results = pool.map(crank_nicol, times)


## closed as off-topic by Kirill, GoHokies, Christian Clason, nicoguaro♦Jan 15 at 1:00

This question appears to be off-topic. The users who voted to close gave this specific reason:

If this question can be reworded to fit the rules in the help center, please edit the question.

• You expect to get a sensible value for exp(-1e11)? Sum up the vector of zeros and then divide by this zero result? – LutzL Jan 12 at 11:45
• is that I don't understand why it gives me a vector of zeros when I only change the sigma value? I believe that python should be able to perform the calculationeven if the number is very small (if this is not the case, please explain me the reason why) – Oliver Jan 14 at 20:16
• Because exp(-1e11) has about 1e10 zeros after the decimal point. The floating point representation of that result can only be 0. Do not forget, if sigma is very small, 0.5/sigma**2 is very large. – LutzL Jan 14 at 20:24