Operations and differentiation with symmatrix in MATLAB... possible bug?

Hi

I have found matlab to behave unexpectedly when using the symmatrix. I dont understand it, maybe i misunderstood the way operators work, or maybe its a bug?

Any help would be greatly appreciated! Thank you!

I declare some symbolic objects using symmatrix

clear all
n = 2;
A = symmatrix('A', n);
lambda = symmatrix('lambda', [n,1]);
g = symmatrix('g', [n,1]);
phi = symmatrix('phi', [n,1]);
l = symmatrix('l', [n,1]);
item = symmatrix('l', [n,1]);

I declare 2 almost identical equations

arga = (l.'*(((lambda).*(A * g)  + (phi.'*lambda) * (A * g)))); 
argb = (l.'*(((lambda).*(item)   + (phi.'*lambda) * (A * g))));

arga and argb should be the interchangable, as both item and A*g are 2x1 parameter vectors.

I differenciate wrt lambda:

sola = diff(arga, lambda) 
solb = diff(argb, lambda)

I get

sola =
 
l.'*(symmatrix(eye(2)) .* A*g + kron(phi.', A*g))
 
 
solb =
 
l.'*(kron(phi.', A*g) + symmatrix(eye(2)) .* l)

ProbIem 1: I cannot evaluate the sum in the parenthesis of sola due to dimensionality issues. the probelm doesnt appear for solb.

Problem 2: Running the numerical example below, sola and solb give two results different whereas this should not be the case

l=[1.2;1.1] ;
g=[1.9876;1.88] ;
phi=[1.0987;1.5192] 
A=[132 123;1222 124] ;
item=A * g;
% sola, having removed only the symmatrix comand
l.'*((eye(2)) .* A*g + kron(phi.', A*g))
% solb, having removed only the symmatrix comand
l.'*(kron(phi.', A*g) + (eye(2)) .* l)

I get

ans =
1.0e+03 *
4.4392    5.9196

ans =
1.0e+03 *
3.8694    5.3495

Answer 1

[Paul]

Hi dominik,

I'm not sure what Problem 1 is refering to.

I think Problem 2 can be explained as follows:

clear all
n = 2;
A = symmatrix('A', n);
lambda = symmatrix('lambda', [n,1]);
g = symmatrix('g', [n,1]);
phi = symmatrix('phi', [n,1]);
l = symmatrix('l', [n,1])

l = 

item = symmatrix('l', [n,1])

item = 

We see that the Matlab variables l and item both reference the same symbolic variable.

arga = (l.'*(((lambda).*(A * g)  + (phi.'*lambda) * (A * g))))

arga = 

argb = (l.'*(((lambda).*(item)   + (phi.'*lambda) * (A * g))))

argb = 

In argb, we see that the symbolic l shows up in the expression on the right.

sola = diff(arga, lambda)

sola = 

solb = diff(argb, lambda)

solb = 

l=[1.2;1.1] ;
g=[1.9876;1.88] ;
phi=[1.0987;1.5192] ;
A=[132 123;1222 124] ;

Here, item is set to A*g, but item is not a variable in the expression for solb. So setting item here has no affect on the results and because l ~= A*g we get different results for sola and solb.

%item=A * g;
% sola, having removed only the symmatrix comand

Also, I think A*g must be enclosed in parentheses.

l.'*((eye(2)) .* (A*g) + kron(phi.', A*g))

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

% solb, having removed only the symmatrix comand
l.'*(kron(phi.', A*g) + (eye(2)) .* l)

ans = 1x2

1.0e+03 * 3.8694 5.3495

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Maybe the intent was for item to be its own variable

%clear all
n = 2;
A = symmatrix('A', n);
lambda = symmatrix('lambda', [n,1]);
g = symmatrix('g', [n,1]);
phi = symmatrix('phi', [n,1]);
l = symmatrix('l', [n,1]);

item = symmatrix('Q', [n,1])

item = 

arga = (l.'*(((lambda).*(A * g)  + (phi.'*lambda) * (A * g)))); 
argb = (l.'*(((lambda).*(item)   + (phi.'*lambda) * (A * g))));
sola = diff(arga, lambda) 

sola = 

solb = diff(argb, lambda)

solb = 

l=[1.2;1.1] ;
g=[1.9876;1.88] ;
phi=[1.0987;1.5192] ;
A=[132 123;1222 124] ;
%item=A * g;
Q = A * g;

% sola, having removed only the symmatrix comand
l.'*((eye(2)) .* (A*g) + kron(phi.', A*g))

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

% solb, having removed only the symmatrix comand
l.'*(kron(phi.', A*g) + (eye(2)) .* Q)

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

[dominik]

Dear Paul,

