MATLAB: My matrix should be symmetric but isn’t by a ridiculous margin

matrix manipulationoptimization

Given the following matrices:
Phi = [...
    1,89393939393939e-06	0	0	0;
    3,78787878787879e-06	1,89393939393939e-06	0	0;
    5,68181818181818e-06	3,78787878787879e-06	1,89393939393939e-06	0;
    7,57575757575758e-06	5,68181818181818e-06	3,78787878787879e-06	1,89393939393939e-06;
    9,46969696969697e-06	7,57575757575758e-06	5,68181818181818e-06	3,78787878787879e-06;
    1,13636363636364e-05	9,46969696969697e-06	7,57575757575758e-06	5,68181818181818e-06;
    1,32575757575758e-05	1,13636363636364e-05	9,46969696969697e-06	7,57575757575758e-06;
    1,51515151515152e-05	1,32575757575758e-05	1,13636363636364e-05	9,46969696969697e-06;
    1,70454545454545e-05	1,51515151515152e-05	1,32575757575758e-05	1,13636363636364e-05;
    1,89393939393939e-05	1,70454545454545e-05	1,51515151515152e-05	1,32575757575758e-05];
lamdba_matrix = 1e-5 * eye(4);
alpha_matrix = 1e-8 * eye(10);
H = 2*(Phi'*alpha_matrix*Phi + lambda_matrix);
H is intended to be used with quadprog. Doing so, I get the following warning: 
Warning: Your Hessian is not symmetric. Resetting H=(H+H')/2.
Well, I got suspicious because it should be symmetric according to the formula. So I investigated:
K>> H == H'
ans =
  4×4 logical array
   1   0   1   0
   0   1   0   1
   1   0   1   1
   0   1   1   1
   
K>> H(4, 1) - H(1, 4)
ans =
   3.0815e-33
What is going on here ? Where does this ridiculously small error come from ? Should I do H = (H + H')/2 to suppress the warnings ? Would it change the result of quadprog ?

Best Answer

Welcome to the world of floating point arithmetic. This is typical and not unusual behavior.
https://www.mathworks.com/matlabcentral/answers/57444-faq-why-is-0-3-0-2-0-1-not-equal-to-zero
You cannot expect floating point arithmetic, particularly matrix multiplication, to maintain properties like this. You must have algorithms in place in your code to deal with it.
In your particular case, if you rearrange your logic a bit you can manage to get MATLAB to call a symmetric matrix multiply routine instead of the generic matrix multiply routine, thus preserving the symmetry:
H = 2*(Phi'*Phi*1e-8 + lamdba_matrix);
Here MATLAB sees the Phi'*Phi matrix multiply and sees that the operands are the same, so it knows the result is symmetric and can call a faster symmetric matrix multiply routine. This will run faster (about 1/2 the work) and preserve symmetry for that operation.
Related Question