MATLAB: Single precision matrix multiplication

hermitiansingle precision

Script:
M = 10000;  N = 100;
A = round(10000 * complex(randn(M, N), randn(M, N)));
B = single(A') * single(A);
C = B' - B;
max(abs(C(:)))  %// ==> answer is non-zero
A = single(A);
B = A' * A;
C = B' - B;
max(abs(C(:)))  %// ==> answer is zero
Not sure what caused the difference? Thanks!

Best Answer

The Matlab interpreter runs a special dedicated "ctranspose-multiply" function to evaluate expressions of the form,
P'*Q
It does not do any explicit transposition of P, thus saving time and memory. This special function likely does a pre-check for the special case when P and Q are the same variable. In this case, it knows that it only has to compute the lower (or upper) triangular part of the resulting matrix, and can just copy the conjugate of that to the upper(lower) triangle. Thus, the output of A'*A will always be perfectly Hermitian.
However, for the expression,
single(P')*single(Q)
no special functions are used. Every operation in this expression (single,ctranpose,mtimes) is done separately and explicitly, and can generate floating point errors that may not be Hermitian across the final matrix.
Related Question