thank you for your extremely swift and comprehensive response.

I indeed had a typo in my code, but you figured that out yourself already... Sorry, and thank you!

The second half of your answer regards the problem which I did mean to post. And indeed, if you enclose A*g in parentheses as you suggest, the problem is solved.

However, why do I have to manually add parentheses? Why are those parentheses not in the symbolic expression that matlab displays? Here is what I get in 2023a, as a screenshot:

The parentheses are missing (marked in red).

Are symbolic matrix operators (their order of execution) not consistent with numeric matrix operators? (I only used symbolic scalars before and i cant recall having noted the same issue).

And if this is the case, is there a way to translate the symbolic expression into the equivalent numeric expression?

Problem 1, which i hadn't explained well enough, is bascially the same as problem 2. If I excecute the highlighted expression in the comand window, right after having execuded the code in thescreeshot, i get an error of matching dimensions. I.e. matlab doesnt accept its own symbolic output (without brackets) as symbolic input.

Thank you very much!

Best

Dominik

PS: How do you get this beautiful math font output?

[Paul]

Hi dominik,

How do you get this beautiful math font output?

Use LiveScript, which I think is also what's being used here on Answers.

Are symbolic matrix operators (their order of execution) not consistent with numeric matrix operators?

Not that I know of, and I would be absolutely shocked if that were the case.

Let's start with the code

n = 2;
A      = symmatrix('A', n);
lambda = symmatrix('lambda', [n,1]);
g      = symmatrix('g',   [n,1]);
phi    = symmatrix('phi', [n,1]);
l      = symmatrix('l',   [n,1]);
item   = symmatrix('Q',   [n,1]);

arga = (l.'*(((lambda).*(A * g)  + (phi.'*lambda) * (A * g)))); 
argb = (l.'*(((lambda).*(item)   + (phi.'*lambda) * (A * g))));

sola = diff(arga, lambda) 

sola = 

solb = diff(argb, lambda)

solb = 

BTW, could the display of the Kronecker term also be ambiguous? It wouldn't be in Matlab because kron is a function, not an operator. But in math, which the display is trying to show, is it always the the case that (regular) product has precedence over the tensor product? Or are they both considered multiplication with precedence from left to right in the usual way?

However, why do I have to manually add parentheses? Why are those parentheses not in the symbolic expression that matlab displays?

No idea. Clearly, these two expressions are different

e1 = symmatrix(eye(2)).*(A*g);
e2 = symmatrix(eye(2)).*A*g;
size(e1)

ans = 1x2

2 2

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

size(e2)

ans = 1x2

2 1

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Interestingly, even though the LiveScript displays of e1 and e2 are different, both expressions are the same using the normal, Matlab operator precedence rules (assuming that bull's eye symbol means .* , though I doubt that's standard math notation, is .* a formal math operator?)

e1,e2

e1 = 

e2 = 

