Why does MATLAB's built-in "chol" yield a large error compared to my written algorithm

I wish to know why MATLAB's built-in function

chol

yields a large error on a matrix

S=randi([1,100],[25,25]);
A=S'*S;
[B]=chol(A)

Executing the above code and calculating the norm produced a massive error :

norm(A-B*B')
>>> 1.4697e+06

How could this be justified knowing that I wrote a modified Cholesky decomposition algorithm that producd an error of magnitude :

9.4093e-11

which is way less than MATLAB's built-in "chol". The reason I am asking is because I thought MATLAB have their matrix operations algorithms optimized using LAPAK which is written in fortan.

Answer 1

[Chunru]

S=randi([1,100],[25,25]);
A=S'*S;
B=chol(A)

B = 25×25

288.6191 181.7724 223.0310 217.7576 189.2702 192.6899 233.5188 254.8757 169.3062 217.2621 256.0676 176.3951 217.1374 222.3657 213.8285 226.2532 192.8077 161.1501 224.4411 202.2978 224.9574 205.7036 223.2388 251.2793 253.9125 0 204.5746 120.0448 103.2272 90.0400 82.0888 11.3881 90.7914 70.5014 97.8750 88.3051 92.2364 54.4418 79.3600 118.8069 67.7181 142.8323 77.1765 109.6460 96.3748 122.9378 102.5677 125.8188 81.1555 88.4944 0 0 174.0472 62.6923 28.3662 49.8936 12.6804 4.9974 75.0782 17.1818 -3.8815 95.6805 48.0570 54.1623 47.6852 27.7242 31.2653 8.8936 48.7813 93.1416 25.3925 81.0642 51.3075 10.0987 13.9983 0 0 0 145.9022 -5.3799 49.5644 -34.8317 -51.6618 28.2230 -1.4599 7.0883 30.2866 29.5747 12.0595 9.0377 18.0409 1.9285 -21.4180 11.6213 -22.8553 13.7036 44.0815 -22.5662 34.3961 -46.1678 0 0 0 0 128.4874 72.2299 28.1380 46.1919 -7.8950 48.8909 50.8244 31.4704 9.1086 -24.1135 66.5609 49.2757 44.6668 40.6703 -12.2016 -3.6623 28.8152 -6.6083 45.1096 11.0880 31.1816 0 0 0 0 0 138.1154 20.0179 51.1375 73.2684 79.1878 44.3817 45.9310 25.4732 21.3693 1.3095 1.2726 6.2278 43.4123 1.2250 42.7116 -3.3906 -2.8876 -9.0278 29.1820 26.1871 0 0 0 0 0 0 113.5595 38.3018 6.9518 12.4341 25.6286 22.5409 26.2495 28.0526 12.9893 57.0123 4.2127 72.7048 44.5154 45.7824 56.6803 12.6413 45.1525 3.0047 40.2931 0 0 0 0 0 0 0 146.7843 -30.7489 -9.4929 48.5814 4.7251 49.9103 -34.9212 -4.6244 9.3329 5.0302 0.3031 -12.0241 47.3106 31.7036 -46.4708 -21.7946 34.8589 -61.8425 0 0 0 0 0 0 0 0 114.3953 -9.8188 -35.7179 11.8700 11.9243 -1.5256 1.1218 31.9549 43.0921 33.9370 24.7113 -0.1008 68.0062 38.8325 28.7445 -5.9466 33.6106 0 0 0 0 0 0 0 0 0 84.8655 -59.5429 -37.2827 -75.1867 -18.7067 32.1745 19.9537 45.7738 -46.2919 -10.1079 -41.5523 33.1476 3.6979 -23.8847 5.1887 -4.0784

For Cholesky factorization, we have A = B'*B (not B*B'), so we should compute:

norm(A-B'*B)

ans = 7.8958e-11

norm(A-B'*B, 'fro')

ans = 1.7493e-10

[Tarek Hajj Shehadi]

Thank you very much for your detailed answer!