I have a linear function $u$ from approximately $\mathbb{R}^{10^6} \to \mathbb{R}^{3\times10^5}$. It has a certain computing cost so I would like to approximate it using a factored matrix représentation. This could be by using the first $k$ component of its SVD (Singular Value Decomposition), or by using a Non Negative Matrix Factorization (NMF) of small rank $k$.
Since the input and output dimension are really huge, it is not possible to actually compute or store the whole matrix $\mathbf{A}$ that would be associated to this linear function.
So my question is: is there a way to actually compute a factored representation (which does not need terabytes of storage) of the matrix $\mathbf{A}$ without ever storing this huge matrix, but by relying only on several evaluations of $u$?
Thank you!