Different behaviour when using symbolic simplification on norm() for scalars and vectors?

I'm building quite a lengthy differential equation (through many substitutions, etc) that involves a norm of the difference between two vectors. I am then trying to use odeToVectorField() and matlabFunction() to convert it into an integration problem. I run into the following error:

Error using symengine
Unable to generate code for piecewise for use in anonymous functions.
Error in sym/matlabFunction>mup2mat (line 431)
    res1 = mupadmex('symobj::generateMATLAB',r.s,ano,spa,splitIt);

I don't think I can write the equation out here because it is far, far too long, but while trying to debug the issue, I found that if I change the following term (which is used several times for different eta and a in constructing the equation):

    d_eta = norm(eta-a); %Both eta and a are defined as real vectors, e.g. syms eta [N 1] real

into

    d_eta = ( ((eta(1)-a(1))^2) + ((eta(2)-a(2))^2) ) ^ (0.5);

then matlabFunction() is able to finish execution and the subsequent integration works fine (well, at least it runs, I'm still trying to come up with a verification case). This is obviously not an ideal solution as I'd like to be able to handle the n-dimensional case, so I'm trying to understand what is happening under the hood in the norm function.

I think the abs() function that is used in the norm() function is what is causing the issue - MATLAB is trying to create piecewise functions for it. What is the best way to have MATLAB understand that eta and a are real, and it can remove the abs() function?

For instance, why does the following work in simplifying norm:

>> syms a b real
>> c = norm(a-b)
c =
(abs(a - b)^2)^(1/2)
>> simplify(c)
ans =
abs(a - b)

But not this?:

>> syms a [2 1] real
>> syms b [2 1] real
>> c = norm(a-b)
c =
(abs(a1 - b1)^2 + abs(a2 - b2)^2)^(1/2)
>> simplify(c)
ans =
(abs(a1 - b1)^2 + abs(a2 - b2)^2)^(1/2)

However, this works:

>> syms a [2 1] real
>> syms b [2 1] real
>> syms c
>> c = norm(a-b)
c =
(abs(a1 - b1)^2 + abs(a2 - b2)^2)^(1/2)
>> rewrite(c,"sqrt")
ans =
((a1 - b1)^2 + (a2 - b2)^2)^(1/2)

Would it be "safe" to use rewrite to overcome this issue? Are there any nuanced differences between these (particularly which could impact the subsequent integration)?

Answer 1

[Torsten]

If you want to define the norm of the difference of two real-valued n-dimensional vectors, use

norm_x_minus_y = sqrt(sum((x-y).^2))

Then you can avoid the trouble with abs.

I don't know why MATLAB uses abs although you defined a and b as real.

[Imran Khan]

Thanks for your reply @Torsten. I know that I can proceed in a number of ways to workaround this, including just defining the vector norm explicitly like you have done (or even component-wise, because I know the max number of dimensions I'm likely to have). But my question was more about why this somewhat strange behaviour is happening, and was wondering what the root cause was - is this a bug? Or is the symbolic toolbox not able to perform a simplification keeping the real assumption "in mind" because the norm function is a matlab defined function (but seemingly the rewrite function does let it do that)? Or perhaps some other nuanced reason.

[Torsten]

But my question was more about why this somewhat strange behaviour is happening, and was wondering what the root cause was - is this a bug?

It's not a bug because using abs is not wrong. I don't know why there is no simplification although you specified a and b as real. I suggest you ask MATLAB support - only the developers of the symbolic toolbox know.

https://uk.mathworks.com/support/contact_us.html

But I think there are many other cases in symbolic computations where you will find "simplifications" where the toolbox doesn't.

[Imran Khan]

Well yes, I suppose you're right, it's not strictly speaking a bug because its not wrong. Interesting case of MATLAB not doing anything wrong, and the coder not doing anything wrong, but still an error occuring (which can happen! :-)

Thanks for your suggestion - I'll reach out to MATLAB support and see if they have any feedback about it.

Answer 2

[Walter Roberson]

Error using symengine

Unable to generate code for piecewise for use in anonymous functions.

The work-around for that is to use the 'file' option to request that the generated code be written to a file. The piecewise() call will then be converted to if/elseif

Note: doing this has the disadvantage that the code cannot be vectorized over the expression that is piecewise()'d over -- all the inputs that go into making up that expression must be scalar. The generated code is a simple if/elseif and does not use logical indexing to selectively operate on portions of the expression.

Also, if you use the 'file' option then be careful about the 'optimize' option. When you use 'file', the default 'optimize' is true. Unfortunately, historically 'optimize' was broken (I don't know if it has been repaired yet.) I am fairly sure it 'optimize' was broken in the R2022a time-frame.

[Imran Khan]

I didn't even know that matlabFunction() can be used to output the function to file - though the doc does mention it now that I looked. Thanks for the suggestion @Walter Roberson and the caution around the function's vectorization!