And they are displayed exactly the same using pretty (I don't know if there is a way to get the normal "screen output" here on Answers.

pretty(e1)

I  .* A g
 2

pretty(e2)

I  .* A g
 2

is there a way to translate the symbolic expression into the equivalent numeric expression?

There are a couple of options.

One option is to use subs to substitute numeric values for symbolic variables if you just want the answer

lnum   = [1.2;1.1] ;
gnum   = [1.9876;1.88] ;
phinum = [1.0987;1.5192] ;
Anum   = [132 123;1222 124] ;
subs(sola,{A l g phi},{Anum, lnum,gnum,phinum})

ans = 

Use vpa to get the numbers, but still in sym world

vpa(ans)

ans = 

Use double to convert to numeric

double(ans)

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

If don't care about overwriting variables, then

l   = [1.2;1.1] ;
g   = [1.9876;1.88] ;
phi = [1.0987;1.5192] ;
A   = [132 123;1222 124] ;
double(subs(sola)) % no other arguments needed

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Another option is to convert sola to a matlabFunction (doesn't work with symmatrix, so need to convert to sym)

ha = matlabFunction(symmatrix2sym(sola))

ha = function_handle with value:
    @(A1_1,A1_2,A2_1,A2_2,g1,g2,l1,l2,phi1,phi2)[l1.*(A1_1.*g1+A1_2.*g2+phi1.*(A1_1.*g1+A1_2.*g2))+l2.*phi1.*(A2_1.*g1+A2_2.*g2),l2.*(A2_1.*g1+A2_2.*g2+phi2.*(A2_1.*g1+A2_2.*g2))+l1.*phi2.*(A1_1.*g1+A1_2.*g2)]

ha(Anum(1,1),Anum(1,2),Anum(2,1),Anum(2,2),gnum(1),gnum(2),lnum(1),lnum(2),phinum(1),phinum(2))

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

It appears that the actual processing of sola is correct based on comparison to solb. Are you able to verify that the expression for sola (with the parentheses on A*g) is, in fact, the correct derivative for your expression? I have no reason to doubt that it is.

Perhaps the only problem is how sola is displayed, though I can see that being a pretty serious problem if you want to copy/paste the expression for sola from the command window.

[dominik]

Dear Paul,

thanks again for taking the time to answer. What you wrote at the very end goes to the core of my problem. Let me give you some background. I want to use symbolic math to derive a system of (matrix) equations that i will later solve numercially. In the above code sola or solb represent this system of equations.

I want my system to be evaluated using numercial matrix operations for computational efficiency reasons, since n will be very large. Hence converting symmatrix to sym and using matlabfunction wont do the trick.

So I was hoping to be able to copy-paste from the comand window, or equivalently convert the symbolic expression to char and then process it further into a code that evaluates the equation system. This is basically a workaround for the shortcoming that matlabfunction does not accept symmatrix. But now it turns out that copy-pasting the symbolic output and using it for numeric calculations is not feasible either.

Got to run now, but will get back to this later. Thanks again!

[Paul]

Perhaps it would be worthwhile to file a bug report with a simple example to illustrate the problem with the display to the screen. Maybe use the e1 and e2 example from above. If you do file a bug report, please add a comment on this thread, or to the original question, with their response.

In the near term, all might not be lost. Starting with the code:

n = 2;
A      = symmatrix('A', n);
lambda = symmatrix('lambda', [n,1]);
g      = symmatrix('g',   [n,1]);
phi    = symmatrix('phi', [n,1]);
l      = symmatrix('l',   [n,1]);
item   = symmatrix('Q',   [n,1]);
arga = (l.'*(((lambda).*(A * g)  + (phi.'*lambda) * (A * g)))); 
argb = (l.'*(((lambda).*(item)   + (phi.'*lambda) * (A * g))));
sola = diff(arga, lambda);

Now, if we know the variables in sola, we can make matlabFunction accept matrix arguments, like so

funcsola = matlabFunction(symmatrix2sym(sola),'Vars',{symmatrix2sym(A), symmatrix2sym(g), symmatrix2sym(phi), symmatrix2sym(l)})

funcsola = function_handle with value:
    @(in1,in2,in3,in4)[in4(1,:).*(in1(1).*in2(1,:)+in1(3).*in2(2,:)+in3(1,:).*(in1(1).*in2(1,:)+in1(3).*in2(2,:)))+in4(2,:).*in3(1,:).*(in1(2).*in2(1,:)+in1(4).*in2(2,:)),in4(2,:).*(in1(2).*in2(1,:)+in1(4).*in2(2,:)+in3(2,:).*(in1(2).*in2(1,:)+in1(4).*in2(2,:)))+in4(1,:).*in3(2,:).*(in1(1).*in2(1,:)+in1(3).*in2(2,:))]

lnum   = [1.2;1.1] ;
gnum   = [1.9876;1.88] ;
phinum = [1.0987;1.5192] ;
Anum   = [132 123;1222 124] ;

funcsola(Anum,gnum,phinum,lnum)

ans = 1x2

1.0e+03 * 4.4603 8.2765

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Unfortunately, the signature of funcsola is opaque in terms of its variable names. There might be a way to modify that, I'm not sure.

If more flexibility is needed to generate funcsola, like if you don't know its variables a priori, there may be other options, but I'm not quite sure what you're looking for.

[dominik]

Dear Paul, thanks again! I will hand in a bug report.

What you suggest as a short term solution generates a function that is suitable for numeric evaluation. But the problem is that this code relies on elementwise calculations of the result vector, which should be slow for large n. Furethermore, n may not be ex ante fixed, and i wouldnt want to have to rerun the symbolic code eevry time i change n.

[Paul]

To make sure I understand correctly ....

sola is derived as symmatrix with n = 2 with the idea that the expression for sola is vaiid for any other value of n, and then we want to define a Matlab expression for sola that can be evaluated for numerical values of the variables using standard matrix operations with the variables being of dimension for any value of n?

[dominik]

exactly

[Paul]

AFAIK there's no way to do that (but will keep thinking about it), maybe someone else will chime in with workable approach.

It would be cool to be able to use symmatrix with symbolic dimensions, and be able to matlabFucntion the same, something like

syms n m integer
A = symmatrix('A',[n n]);
B = symmatrix('B',[n m]);
C = A*B;
cfunc = matlabFunction(C);

Maybe you should consider submitting an enhancement request. They do improve functionality, e.g., diff on symmatrix arrived in r2023b.

Did you get a reponse on your bug report?