spectral decomposition of a bivariate function Now I have a function $f=f(x,y)$, smooth and symmetric(i.e. $f(x,y)=f(y,x)$ everywhere), with arguments defined on a compact set: $(x,y)\in[0,1]\times[0,1]$.
I'd wish to know if $f$ can be expanded as $f(x, y) = \sum_{l=1}^\infty \lambda_l \xi_l(x)\xi_l(y)$, where $\xi_l(\cdot)$ are functions that are norm 1 and orthogonal to each other -- here I define the needed inner product between two functions $u, v$ on $[0,1]$ to be $<u, v> := \int_0^1 u(t)v(t)dt$.
A reasonably easy to follow reference will be good enough for my question. Just I didn't know what the correct key words to use -- I tried "functional spectral theory" etc and got something too abstract for me to understand(and they were unnecessarily general for my purpose of application). Thank you so much in advance!
 A: For simplicity, assume $f(x,y)=f(y,x)$ is a real function. Because $f$ is in $L^{2}([0,1]\times[0,1])$, then the integral operator $K$ given by
$$
          Kg = \int_{0}^{1}f(x,y)g(y)\,dy
$$
is a selfadjoint Hilbert-Schmidt integral operator on $L^{2}[0,1]$. So there is an orthonormal basis $\{\varphi_{n}\}_{n=1}^{\infty}$ consisting of real eigenfunctions of $K$ with corresponding real eigenvalues of $K$. The eigenfunctions with non-zero eigenvalues are smooth because they inherit smoothness from $k$:
$$
             \varphi_{n}(x)=\frac{1}{\lambda_{n}}\int_{0}^{1}f(x,y)\varphi_{n}(y)\,dy.
$$
To show that $\{ \varphi_{n}(x)\varphi_{m}(y)\}_{n,m=1}^{\infty,\infty}$ is an orthonormal basis of $L^{2}([0,1]\times[0,1])$, notice that this set is an orthonormal subset of $L^{2}([0,1]\times[0,1])$, and Parseval's equality for $L^{2}[0,1]$ gives Parseval's equality for every $g \in L^{2}([0,1]\times[0,1])$:
$$
\begin{align}
    &\sum_{n,m}\left|\int_{0}^{1}\int_{0}^{1}g(x,y)\varphi_{n}(x)\varphi_{n}(y)\,dx\,dy\right|^{2} \\
   & = \sum_{n}\sum_{m}\left|\int_{0}^{1}\left(\int_{0}^{1}g(x,y)\varphi_{n}(x)\,dx\right)\varphi_{m}(y)\,dy\right|^{2} \\
   & = \int_{0}^{1}\sum_{n}\left|\int_{0}^{1}g(x,y)\varphi_{n}(x)\,dx\right|^{2}\,dy \\
   & = \int_{0}^{1}\int_{0}^{1}|g(x,y)|^{2}\,dx\,dy.
\end{align}
$$
So, the following sum converges in $L^{2}([0,1]\times[0,1])$:
$$
\begin{align}
        K(x,y) = \sum_{n,m}\left[\int_{0}^{1}\int_{0}^{1}f(x,y)\varphi_{n}(y)\varphi_{m}(x)\,dydx\right]\varphi_{n}(x)\varphi_{m}(y)
\end{align}
$$
Because $\varphi_{n}$ is an eigenfunction of $K$ with eigenvalue $\lambda_{n}$, then
$$
\begin{align}
  K(x,y)
    & = \sum_{n,m}\left[\int_{0}^{1}\left[\int_{0}^{1}k(x,y)\varphi_{n}(y)\,dy\right]\varphi_{m}(x)\,dx\right]\varphi_{n}(x)\varphi_{m}(y) \\
    & = \sum_{n,m}\left[\lambda_{n}\int_{0}^{1}\varphi_{n}(x)\varphi_{m}(x)\,dx\right]\varphi_{n}(x)\varphi_{m}(y) \\
    & = \sum_{n}\lambda_{n}\varphi_{n}(x)\varphi_{n}(y).
\end{align}
$$
All functions $\varphi_{n}$ in the final sum are smooth because $0$ eigenvalue terms drop out. The convergence of the sum is in $L^{2}([0,1]\times[0,1])$.
