I have two matrix $A$ and $B$ and $fl(A), fl(B)$ denote the stored version them, respectively. Let $fl(AB)$ be the stored version of the product of $A$ and $B$. Is it true that $fl(AB) = fl(A) * fl(B)$?
Intuitively, I don't think this is the case since the right hand side must have more error, correct?
But I'm not sure how to show this.
There are usually two meanings of $\mathrm{fl}(\cdot)$. The first one is the representation of an object (matrix, vector) in the finite precision arithmetic. So for an $A\in\mathbb{R}^{m\times n}$, $\mathrm{fl}(A)\in\mathbb{R}^{m\times n}$ is what you call the "stored version" of $A$. But when you write $\mathrm{fl}(AB)$, it normally means the computed result of the product of $A$ and $B$, which can be quite different from $\mathrm{fl}(A)\mathrm{fl}(B)$.
Consider a matrix $$ A=B=\pmatrix{10^{-2} & 10^2 \\ 10^1 & 10^0} $$ and the 2-digit floating point arithmetic with base 10. We have $\mathrm{fl}(A)=\mathrm{fl}(B)=A=B$ (each mantisa contains only a single digit) but $$ \mathrm{fl}(AB)=\pmatrix{10^3 & 10^2 \\ 10^1 & 10^3}\neq \pmatrix{10^3 & 1.01\cdot 10^2\\1.01\cdot 10^1&1.001\cdot 10^3}=AB. $$
In general, usual implementations of the matrix product (it is not true, e.g., for the Strassen algorithm) give for $A\in\mathrm{R}^{m\times n}$ and $B\in\mathrm{R}^{n\times p}$ that $$ |\mathrm{fl}(AB)-AB|\leq\gamma_n|A||B|, \quad \gamma_n:=\frac{nu}{1-nu}, $$ where $u$ is the unit roundoff. See, e.g., Section 3.5 here.