Solving a double integral equation for a third var

问题:

I've been trying to solve a double integral equation for a third variable in matlab.

An example:

At first, I tried to solve it symbolically (when k=1) as below:

syms x y h
F = @(x,y,h) 2*x*y+4*h; 
f = @(x) x/2;
solve(int(int(F(x,y,h)*f(x),x,0,3)*f(y),y,0,1)-3, h)

The code gives the right answer, i.e. h=2/3 when k=1.

For my real problem, however, the functions F and f are so much more complex. And when I applied the same code above with the complex F and f, the matlab does not solve it in an appropriate time. Actually I do not know if it ever solves - I have let the code run for 30 mins and forced it to terminate. I'll have to pursue further this route but for now, I'm trying to solve it numerically. The code below is what I have tried:

F = @(x,y,h) 2.*x.*y+4.*h; 
f = @(x) x./2;
g1 = @(y,h) integral(@(x) F(x,y,h).*f(x),0,3)
g2 = @(h) integral(@(y) g1(y,h).*f(y),0,1)-3
bsolve = fsolve(g2,0)

Why does this code give me the wrong answer of 0.5833?

回答1:

From the documentation for integral:

For scalar-valued problems the function Y = FUN(X) must accept a vector argument X and return a vector result Y, the integrand function evaluated at each element of X. For array-valued problems (see the 'ArrayValued' option below) FUN must accept a scalar and return an array of values.

Because you're nesting calls to integral, the outer integral is passing vectors into the inner call. To fix this, the inner integral, g1 must configured to operate on scalars:

g1 = @(y,h)integral(@(x)F1(x,y,h).*f1(x),0,3,'ArrayValued',true)

The code will now return 2/3 in floating point. You should probably also use fzero instead of fsolve for this univariate root-finding problem.