piecewise for use in anonymous functions
Mostrar comentarios más antiguos
Hello,
clc,clear all,close all;
opengl('save', 'software')
syms beta theta r ;
B = 10;
H = 10;
C = 10;
L1 = 4;
Ri = B/2 + L1;
Rf = C + B/2 + L1;
Rbeta0 = Ri;
theta1 = atan((H/2)/(B/2+L1));
theta2 = pi - atan((H/2)/(B/2-L1));
Rt0 = piecewise((0<=theta) & (theta<theta1),(B/2+L1)/cos(theta),(theta1<=theta) & (theta<theta2),(H/2)/sin(theta),(theta2<=theta) & (theta<=pi),(B/2-L1)/cos(pi-theta));
Rmaxtb0 = (Rt0*Rbeta0)/Ri;
Rmaxtb00 = matlabFunction(Rmaxtb0)
Vmbeta0 = (2*(pi^2)*((B+(2*L1))^2)) / ((((2*L1*pi)+(pi*B))^2));
Vbeta0 = Vmbeta0 * (1 - ((r^2) / (Rmaxtb0^2)));
Vbeta00 = Vbeta0 * r;
Vbeta000 = matlabFunction(Vbeta00,'Vars',[theta r]);
Wsigmat = 2 * integral2(Vbeta000,0,pi,0,Rmaxtb00,'AbsTol',1e-5, 'RelTol',1e-5)
Error using symengine
Unable to generate code for piecewise for use in anonymous functions.
As you can see in my piecewise function there is 'theta' in the limits which is defined as 'syms' and is not defined prior
and I think it couldn't be written as if/else!! (I don't know if it can)
Is there any way to define a piecewise function as an anonymous function?
or any way to use an alternative for piecewise function??
thanks in advance for any help
Respuesta aceptada
Walter Roberson
el 6 de Ag. de 2020
Rmaxtb00 = matlabFunction(Rmaxtb0, 'File', 'Rmaxtb00.m', 'optimize', false);
The generated function in the .m file will use if/elseif to code the piecewise()
However, the generated if/elseif code will not be vectorized, so when you use it with integral() you need to use 'ArrayValued', true so that integral() does not pass in a vector of values to be evaluated.
8 comentarios
Pooyan Jafari
el 6 de Ag. de 2020
Walter Roberson
el 6 de Ag. de 2020
ArrayValued is for integral() not integral2()
Pooyan Jafari
el 8 de Ag. de 2020
Editada: Pooyan Jafari
el 8 de Ag. de 2020
Walter Roberson
el 9 de Ag. de 2020
The most immediate problem is that the way you have defined Rt0 with piecewise, when diff(Vbeta,theta) is taken, the derivative is not defined at 0, so you get a NaN being generated (because piecewise gives nan for any value out of range.)
That is something you can define away by changing to
Rt0 = piecewise((-pi/2 < theta) & (theta<theta1),(B/2+L1)/cos(theta),(theta1<=theta) & (theta<theta2),(H/2)/sin(theta),(theta2<=theta) & (theta<=pi),(B/2-L1)/cos(pi-theta));
This extends the first case across 0, allowing the derivative to exist at 0 exactly.
That done, you still have problems. The (extended) integral does not converge at r = 0, theta = 0, beta = 0.
temp = epsmaxxx(sym('r'),0,0);
limit(temp, r, 0, 'left')
limit(temp, r, 0, 'right')
MATLAB will not be able to find those two limits, but if you put the expression through another package such as Maple, you will find that the left limit is -infinity and the right limit is +infinity .
The message about reaching the maximum number of intervals is generated while integrating with respect to r, when beta and theta are both small and positive, such as
theta: 1.71978371262361e-06
beta: 8.59891856311806e-07
By that point, the function being integrated has reached about 2E+17 and getting larger rapidly as r gets smaller.
Because it is increasing in value so quickly, and the values are all positive at the time the problem is detected, we know that this is not just a removable discontinuity being dealt with, and that exact way that we extended Rt0 to be continuous across 0 is not the cause of the problem (because theta is never negative in the calculation.)
So, I don't think you can caclulate the integral -- not unless you put in non-zero boundary conditions.
Pooyan Jafari
el 9 de Ag. de 2020
Editada: Pooyan Jafari
el 9 de Ag. de 2020
Pooyan Jafari
el 9 de Ag. de 2020
Walter Roberson
el 9 de Ag. de 2020
test = diff(Vbeta,theta)
test =
piecewise(theta in Dom::Interval(0, 6257646083598987/9007199254740992), etc
Notice that the Interval given starts with 0. We are looking at MuPad internal notation here, and in that notation there is a difference between the scalar 0 and the "list" [0] . MATLAB treats 0 and [0] the same (except that the [] involves an extra call). However MuPad does not treat them the same: 0 is DOM::Numeric and [0] is DOM::List which is a data structure. In the context of DOM:Interval, an endpoint that numeric is Open Interval, and an endpoint that is DOM::List is closed on that side.
Thus the closed interval endpoint in Vbeta includes 0 and you can subs(Vbeta, theta, 0) and get something useful. However when you differentiate the fact that Vbeta is not defined before 0 means that the derivative at 0 exactly is undefined. When you
subs(diff(Vbeta, theta), theta, 0)
then you get nan because 0 is outside any interval set up by the piecewise that remains after the derivative.
since the denominator cos(theta) is nonzero at 0, then it is not required that we introduce a singularity at 0 in the derivative; we just extend the piecewise slightly so that the derivative exists in theory at 0, and proceed knowing that we will be able to evaluate the derivative at 0 exactly in the piecewise without triggering nan.
Walter Roberson
el 9 de Ag. de 2020
"I think because Rmaxtb is in the denominator of Vbeta, there shouldn't be any problem regarding differentiation (I am not sure);"
You are missing that when you calculate a derivative then for any given answer to considered correct, it must be valid on every compact infinitesimal interval near the point. However when you have a function with a closed boundary and define the value outside the interval to be nan (because that is what piecewise does), then although the derivative might be correct on "one side" of the boundary, on the other side of the boundary the derivative is undefined. There is no non-zero infinitesimal that you can go negative to, and at which a hypothetical derivative is defined.
Can we potentially match endpoints between two piecewise that are defined on different sides of the boundary, piecewise one way and piecewise the other way agreeing on a value where they meet? That is a valid approach for some cases, but this case says that to the left (negative) the function value is undefined (nan). No possible meeting of functions.
subs theta=0 into the diff(Vbeta, theta) and you will get nan, for the reasons described. So you cannot integrate the triple integral on any interval that includes theta = 0.
... not unless you extend the piecewise in a way that is consistent, such as just moving the lower boundary from zero to something negative but greater than π/2
Más respuestas (0)
Categorías
Más información sobre Mathematics en Centro de ayuda y File Exchange.
el 5 de Ag. de 2020
el 9 de Ag. de 2020
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
