How to evaluate a symbolic expression having `max` and `diff`?

I have calculated the jacobian of two functions where variables are x1, x2, x3.

The jacobian is as follows-

JacobianF =
 
[                                            diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) + 96*pi*cos(6*pi*x1)*(x3 + sin(6*pi*x1)) + 160*3^(1/2)*pi^2*cos(6*pi*x1)*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3) + 1,                                                                                                                                                                                        0, 16*x3 + 16*sin(6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x3) + (80*3^(1/2)*pi*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3))/3]
[diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) - 96*pi*cos((2*pi)/3 + 6*pi*x1)*(x2 - sin((2*pi)/3 + 6*pi*x1)) - 240*2^(1/2)*pi^2*cos((2*pi)/3 + 6*pi*x1)*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2) - 1, 16*x2 - 16*sin((2*pi)/3 + 6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x2) + 40*2^(1/2)*pi*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2),                                                                                                                                                                      0]

Now, I need to evaluate this JacobianF at

X = [0.2703    0.6193    0.9370];

where X(1) is x1 and so on.

To evaluate this JacobianF, I have used the following code-

Var_List = sym('x', [1, 3]);
df=double(subs(JacobianF, Var_List, X));

However, I get the following error. What is the cause of this error? How to resolve it and calculate the JacobianF at the specified position?

Error using symengine
Unable to convert expression containing remaining symbolic function calls into double array. Argument must be
expression that evaluates to number.

Error in sym/double (line 872)
Xstr = mupadmex('symobj::double', S.s, 0);

[Dyuman Joshi]

"What is the cause of this error?"

You can't use double() on a symolic expression -

syms x1 x2 x3

JacobianF = [diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) + 96*pi*cos(6*pi*x1)*(x3 + sin(6*pi*x1)) + 160*3^(1/2)*pi^2*cos(6*pi*x1)*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3) + 1,                                                                                                                                                                                        0, 16*x3 + 16*sin(6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x3) + (80*3^(1/2)*pi*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3))/3;
diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) - 96*pi*cos((2*pi)/3 + 6*pi*x1)*(x2 - sin((2*pi)/3 + 6*pi*x1)) - 240*2^(1/2)*pi^2*cos((2*pi)/3 + 6*pi*x1)*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2) - 1, 16*x2 - 16*sin((2*pi)/3 + 6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x2) + 40*2^(1/2)*pi*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2),                                                                                                                                                                      0]

JacobianF = 

X = [0.2703    0.6193    0.9370];

df= subs(JacobianF, [x1 x2 x3], X)

df = 

[Rounak Saha Niloy]

I understand.

But I need to have the evaluated values of each component of JacobianF as a matrix. How can I get that in this case?

[Dyuman Joshi]

As for why the expression is not evaluated to a numeric value, I am not so sure.

If I had to guess, I'd say it is because you have assumed x1, x2 and x3 to be not-real i.e. complex, yet you have provided a real value.

But then, the Sym engine is unable to evaluate the expression for complex numbers as well -

syms x1 x2 x3

