'integral' with piecewise expressions

I have the following expression for the angle Beta:

L = 1;
syms y 
Beta = piecewise(0<=y<L/2, pi/2, ...
    L/2 <=y <L*(3/4),pi/2+0.3491, ...
    (3/4)*L <=y <=L,pi/2-0.3491);

i want to compute this integral:

    a = pi/4;

    fun_f=@(y) (1./(50*(1+y)*(cos(a)*cos(Beta)+sin(a)*sin(Beta))+470));
     
    int= integral(fun_f,0,L)

Error using integralCalc/finalInputChecks
Input function must return 'double' or 'single' values. Found 'sym'.

Error in integralCalc/iterateScalarValued (line 315)
finalInputChecks(x,fx);

Error in integralCalc/vadapt (line 132)
[q,errbnd] = iterateScalarValued(u,tinterval,pathlen);

Error in integralCalc (line 75)
[q,errbnd] = vadapt(@AtoBInvTransform,interval);

Error in integral (line 87)
Q = integralCalc(fun,a,b,opstruct);

Does anybody know how to calculate the integral with using integral, but with Beta defined as a piecewise function?

I looked for documentation for piecewise but I don't really know if it's needed to use it. I'd prefer not to.

If you need further information i'll be happy to provide it in order to solve my problem

Answer 1

[Torsten]

L = 1;
Beta = @(y,L) pi/2.*((0<=y) & (y<L/2)) + (pi/2+0.3491).*((L/2<=y) & (y<3/4*L)) + (pi/2-0.3491).*((3/4*L<=y) & (y<=L));
figure(1)
plot((0:0.01:L),Beta(0:0.01:L,L))

a = pi/4;
fun_f=@(y,L) (1./(50*(1+y).*(cos(a)*cos(Beta(y,L))+sin(a)*sin(Beta(y,L)))+470));  
figure(2)
plot((0:0.01:L),fun_f(0:0.01:L,L))

format long
int_value= integral(@(y)fun_f(y,L),0,L)

int_value = 
   0.001918690636037

int_value_improved = integral(@(y)fun_f(y,L),0,L/2) + integral(@(y)fun_f(y,L),L/2,3*L/4) + integral(@(y)fun_f(y,L),3/4*L,L)

int_value_improved = 
   0.001918689980576

[anto]

I tried implementing it on my code: but it gives me the error : 'Incorrect number or types of inputs or outputs for function 'cos'.'

LATO = 1;
intersezione_34LATO= (3/4)*LATO;
i = 1;
p = [0.125000000000000	0];
L_AB= 0.625000000000000;
alpha_1 = pi/2;
        Beta = @(l) (pi/2+0.3491).*(LATO/2<=(p(i,1)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i)) <intersezione_34LATO)) + (pi/2-0.3491).*(intersezione_34LATO<=(p(i,2)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i))<= LATO));

        fun_f=@(l) (1./(cos(Beta(l))));

        tau_f(i) = integral(@(l)fun_f(l),0,L_AB(i),'AbsTol',0,'RelTol',1e-20,'ArrayValued',false); 

Incorrect number or types of inputs or outputs for function 'cos'.

Error in solution (line 9)
fun_f=@(l) (1./(cos(Beta(l))));

Error in solution (line 11)
tau_f(i) = integral(@(l)fun_f(l),0,L_AB(i),'AbsTol',0,'RelTol',1e-20,'ArrayValued',false);

Error in integralCalc/iterateScalarValued (line 314)
fx = FUN(t);

Error in integralCalc/vadapt (line 132)
[q,errbnd] = iterateScalarValued(u,tinterval,pathlen);

Error in integralCalc (line 75)
[q,errbnd] = vadapt(@AtoBInvTransform,interval);

Error in integral (line 87)
Q = integralCalc(fun,a,b,opstruct);

I even tried with the syntax for Beta you suggested:

LATO = 1;
intersezione_34LATO= (3/4)*LATO;
i = 1;
p = [0.125000000000000	0];
L_AB= 0.625000000000000;
alpha_1 = pi/2;
   
        
        Beta = @(l,p,alpha_1,LATO) (pi/2+0.3491).*(LATO/2<=(p(i,1)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i)) <intersezione_34LATO)) + (pi/2-0.3491).*(intersezione_34LATO<=(p(i,2)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i))<= LATO));

                fun_f=@(l) (1./(cos(Beta(l,p,alpha_1,LATO))));
                
                  tau_f(i) = integral(@(l)fun_f(l),0,L_AB(i),'AbsTol',0,'RelTol',1e-20,'ArrayValued',false);

