Why Does fplot() Show a Phantom Pole?

Example of fplot() showing a pole where clearely one does not exist. I know I can get rid of the dased line with the ShowPoles option, but that would eliminate all of the vertical lines, even for actual poles should there be any. Any idea why fplot() can't handle such a seemingly simple function?

syms t real
s(t) = piecewise(t<-1,0, t>2,0, exp(-abs(t)))

s(t) = 

figure;
fplot(s(t),[-3 3])

[Walter Roberson]

Interestingly, plotting [-2.89 2.9] does not produce the line, but [-2.9 2.9] does.

[Paul]

That is weird. I wonder why the x-limits of the plot influence whehter or not some underlying algorithm finds a pole between those limits.

Unfortunately, this problem is coming up for me quite a bit in some work I've been doing lately and I'm getting tired of fighting with fplot(). In some cases, I've gone to using the matlabFunction/plot() sequence, which is o.k. when it's possible (and it's not possible with piecewise, at least in 2019a), but it's annoying, especially for plots of functions that do have singularities.

[Walter Roberson]

The trick for using matlabFunction() with piecewise is:

1. ask matlabFunction() to write to a file using the 'File' option; and then
2. Use arrayfun() in invoking the resulting function, so that at the point it would need to do a piecewise(), only a single value would be processed
3. As you are using R2019a time frame, when you use matlabFunction() with 'File' option, be sure to use 'Optimize', false . The optimization phase malfunctioned in that time frame, returning the results of intermediate computations instead of final solutions. (Unfortunately, the optimization phase got worse in some of the later releases.)

[Paul]

I was wondering if the File option might do the trick, as it seems there would be more flexibility there than with an anonymous function. Thanks for pointing this out. And I'm soon to move to 2020b, so I'll keep in mind the warning about not using Optimize.

[Walter Roberson]

matlabFunction() challenges:

• R2020a onward: if you matlabFunction() something that has an int() in it, then the variable of integration should be "bound" within the int(), effectively variable "local" to that expression, and should not appear outside the equivalent integrate() unless possibly as a named parameter for an anonymous function. But in R2021b, sub-expression optimization can be fooled into thinking it is an input variable that can be moved outside of the int() .

For example, if you had int(4*x^2 + x*cos(x^2), x, -pi, pi) then if optimization is on in matlabFunction(), it might decide to create

t1 = x.^2; output = integrate(@(x) 4*t1 + x*cos(t1), x, -pi, pi)

which is obviously wrong.

• matlabFunction() does not work on vpaintegral() (which was introduced in R2016b)
• matlabFunction() does not work for the int() 'hold' option (which was introduced in R2019b)
• matlabFunction() R2021a fails with 'file' option if the problem involves more than one variable (this was perhaps fixed in an Update, as it was an easy fix)
• matlabFunction() R2021a with 'file' and 'optimize' can produce very inefficient code for expressions whose optimization calls for more than 500 temporary variables
• matlabFunction() R2021a with 'file' and 'optimize' can generate incorrect code, in a situation that is especially triggered by piecewise() with at least two cases when total length of all expressions is more than about 128 characters, or if some of the outputs are pure numeric instead of being symbolic expressions, or even without piecewise() if arrays of outputs are being created and some of the entries are numeric but not all are
• matlabFunction() R2020b (and possibly earlier): order of variables inserted as parameters is not as documented
• matlabFunction() R2020b (and possibly earlier): fails to generate the correct size of output in certain cases that mathematically evaluate to constants
• matlabFunction() R2020a, R2020b: int() or symsum() raised to power will generate error
• matlabFunction() R2020a: fails on besselh()
• matlabFunction() R2019a: can fail on piecewise() when 'optimize' is set, returning intermediate caclulations (this problem was fixed at some point.)

[Paul]

Wow. Have these all been reported as bugs?

Do you have examples at hand for either of these two:

• matlabFunction() R2020b (and possibly earlier): order of variables inserted as parameters is not as documented
• matlabFunction() R2020b (and possibly earlier): fails to generate the correct size of output in certain cases that mathematically evaluate to constants

If not, please don't spend any effort on it.

[Walter Roberson]

@Paul

I pulled all of those from bug reports I have sent to Mathworks.

The first one, my report was:

===

matlabFunction 'vars' documents "By default, when you convert symbolic expressions, the order is alphabetical"

However, that is not the case.

matlabFunction(str2sym('A+A1+A2+A10+a+a1+a2+a10+B+b+x+y+z+X+Y+Z+A1_1+A1_10+A10_1+A10_10+A2_1+A1_2')) 

ans = function_handle with value:
    @(A,A1,A2,A10,A10_1,A10_10,A1_1,A1_2,A2_1,A1_10,B,X,Y,Z,a,a1,a2,a10,b,x,y,z)A+A1+A2+A10+A10_1+A10_10+A1_1+A1_2+A2_1+A1_10+B+X+Y+Z+a+a1+a2+a10+b+x+y+z

ans =

function_handle with value:

@(A,A1,A2,A10,A10_1,A10_10,A1_1,A1_2,A2_1,A1_10,B,X,Y,Z,a,a1,a2,a10,b,x,y,z)A+A1+A2+A10+A10_1+A10_10+A1_1+A1_2+A2_1+A1_10+B+X+Y+Z+a+a1+a2+a10+b+x+y+z

So all capital letters are before all lower-case letters, and the first level sorting is according to the first letter.

For the variables that have numbers without underscore, the sort is numeric rather than alphabetic (alphabetic would require that A10 be before A2 because '1' is before '2'.) But notice that A10_1 is before A1_1 and that A1_10 is after A2_1 and A10_10 is before A1_10

I think you would be hard-pressed to call it "alphabetical order".

You can potentially explain why A10_1 might appear before A2_1 but you cannot use alphabetical order to explain why A1_10 appears after A2_1

[Walter Roberson]

For the second of those:

===

Reference:

https://www.mathworks.com/matlabcentral/answers/694280-error-calculating-an-integral-input-function-must-return-double-or-single-values-found-sym#answer_578865

Suppose you have a symbolic formula that you want to turn into an anonymous function, such as

syms x 
A = randi([-2 2], 1, 5); 
for K = 1 : 5; F{K} = matlabFunction(A(K)*x, 'vars', {x}); end 
for K = 1 : 5; result(K) = integral(F{K}, 0, 1); end 

This works fine, calculating the integrals as needed. Unless, that is, that one of the A entries happens to be 0.

Because then A(K)*x -> 0*x -> 0 (symbolic) and

matlabFunction(sym(0), 'Vars', {x}) 

will generate

@(x) 0.

which is a problem because integral() requires that the result of executing the function on a vector of values is to return a vector the same size, but the result of executing @(x) 0. on a vector of values is to return the scalar constant 0.

The workaround is to know that use 'ArrayValued', true in integral -- which also reduces efficiency for all the integrations. If your symbolic expression does not use int() or piecewise() then matlabFunction will probably be able to vectorize the calculation, and that can be important for calculation purposes.

The improvement needed for matlabFunction is that functions of a single variable but constant value, should automatically use

zeros(size(Variable)) or constant*ones(size(Variable)) .

If the user did not specify 'vars' and the expression has no variables, then returning the scalar is perhaps still warranted.

If the function has more than one variable and the user specified 'vars', with {} with a single entry, and that single entry is a row vector, then matlabFunction is going to vectorize along columns, so the return should be constant*ones(size(Variable,1),1) . Likewise if the single entry is a column vector then matlabFunction is going to vectorize along rows, so the return should be

constant*ones(1,size(Variable,2)) .

Other 'vars' with cell cases lead to ambiguity about what to vectorize over, I think.

If the function has more than one variable and the user specified 'vars' without a cell, or if the user did not specify 'vars', then you have potential ambiguity about what size to return. For example, if the code would have normally been

@(x,y) constant .* x .* y

except that the constant was 0 so matlabFunction received only a 0 to know about, then it could be incorrect to use either x or y as the size, because of implicit expansion. You could make guesses such as

@(x,y) constant .* ones(max(size(x), size(y)))

but by that point it might be best to simply add an additional option indicating which variable to extract the size information from.

[Walter Roberson]

The second one about incorrect size also refers to some of the details about the problems when using piecewise() with to a File with optimization turned on, which is a messy complicated one.

Answer 1

[Chunru]

The default number of t points (23) is not big enough to show the details around t=0. Change it to a larger value.

syms t real
s(t) = piecewise(t<-1,0, t>2,0, exp(-abs(t)))

s(t) = 

fplot(s(t),[-3 3], 'MeshDensity', 200)

[Paul]

I certainly like this approach more than using ShowPoles, and it has the added bonus of the plot showing the full range of the function up to s(0) = 1. Still makes me wonder how fplot() determines pole locations.

[Chunru]

it "adaptively" setting the evaluation points based on the defaults. I am afraid there is no fool-proof methods for determining the pole locations for all possible scenarios.

[Paul]

MeshDensity seems to be working well for the cases I'm interested in. Thanks.