Why does abs cause bad matlabfunction evaluation?

4 Ansichten (letzte 30 Tage)
Alec Jacobson
Alec Jacobson am 9 Dez. 2023
Bearbeitet: Catalytic am 9 Dez. 2023
If I run
x = sym('x','real');
y = sym('y','real');
f = abs(((x - 1)^2 + y^2)^(1/2) - 1)^2 + ((x^2 + y^2)^(1/2) - 1)^2;
d2fdx2 = diff(diff(f,x),x);
[eval(subs(subs(d2fdx2,x,1),y,1)) feval(matlabFunction(d2fdx2),1,1)]
Then I see
ans =
1.2929 NaN
If I remove the abs then I get the correct evaluation for both ([1.2929 1.2929]).
I'm wondering why this abs is causing trouble in matlabFunction . (in this minimal example it's of course easy to remove)
  2 Kommentare
Alec Jacobson
Alec Jacobson am 9 Dez. 2023
Symbolically we know that if x and y are real then ((x - 1)^2 + y^2)^(1/2) is real, but perhaps the symbolic toolbox is not acknowledging this?
Dyuman Joshi
Dyuman Joshi am 9 Dez. 2023
It is acknowledging it -
syms x y real
assumptions
ans = 
eqn = ((x - 1)^2 + y^2)^(1/2)
eqn = 
isAlways(eqn>=0)
ans = logical
1

Melden Sie sich an, um zu kommentieren.

Antworten (3)

Catalytic
Catalytic am 9 Dez. 2023
Bearbeitet: Catalytic am 9 Dez. 2023
The real problem appears to be that f is not fully simplified, and therefore the Symbolic Toolbox tries to naively apply the Chain rule to sub-expressions abs(z)^2. Even though abs(z)^2=z^2 is differentiable for real z, the chain rule isn't applicable because the nested function abs(z) is not continuously differentiable.
This causes the Toolbox to express the derivatives using Dirac deltas which evaluate to NaN at x=y=1 even though that is fake calculus -
x = sym('x','real');
y = sym('y','real');
f = abs(((x - 1)^2 + y^2)^(1/2) - 1)^2 + ((x^2 + y^2)^(1/2) - 1)^2
f = 
d2fdx2 = diff(diff(f,x),x) %Dirac delta in 6th term.
d2fdx2 = 
On the other hand, when f is properly simplified, the Chain rule works properly, because there are no abs(z)^2 sub-expressions -
fsimp=simplify( abs(((x - 1)^2 + y^2)^(1/2) - 1)^2 ) + ((x^2 + y^2)^(1/2) - 1)^2
fsimp = 
d2fdx2 = diff(diff(fsimp,x),x) %no Diracs
d2fdx2 = 
eval(subs(d2fdx2,[x,y],[1,1]))
ans = 1.2929

Walter Roberson
Walter Roberson am 9 Dez. 2023
You are taking the derivative of abs(), but abs() does not have a continuous derivative.
You have the sub-expression abs(((x - 1)^2 + y^2)^(1/2) - 1)^2 but as x and y both approach 1, then that goes to abs(0) and the derivative of abs at 0 is not continuous.
Depending on the order of substitution, as you partially evaluate the derivative, the Symbolic Toolbox might turn part of the expression into a dirac delta. But when that dirac delta is evaluated at 0 then the result is infinity; that infinity happens to be multiplied by a 0 in that case, and that leads to 0*inf which is NaN.
[eval(subs(subs(d2fdx2,x,1),y,1)) feval(matlabFunction(d2fdx2),1,1)]
There is no documented meaning for eval() of a symbolic expression. You will not find eval() of a symbolic expression described anywhere in doc .
Symbolic expressions are represented to users in a human-readable form rather than as MATLAB expressions or as expressions in the internal symbolic computing language MuPAD. eval() of them can fail.

Matt J
Matt J am 9 Dez. 2023
Bearbeitet: Matt J am 9 Dez. 2023
It is not the matlabFunction version that is evaluating things badly. It is the first version, with the double application of subs. The comparison that you should be doing is,
x = sym('x','real');
y = sym('y','real');
f = abs(((x - 1)^2 + y^2)^(1/2) - 1)^2 + ((x^2 + y^2)^(1/2) - 1)^2;
d2fdx2 = diff(diff(f,x),x);
[eval(subs(d2fdx2,[x,y],[1,1])) feval(matlabFunction(d2fdx2),1,1)]
ans = 1×2
NaN NaN
which is the correct result.
  2 Kommentare
Alec Jacobson
Alec Jacobson am 9 Dez. 2023
Oh, I didn't know it was an issue to apply `subs` twice. Why is that?
Matt J
Matt J am 9 Dez. 2023
Bearbeitet: Matt J am 9 Dez. 2023
When subs() is used, a ggeneric simplification is triggered, which is not valid in all cases. Consider the simpler example,
syms x y
expr=(x-1)/(y-1)
expr = 
This expression should be NaN at x=y=1. However, if you first evaluate at x=1, giving the expression 0/(1-y), the symbolic engine simplifies this to its generic value, which is 0,
step1=subs(expr,x,1)
step1 = 
0
and so the second substitution step of course results in 0 as well.
step2=subs(step1,y,1)
step2 = 
0
Conversely,
subs(expr,[x,y],[1,1])
Error using sym/subs
Division by zero.

Melden Sie sich an, um zu kommentieren.

Produkte


Version

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by