Odeint to solve Mukhanov equation

[scidam]

As far as I understood, the problem is in lines 38-41. You can resolve it, e.g. assuming that ep = H,

solutions = []
for ep, H in zip(range(len(t)), range(len(t))):    
    asol = odeint(sol, Y0, t, args=(k, ep, H), h0 = 0.005)  
    solutions.append(asol)

Another approach would be to use dictionary to store results:

solutions = dict()
ep_values = [1,2,3]
H_values =[4,3,1]
for ep in ep_values:
    for H in H_values:
        asol = odeint(sol, Y0, t, args=(k, ep, H), h0 = 0.005)  
        solutions.update({(ep, H): asol})

But this will require to restructure the rest part of your code.