0

I encountered some numerical questions when running simulation on MatLab. Here please find the questions:

  1. I found that A*A' (a matrix times its transpose) is not guaranteed to be symmetric in MatLab. Can I know what is the reason? And because I will have A*C*A', where C is a symmetric matrix, and I would like to keep A*C*A' as symmetric. Is there any method to fix the numerical difference created by the transpose operation?

  2. I implemented a for loop in Matlab to compute a set of matrices. Small numerical difference (around 10^(-10)) in each round accumulates to the next run, and it finally diverges after around 30 rounds. Is there any method to fix small error in each run and do not affect the result at the same time.

Thank you for reading my questions!

user2874626
  • 13
  • 4
  • 2
    As you noticed, SO doesn't do LaTeX. Please use code (enclosed in backpacks) to talk about computations. Also, please ask only one question in each post. The first question is explained by [this post](https://stackoverflow.com/q/686439/7328782), it is not the transposition that causes differences, but the matrix multiplication -- the order of operations is different for the upper and lower triangles, and therefore the rounding error is different too. – Cris Luengo Oct 21 '21 at 07:27
  • The other question I don't understand, you cannot fix the errors, if you knew what the error is, then you wouldn't be making an error. You need to develop algorithms that are robust against rounding errors. That is the cost of using numerical methods. – Cris Luengo Oct 21 '21 at 07:28
  • Thank you for the answers! I will apply code next time. – user2874626 Oct 21 '21 at 13:56
  • You can [edit] your question right now, no need to wait for next time. – Cris Luengo Oct 21 '21 at 14:09

1 Answers1

1

"I found that A*A' (a matrix times its transpose) is not guaranteed to be symmetric in MatLab."

I would dispute that statement as written. The MATLAB parser is smart enough to recognize that the operands of A*A' are the same and call a symmetric BLAS routine in the background to do the work, and then manually copy one triangle into the other resulting in an exactly symmetric result. Where one usually gets into trouble is by coding something that the parser cannot recognize. E.g.,

A = whatever;
B = whatever;
X = A + B;
(A+B) * (A+B)'  <-- MATLAB parser will call generic BLAS routine
X * X'  <-- MATLAB parser will call symmetric BLAS routine

In the first matrix multiply above, the MATLAB parser may not be not smart enough to recognize the symmetry so a generic matrix multiply BLAS routine (e.g., dgemm) could be called to do the work and the result is not guaranteed to be exactly symmetric. But in the second matrix multiply above the MATLAB parser does recognize the symmetry and calls a symmetric BLAS matrix multiply routine.

For the ACA' case, I don't know of any method to force MATLAB to generate an exact symmetric result. You could manually copy one resulting triangle into the other after the fact. I suppose you could also factor C into two parts X*X' and then regroup but that seems like too much work for what you are trying to do.

James Tursa
  • 2,242
  • 8
  • 9
  • Yes, thank you for the answer. I tried to simplify the scenario when described it. The details were after the first statement. The case I encountered is that 'A*C*A'' is not symmetric. I fixed the numerical error by adding the result with its transpose and dividing the matrix sum by 2. – user2874626 Nov 08 '21 at 21:33