-1

I am receiving a ValueError when integrating using scipy.integrate.quad. Here is my simplified code:

import numpy as np
import scipy.integrate as integrate
p = np.arange(0,1,1/1000)
h = lambda p: p**2/2+p*(1-p)
Kl = lambda p: h(p) + 0.02
K = Kl(p)
R = 0.5*h(p) + 0.5*h(1)
Vl = lambda p: np.minimum.reduce([p, K, R])
integrate.quad(Vl, 0, 1)[0]

The last line gives the exception:

ValueError: setting an array element with a sequence.

Can someone please propose how to do this integration properly? Thanks

splinter
  • 3,727
  • 8
  • 37
  • 82
  • Could you explain your `Vl` function? What do you want `np.minimum.reduce([p, K, R])` to do? – DSM Feb 20 '17 at 18:12
  • Thanks @DSM. I would like it to take the element-wise minimum of p, K and R. I think of it as a mathematical function Vl(p,K(p),R(p)) – splinter Feb 20 '17 at 18:15
  • This is the same code as http://stackoverflow.com/questions/42342409/valueerror-when-defining-a-lambda-function-in-python; just a slightly different question. – hpaulj Feb 20 '17 at 20:07

1 Answers1

1

I think this code confuses functions on floats with arrays. The function argument to quad needs to accept a float and return a float. In

Vl = lambda p: np.minimum.reduce([p, K, R])

p is a float (and is not the p array-- it's usually a bad idea to have a module global with the same name as an argument) and K and R are arrays, which isn't what we want (that's what's throwing the error.)

I think you're just looking for

def h(p):
    return p**2/2+p*(1-p)

def K(p):
    return h(p) + 0.02

def R(p):
    return 0.5*h(p) + 0.5*h(1)

def Vl(p):
    return min(p, K(p), R(p))

which gives me

In [177]: integrate.quad(Vl, 0, 1.0)
Out[177]: (0.34689543041336846, 4.8736376714649885e-09)
DSM
  • 342,061
  • 65
  • 592
  • 494