this is part of a for cycle. but the syntax is as you suggested, but for some reasons it still gives:

Incorrect number or types of inputs or outputs for function 'cos'.

Do you know why? Thanks

[Walter Roberson]

Your Beta is returning logical values. You probably have () in the wrong place.

[anto]

The expression for Beta is OK, also the parenthesis in the fun_f expression are OK though

[Torsten]

I doubt that Beta is what you want since the result is of type "logical".

What function do you want to use (in a mathematical notation) ?

But if you think everything is as wanted - here is the result:

LATO = 1;
intersezione_34LATO= (3/4)*LATO;
i = 1;
p = [0.125000000000000	0];
L_AB= 0.625000000000000;
alpha_1 = pi/2;


Beta = @(l,p,alpha_1,LATO) (pi/2+0.3491).*(LATO/2<=(p(i,1)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i)) <intersezione_34LATO)) + (pi/2-0.3491).*(intersezione_34LATO<=(p(i,2)+l.*sin(alpha_1(i)))) & ((p(i,2)+l.*sin(alpha_1(i))<= LATO));

fun_f=@(l) (1./(cos(double(Beta(l,p,alpha_1,LATO)))));
l=0:0.01:L_AB;
plot(l,fun_f(l))

tau_f(i) = integral(@(l)fun_f(l),0,L_AB(i))%,'AbsTol',0,'RelTol',1e-20,'ArrayValued',false)

tau_f = 0.8377

[anto]

OK, thanks a lot. With double it works, have a good day

Answer 2

[Walter Roberson]

Use matlabFunction with the piecewise expression, giving the 'file' option and 'optimize' false. And when you integral specify 'arrayvalued' true

matlabFunction can convert piecewise to if/else but only when writing to file, and the result cannot accept vectors

[Walter Roberson]

format long g
L = 1;
syms y 
Pi = sym(pi);
Beta = piecewise(0<=y<L/2, Pi/2, ...
    L/2 <= y < L*(3/4), Pi/2 + sym(3491)/10^4, ...
    (3/4)*L <= y <=L, Pi/2 - sym(3491)/10^4);

a = Pi/4;

fun_f = (1./(50*(1+y)*(cos(a)*cos(Beta)+sin(a)*sin(Beta))+470))

fun_f = 

     
result_symbolic = int(fun_f, y, 0, L)

result_symbolic = 

result_vpa = vpa(result_symbolic, 16)

result_vpa = 

0.001918689980575797

fun_f_h = matlabFunction(fun_f, 'vars', y, 'File', 'fun_f.m', 'optimize', false)

fun_f_h = function_handle with value:
    @fun_f

result_numeric = integral(fun_f_h, 0, L, 'arrayvalued', true)

result_numeric = 
       0.00191869063603731

dbtype fun_f.m

1     function fun_f = fun_f(y)
2     %FUN_F
3     %    FUN_F = FUN_F(Y)
4     
5     %    This function was generated by the Symbolic Math Toolbox version 9.2.
6     %    14-Jan-2023 22:18:56
7     
8     if ((y < 1.0./2.0) & (0.0 <= y))
9         fun_f = 1.0./((sqrt(2.0).*(y.*5.0e+1+5.0e+1))./2.0+4.7e+2);
10    elseif ((y < 3.0./4.0) & (1.0./2.0 <= y))
11        fun_f = 1.0./((y.*5.0e+1+5.0e+1).*((sqrt(2.0).*cos(pi./2.0+3.491e-1))./2.0+(sqrt(2.0).*sin(pi./2.0+3.491e-1))./2.0)+4.7e+2);
12    elseif ((y <= 1.0) & (3.0./4.0 <= y))
13        fun_f = 1.0./((y.*5.0e+1+5.0e+1).*((sqrt(2.0).*cos(pi./2.0-3.491e-1))./2.0+(sqrt(2.0).*sin(pi./2.0-3.491e-1))./2.0)+4.7e+2);
14    else
15        fun_f = NaN;
16    end

[anto]

Thanks a lot