Double integral with variable limits over implicit function in MATLAB

I'm experiencing the problem because of variable limit and implicit funtion being together.

So let's simplify it to this:

s(y)=y - our "implicit" function

Int [Int(x*s(y)*dy, 1,x)*dx, 1, 2] - our double integral (which equals 9/8).

(So you can even separate in into 2 integral I_small= s(y)dy and I=I_small * x*dx)

All that I figured out:

1) I tried using quad2d (so there is no ploblem with variable limit) - but I can't put the root of implicit function in it . So it wotks for non-implicit function:

 function main
     quad2d(@myfun,1,2,1,@(x)x)
 end
 function value=myfun(x,y)
    value=x.*y;
 end

But for implicit I've tried that - and it doesn't work. I know there is something wrong with this code - matlab doesn't understand that argument "y" in function name and "y" in the function itself are the same. But don't know how to fix it.

function main
     quad2d(@myfun,1,2,1,@(x)x)
 end
 function value=myfun(x,y)
    f=@(s,y)s-y;
    value=x.*fzero(@(s)f(y,s), 0.5);
 end

2) This code solves the opposite I = s(x).*y and I can't understand how to switch x to y because fzero doesnt work if I place y in it instead of x(j)

function main
    quad(@myfun, 0,1)
end
function z=myfun(x)
    z=zeros(size(x));
    f=@(x,s) s-x;
    for j=1:length(x);
         s(j)=fzero(@(s)f(x(j),s), 0.5);
         h=@(y) (s(j).*y);
         z(j)=quad(h,1,x(j));
    end
end

3) I also tried the nested quads, but it only works with constant limits. Can't fiqure it out how instead of Upperlimit should I place @(x)x.

function main
    quad(@(y)y.*first_int(2),1,2)
end
function value=first_int(UpperLimit)
    value=quad(@(x)yfunction(x,1),1,UpperLimit);
end
function value=yfunction(x,l)
    syms y;
    f=@(x,y) l.*x-y;
    for k=1:length(x)
        value(k)=fzero(@(y)f(x(k),y), 0.5);
end

Could you guys help with that?

Solution

The command quad2d (as its modern and better counterpart integral2) requires that the function to be integrated accept matrices as inputs.

The function Z=FUN(X,Y) must accept 2D matrices X and Y of the same size and return a matrix Z of corresponding values.

On the other hand, fzero only solves one equation, not a whole bunch of them at once. So, in order for your implicit function to accept matrix inputs, it needs to be written with a loop:

function value = myfun(x,y)
    f=@(s,y)s-y;
    value = zeros(size(x));   
    for i=1:size(x,1)
        for j=1:size(x,2)
            value(i,j) = x(i,j)*fzero(@(s) f(y(i,j),s), 0.5);
        end
    end
 end

Then quad2d(@myfun,1,2,1,@(x)x) or integral2(@myfun,1,2,1,@(x)x) will work.