Subtraction of two symbolic matrices gives error "Array sizes must match" when the sizes are the same

I am attempting to confirm that two ways of manipulating symbolic variables are equivalent. The variables ogLambda and ogN are the original method that I know is correct. The variables newLambda and newN are the new formulations that I am trying to confirm. The comparison of ogN and newN works as expected, outputting a value of 1 from the isequal call, however ogLambda and newLambda are acting strangely. When I display both variables they appear identical, and when I used isequal and nested for loops to check every row and column index independently the isequal command returns 1 (for example isequal(ogLambda(i,j),newLambda(i,j)) outputs 1). However when comparing both matrices in their entierty isequal returns 0. Even more confusingly I cannot subtract the two martices to compare their difference because when I attempt to I am met with the following error:

Error using -

Array sizes must match.

This is even more confusing because when I call size(ogLambda) and size(newLambda) both outputs are 3 3. If someone could help me to understand what is going on here I would be very grateful as this is related to my thesis which is due in a month.

rng(3)
L=randi([1 10],3,3,4,4);
N=randi([1 10],3,4);
syms xvec [4,1] real
syms xmat [4,4] real
ogLambda=zeros(3,3);
for ind1=1:4
    for ind2=1:4
        ogLambda=ogLambda+L(:,:,ind1,ind2)*xmat(ind1,ind2);
    end
end
ogN=zeros(3,1);
for ind1=1:4
    ogN=ogN+N(:,ind1)*xvec(ind1);
end

newN=N*xvec;
expandxmat=permute(repmat(xmat,[1,1,3,3]),[3,4,1,2]);
newLambda=sum(sum(L.*expandxmat,3),4);

N_equivalence=isequal(ogN,newN)

N_equivalence = logical
   1

Lambda_equivalence=isequal(ogLambda,newLambda)

Lambda_equivalence = logical
   0

size(ogLambda)

ans = 1×2

3 3

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

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

size(newLambda)

ans = 1×2

3 3

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

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

ogLambda-newLambda

Error using - (line 7)
Array sizes must match.

>> ogLambda

ogLambda =

[6*xmat1_1 + 7*xmat1_2 + 8*xmat1_3 + 10*xmat1_4 + 5*xmat2_1 + 9*xmat2_2 + 6*xmat2_3 + 2*xmat2_4 + 3*xmat3_1 + 5*xmat3_2 + 8*xmat3_3 + 10*xmat3_4 + 3*xmat4_1 + 6*xmat4_2 + 2*xmat4_3 + 3*xmat4_4, 6*xmat1_1 + 3*xmat1_2 + xmat1_3 + 5*xmat1_4 + 7*xmat2_1 + 3*xmat2_2 + 7*xmat2_3 + 3*xmat2_4 + 7*xmat3_1 + 8*xmat3_2 + 8*xmat3_3 + xmat3_4 + 10*xmat4_1 + xmat4_2 + 8*xmat4_3 + 4*xmat4_4, 2*xmat1_1 + 3*xmat1_2 + 2*xmat1_3 + xmat1_4 + 6*xmat2_1 + 3*xmat2_2 + 3*xmat2_3 + 5*xmat2_4 + 6*xmat3_1 + 4*xmat3_2 + 7*xmat3_3 + 6*xmat3_4 + 9*xmat4_1 + 3*xmat4_2 + 6*xmat4_3 + 3*xmat4_4]

[ 8*xmat1_1 + 6*xmat1_2 + 5*xmat1_3 + 8*xmat1_4 + xmat2_1 + 6*xmat2_2 + 2*xmat2_3 + 5*xmat2_4 + 5*xmat3_1 + 4*xmat3_2 + 5*xmat3_3 + 5*xmat3_4 + 4*xmat4_1 + 4*xmat4_2 + 2*xmat4_3 + 7*xmat4_4, 9*xmat1_1 + 5*xmat1_2 + xmat1_3 + 3*xmat1_4 + 3*xmat2_1 + 5*xmat2_2 + 6*xmat2_3 + 8*xmat2_4 + 5*xmat3_1 + 7*xmat3_2 + 7*xmat3_3 + xmat3_4 + 7*xmat4_1 + 3*xmat4_2 + 10*xmat4_3 + xmat4_4, 3*xmat1_1 + 3*xmat1_2 + 2*xmat1_3 + 9*xmat1_4 + xmat2_1 + 6*xmat2_2 + 10*xmat2_3 + 3*xmat2_4 + 8*xmat3_1 + 7*xmat3_2 + 4*xmat3_3 + 9*xmat3_4 + 4*xmat4_1 + 7*xmat4_2 + xmat4_3 + 8*xmat4_4]

[ 3*xmat1_1 + 4*xmat1_2 + 6*xmat1_3 + 5*xmat1_4 + 5*xmat2_1 + 3*xmat2_2 + 10*xmat2_3 + xmat2_4 + 3*xmat3_1 + 8*xmat3_2 + 9*xmat3_3 + xmat3_4 + 10*xmat4_1 + 5*xmat4_2 + 8*xmat4_3 + 9*xmat4_4, 9*xmat1_1 + 5*xmat1_2 + 2*xmat1_3 + 2*xmat1_4 + 7*xmat2_1 + 3*xmat2_2 + 8*xmat2_3 + 7*xmat2_4 + 2*xmat3_1 + 7*xmat3_2 + 10*xmat3_3 + 10*xmat3_4 + 10*xmat4_1 + 10*xmat4_2 + 9*xmat4_3 + 9*xmat4_4, xmat1_1 + 5*xmat1_2 + 2*xmat1_3 + 7*xmat1_4 + 6*xmat2_1 + xmat2_2 + 9*xmat2_3 + 7*xmat2_4 + 4*xmat3_1 + 4*xmat3_2 + 6*xmat3_3 + 2*xmat3_4 + xmat4_1 + 7*xmat4_2 + 5*xmat4_3 + 2*xmat4_4]