JacobianF = [diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) + 96*pi*cos(6*pi*x1)*(x3 + sin(6*pi*x1)) + 160*3^(1/2)*pi^2*cos(6*pi*x1)*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3) + 1,                                                                                                                                                                                        0, 16*x3 + 16*sin(6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x3) + (80*3^(1/2)*pi*sin((3^(1/2)*pi*(20*x3 + 20*sin(6*pi*x1)))/3))/3;
diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x1) - 96*pi*cos((2*pi)/3 + 6*pi*x1)*(x2 - sin((2*pi)/3 + 6*pi*x1)) - 240*2^(1/2)*pi^2*cos((2*pi)/3 + 6*pi*x1)*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2) - 1, 16*x2 - 16*sin((2*pi)/3 + 6*pi*x1) + diff(max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real')), x2) + 40*2^(1/2)*pi*sin((2^(1/2)*pi*(20*x2 - 20*sin((2*pi)/3 + 6*pi*x1)))/2),                                                                                                                                                                      0];

X = [0.2703    0.6193    0.9370];
Y = X + X*i;

%Real
df1 = subs(JacobianF, [x1 x2 x3], X);
vpa(df1)

ans = 

%Complex
df2 = subs(JacobianF, [x1 x2 x3], Y);
vpa(df2)

ans = 

[Rounak Saha Niloy]

x1, x2 and x3 are assumed as real and the values provided are also real.

[Walter Roberson]

syms x real
F = diff(max([0,x]),x)

F = 

subs(F, x, -1)

ans = 

double(ans)

Error using symengine
Unable to convert expression containing remaining symbolic function calls into double array. Argument must be expression that evaluates to number.

Error in sym/double (line 761)
Xstr = mupadmex('symobj::double', S.s, 0);

[Walter Roberson]

syms x real
F = diff(piecewise(x <= 0, 0, x),x)

F = 

subs(F, x, -1)

ans = 

0

double(ans)

ans = 0

So your problem is the diff(max) and you can avoid it by rewriting as piecewise.

[Rounak Saha Niloy]

Yes, I am rewriting the problems as piecewise. I will get back once evaluated and verified. Thanks again.

[Rounak Saha Niloy]

I have sorted the max function. But I am having trouble with min functions.

I have the following function which needs to be converted-

function Output = b_flat(y,A,B,C)
    Output = A+min(0,floor(y-B))*A.*(B-y)/B-min(0,floor(C-y))*(1-A).*(y-C)/(1-C);
    Output = round(Output*1e4)/1e4;
end

I have converted it as-

function Output = b_flat(y,A,B,C)
    Output = A+piecewise(0<=floor(y-B),0,floor(y-B))*A.*(B-y)/B-piecewise(0<=floor(C-y),0,floor(C-y))*(1-A).*(y-C)/(1-C);
    Output = round(Output*1e4)/1e4;
end

But this gives me the following error. What might be causing the issue?

Output = A+piecewise(0<=floor(y-B),0,floor(y-B))*A.*(B-y)/B-piecewise(0<=floor(C-y),0,floor(C-y))*(1-A).*(y-C)/(1-C);

Error using mupadengine/feval_internal

First argument must be a condition.

Error in sym/piecewise (line 49)

pw = feval_internal(symengine, 'piecewise', lists{:}, 'ExclusiveConditions');

Error in CalculateJacobianF>b_flat (line 364)

Output = A+piecewise(0<=floor(y-B),0,floor(y-B))*A.*(B-y)/B-piecewise(0<=floor(C-y),0,floor(C-y))*(1-A).*(y-C)/(1-C);

[Walter Roberson]

You would encounter that problem if both y and B or both y and C are numeric. piecewise() is only defined when the condition is symbolic

[Rounak Saha Niloy]

in this case, y is symbolic and B,C are numeric.

[Torsten]

Use

min(x,0) = 0.5*(x-abs(x))

as I used

max(x,0) = 0.5*(x+abs(x))

below.

[Walter Roberson]

Looks like it works for me when y is symbolic.

syms y
b_flat(y, 1, 2, 3)

ans = 

b_flat(y, -10, 5, 17)

ans = 

function Output = b_flat(y,A,B,C)
    Output = A+piecewise(0<=floor(y-B),0,floor(y-B))*A.*(B-y)/B-piecewise(0<=floor(C-y),0,floor(C-y))*(1-A).*(y-C)/(1-C);
    Output = round(Output*1e4)/1e4;
end

Answer 1

[Walter Roberson]

The derivative of max() is not generally defined.

You would probably have more success if you defined in terms of piecewise() instead of in terms of max()

[Dyuman Joshi]

I guess the Sym engine does not have the ability to recognise that the definition of max() can be broken into a piecewise definition, than a derivative can be calculated.

I wonder if that is possible to implement or not.

Answer 2

[Torsten]

Use

max(x,0) = 0.5*(abs(x)+x)

for real x.

syms x1
f1 = max([0, (7*sin(4*pi*x1))/10], [], 2, 'omitnan', ~in(x1, 'real'))

f1 = 

f2 = 0.5*(abs(7*sin(4*pi*x1)/10)+7*sin(4*pi*x1)/10)

f2 = 

figure(1)
hold on
fplot(f1,[-0.5 0.25])
fplot(f2,[-0.5 0.25])
hold off

df1 = diff(f1,x1)

df1 = 

df2 = diff(f2,x1)

df2 = 

figure(2)
%fplot(df1,[-0.5 0.25])
fplot(df2,[-0.5 0.25])