How to find double integral in Matlab to the following code

I have the following code, I have a problem to solve it, please i need help

syms s t
a=0.1;b=0.4;
Rs=-log(s);
Ws=vpa([int(log(s)*log(s)*Rs,s,a,b), -int(log(s)*Rs,s,a,b), -2*int(log(s)*cos(2*pi*s)*Rs,s,a,b); ...
    -int(log(s)*Rs,s,a,b), int(Rs,s,a,b), 2*int(cos(2*pi*s)*Rs,s,a,b); ...
    -2*int(log(s)*cos(2*pi*s)*Rs,s,a,b), 2*int(cos(2*pi*s)*Rs,s,a,b), 4*int(cos(2*pi*s)*cos(2*pi*s)*Rs,s,a,b)]);
W11=inv(Ws);
w1 = W11(:,1);
g(s)=[-log(s);1;2*cos(2*pi*s)];
g(t)=[-log(t);1;2*cos(2*pi*t)];
Gs=((w1)'*g(s));
Gt=((w1)'*g(t));
fin = @(s,t) Gs.*Gt.*Rs.*Rt.*(min(1-s,1-t)-((1-s).*(1-t)));
MM=integral2(fin,0.1,0.4,0.1,0.4)

The problem in the min(1-s,1-t)-((1-s).*(1-t))

Thank you in advance

[Walter Roberson]

g(s)=[-log(s);1;2*cos(2*pi*s)];
g(t)=[-log(t);1;2*cos(2*pi*t)];

That redefines g using a different variable name but otherwise the same formula

Gs=((w1)'*g(s));
Gt=((w1)'*g(t));

You seem to expect them to be different ?

When you define a symbolic function, it is like a regular function, in that the variable you use will be replaced by whatever is passed to it, so you could define g just once but call it with s and t... however the only difference between the results would be in the variable used in the formula.

[Walter Roberson]

You do not define Rt

[norh hameed]

Thank you for your comment, Rt=-log(t). But the problem in fin = @(s,t) Gs.*Gt.*Rs.*Rt.*(min(1-s,1-t)-((1-s).*(1-t)))

MM=integral2(fin,0.1,0.4,0.1,0.4)

I cant calculate this integral.

Many thanks

Answer 1

[Walter Roberson]

syms s t v
Q = @(v) sym(v);
a = Q(0.1);
b = Q(0.4);
Rs = -log(s);
Rt = -log(t);
Pi = Q(pi);
Int = @vpaintegral;
Ws = ([Int(log(s)*log(s)*Rs,s,a,b), -Int(log(s)*Rs,s,a,b), -2*Int(log(s)*cos(2*Pi*s)*Rs,s,a,b); ...
    -Int(log(s)*Rs,s,a,b), Int(Rs,s,a,b), 2*Int(cos(2*Pi*s)*Rs,s,a,b); ...
    -2*Int(log(s)*cos(2*Pi*s)*Rs,s,a,b), 2*Int(cos(2*Pi*s)*Rs,s,a,b), 4*Int(cos(2*Pi*s)*cos(2*Pi*s)*Rs,s,a,b)]);
W11 = inv(Ws);
w1 = W11(:,1);
g(v) = [-log(v);1;2*cos(2*Pi*v)];
Gs=((w1)'*g(s));
Gt=((w1)'*g(t));

Min = @(a,b) piecewise(a<=b, a, b);

fin = Gs.*Gt.*Rs.*Rt.*(Min(1-s,1-t)-((1-s).*(1-t)));

MM = Int(Int(fin,s,Q(0.1),Q(0.4)),t,Q(0.1),Q(0.4));

disp(MM)

Your Gs and Gt are separable and could be calculated independently if not for the (Min(1-s,1-t)-((1-s).*(1-t))) part .

[norh hameed]

My Matlab is R2016a. When I apply the above code I got the following error

Undefined function or variable 'vpaintegral'

Undefined function or variable 'piecewise'.

Thank you

[Walter Roberson]

vpaintegral() was introduced in R2016b. If you could upgrade then you could get a significant performance improvement for this operation.

IsOld = verLessThan('matlab', '9.1');    %before R2016b

syms s t v
Q = @(v) sym(v);
a = Q(0.1);
b = Q(0.4);
Rs = -log(s);
Rt = -log(t);
Pi = Q(pi);

if IsOld
    Int = @int;
    symLIST = @(varargin)feval(symengine,'DOM_LIST',varargin{:});
    piecewise = @(cond, val1, val2) feval(symengine, 'piecewise', symLIST(cond, val1), symLIST(~cond, val2));
else
    Int = @vpaintegral;
end

Ws = ([Int(log(s)*log(s)*Rs,s,a,b), -Int(log(s)*Rs,s,a,b), -2*Int(log(s)*cos(2*Pi*s)*Rs,s,a,b); ...
    -Int(log(s)*Rs,s,a,b), Int(Rs,s,a,b), 2*Int(cos(2*Pi*s)*Rs,s,a,b); ...
    -2*Int(log(s)*cos(2*Pi*s)*Rs,s,a,b), 2*Int(cos(2*Pi*s)*Rs,s,a,b), 4*Int(cos(2*Pi*s)*cos(2*Pi*s)*Rs,s,a,b)]);
W11 = inv(Ws);
w1 = W11(:,1);
g(v) = [-log(v);1;2*cos(2*Pi*v)];
Gs=((w1)'*g(s));
Gt=((w1)'*g(t));

Min = @(a,b) piecewise(a<=b, a, b);

fin = Gs.*Gt.*Rs.*Rt.*(Min(1-s,1-t)-((1-s).*(1-t)));

if IsOld
    fin = vpa(fin);
end

MM = Int(Int(fin,s,Q(0.1),Q(0.4)),t,Q(0.1),Q(0.4));

MMd = double(MM);

disp(MMd)

[norh hameed]

Thank you so much, but it takes a long time to apply

[Walter Roberson]

Yes it does. If you were to upgrade to even just one release after what you have now, you would get much faster results.

You could probably improve performance by splitting into two cases at the boundary of the two expressions being min(). If you were to do that then you could use matlabFunction() and integral2()

[norh hameed]

Thank you so much.