>> newLambda

newLambda =

[6*xmat1_1 + 7*xmat1_2 + 8*xmat1_3 + 10*xmat1_4 + 5*xmat2_1 + 9*xmat2_2 + 6*xmat2_3 + 2*xmat2_4 + 3*xmat3_1 + 5*xmat3_2 + 8*xmat3_3 + 10*xmat3_4 + 3*xmat4_1 + 6*xmat4_2 + 2*xmat4_3 + 3*xmat4_4, 6*xmat1_1 + 3*xmat1_2 + xmat1_3 + 5*xmat1_4 + 7*xmat2_1 + 3*xmat2_2 + 7*xmat2_3 + 3*xmat2_4 + 7*xmat3_1 + 8*xmat3_2 + 8*xmat3_3 + xmat3_4 + 10*xmat4_1 + xmat4_2 + 8*xmat4_3 + 4*xmat4_4, 2*xmat1_1 + 3*xmat1_2 + 2*xmat1_3 + xmat1_4 + 6*xmat2_1 + 3*xmat2_2 + 3*xmat2_3 + 5*xmat2_4 + 6*xmat3_1 + 4*xmat3_2 + 7*xmat3_3 + 6*xmat3_4 + 9*xmat4_1 + 3*xmat4_2 + 6*xmat4_3 + 3*xmat4_4]

[ 8*xmat1_1 + 6*xmat1_2 + 5*xmat1_3 + 8*xmat1_4 + xmat2_1 + 6*xmat2_2 + 2*xmat2_3 + 5*xmat2_4 + 5*xmat3_1 + 4*xmat3_2 + 5*xmat3_3 + 5*xmat3_4 + 4*xmat4_1 + 4*xmat4_2 + 2*xmat4_3 + 7*xmat4_4, 9*xmat1_1 + 5*xmat1_2 + xmat1_3 + 3*xmat1_4 + 3*xmat2_1 + 5*xmat2_2 + 6*xmat2_3 + 8*xmat2_4 + 5*xmat3_1 + 7*xmat3_2 + 7*xmat3_3 + xmat3_4 + 7*xmat4_1 + 3*xmat4_2 + 10*xmat4_3 + xmat4_4, 3*xmat1_1 + 3*xmat1_2 + 2*xmat1_3 + 9*xmat1_4 + xmat2_1 + 6*xmat2_2 + 10*xmat2_3 + 3*xmat2_4 + 8*xmat3_1 + 7*xmat3_2 + 4*xmat3_3 + 9*xmat3_4 + 4*xmat4_1 + 7*xmat4_2 + xmat4_3 + 8*xmat4_4]

[ 3*xmat1_1 + 4*xmat1_2 + 6*xmat1_3 + 5*xmat1_4 + 5*xmat2_1 + 3*xmat2_2 + 10*xmat2_3 + xmat2_4 + 3*xmat3_1 + 8*xmat3_2 + 9*xmat3_3 + xmat3_4 + 10*xmat4_1 + 5*xmat4_2 + 8*xmat4_3 + 9*xmat4_4, 9*xmat1_1 + 5*xmat1_2 + 2*xmat1_3 + 2*xmat1_4 + 7*xmat2_1 + 3*xmat2_2 + 8*xmat2_3 + 7*xmat2_4 + 2*xmat3_1 + 7*xmat3_2 + 10*xmat3_3 + 10*xmat3_4 + 10*xmat4_1 + 10*xmat4_2 + 9*xmat4_3 + 9*xmat4_4, xmat1_1 + 5*xmat1_2 + 2*xmat1_3 + 7*xmat1_4 + 6*xmat2_1 + xmat2_2 + 9*xmat2_3 + 7*xmat2_4 + 4*xmat3_1 + 4*xmat3_2 + 6*xmat3_3 + 2*xmat3_4 + xmat4_1 + 7*xmat4_2 + 5*xmat4_3 + 2*xmat4_4]

>>

Answer 1

[Paul]

Hi Jacob,

Start of code

rng(3)
L=randi([1 10],3,3,4,4);
N=randi([1 10],3,4);
syms xvec [4,1] real
syms xmat [4,4] real
ogLambda=zeros(3,3);
for ind1=1:4
    for ind2=1:4
        ogLambda=ogLambda+L(:,:,ind1,ind2)*xmat(ind1,ind2);
    end
end
ogN=zeros(3,1);
for ind1=1:4
    ogN=ogN+N(:,ind1)*xvec(ind1);
end
newN=N*xvec;
expandxmat=permute(repmat(xmat,[1,1,3,3]),[3,4,1,2]);

Problem appears to be with this line, where newLambda isn't being returned properly and the display looks like an internal mupad format.

newLambda=sum(sum(L.*expandxmat,3),4)

newLambda = 

This doesn't work either

newLambda = sum(L.*expandxmat,[3 4])

newLambda = 

●

But this does; I don't know why the others don't.

newLambda = sum(sum(L.*expandxmat,4),3)

newLambda = 

N_equivalence=isequal(ogN,newN)

N_equivalence = logical
   1

Lambda_equivalence=isequal(ogLambda,newLambda)

Lambda_equivalence = logical
   1

size(ogLambda)

ans = 1×2

3 3

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

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

size(newLambda)

ans = 1×2

3 3

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

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

ogLambda-newLambda

ans = 

BTW, be careful using isequal. You might really want to use isAlways.

I suggest you file a bug report with Tech Support.

[Jacob]

Very strange, but thank you for the fix!