6 views (last 30 days)

Show older comments

Dear All,

I have a function which depend on a parameter of an integral, like the following one

g(c) = integral(x^2 + c*x + 1)

where the integration range is [0,1] and the integration variable is x. In matlab this function can be defined as a multiple anonymous function

g = @(c) (integral(@(x) (x.^2 + c*x + 1),0,1));

Now, what I need to do is to integrate a function of g(c), let's say the fourth power, over c in the same range [0,1]. Hence I wrote the code

I = integral(@(c) g(c).^4,0,1)

but it doesn't work and the reason seems to be the inner product 'c*x'. Indeed I have got the same error even if I simply do the integral of g: I = integral(g,0,1) I can simply sample the c axis and use trapz routine instead of integral or implement other quadrature schemes, but still I would like to figure out why it does not work in this way and if I can overcome the problem. Any suggestion?

Thanks Nicola

Mike Hosea
on 18 Oct 2014

Edited: Mike Hosea
on 18 Oct 2014

So if your function g works with a scalar value of c, then you need to vectorize it. You'll have the same problem if you do something like

x = linspace(0,1);

plot(x,g(x))

The easiest way to vectorize g is with arrayfun:

gv = @(c)arrayfun(g,c);

then something like @(c)gv(c).^4 will be vectorized properly for plotting, integrating, or whatever.

NICOLA
on 20 Oct 2014

Mike Hosea
on 20 Oct 2014

For technical reasons, nested integration of smooth functions often results in extra accuracy, but nothing comes with out cost. You might want to loosen the tolerances a bit, say add

'AbsTol',1e-5,'RelTol',1e-3

to the integral2 and integral3 calls. Of course, if you are nesting integral3 and integral2, it means that you are ultimately doing the equivalent work of a 5-D integral. Nested adaptive quadrature is often the first-attempted method, since it is easy to try, but you might need to resort to sparse grid or Monte Carlo methods for better speed.

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!