I have asked this question at stackexchange for Computational Science, but have so far no answers (https://scicomp.stackexchange.com/questions/41041/generate-polynomial-basis-through-a-sequence-of-svd). I will instead try my luck here. I will remove the question here if it's considered a duplicate.
I need help to understand how to use the result given by an algorithm for constructing an orthonormal polynomial basis over $L^{2}(X)$, where $X\subset\mathbb{R}^2$, with respect to the inner product $(f,g) = \iint_X f(x,y)\overline{g(x,y)}\,dxdy$. We denote $\mathcal{P}_{n} = \{\text{Polynomials }P(x,y): \text{ total degree of } P \leq n\}\subset L^{2}(X)$. The algorithm is as follows, as described in the paper https://epubs.siam.org/doi/abs/10.1137/110860082 (paywall) or in a previous version of it in https://www.cs.yale.edu/publications/techreports/tr1443.pdf (free)
(1) Consider the vector corresponding to the constant polynomial $1$. Normalize it. The resulting vector provides an orthonormal basis of $\mathcal{P}_{0}$.
(2) Suppose we have constructed an orthonormal basis of $\mathcal{P}_{k}$, consisting of $K$ vectors. Multiply each of the vectors in this basis by $x$ and $y$, respectively, and put the resulting vectors, together with the vectors in the basis of $\mathcal{P}_{k}$, in a matrix having $3K$ columns. Calculate the SVD of the matrix. The SVD yields an orthonormal basis of the column space of the matrix, which is also an orthonormal basis of $\mathcal{P}_{k+1}$.
(3) Apply the procedure described in step (2) for each $k$ in $\{0, 1, \ldots , n − 1\}$ in order to obtain an orthonormal basis of $\mathcal{P}_n$.
My situation: I want to do this numerically, i.e. without symbolic software. We can assume that an SVD with respect to the given inner product is given.
Questions: (1) What does it mean to multiply with x and y in step (2)? The SVD is not performed symbolically. Should I sample $X$ and multiply by those values?
(2) Executing the algorithm above I get a basis expressed as a unitary matrix, how do I use the result to evaluate the resulting orthonormal polynomial basis at a point in $X$?