Scalar product and uniform convergence of polynomials Given two functions $u$ and $v$, you can compute $(u|v) = \int_0^1 (u(t)'-u(t))(v(t)'-v(t)) \,\mathrm{d}t$. It resembles the scalar product $\int_{-1}^1 u(t)v(t) \,\mathrm{d}t$, which leads to Legendre polynomials, with a major difference: while the latter is a scalar product over continuous functions, the former is not, since for $f(t)=e^t$, $(f|f)=0$.
However, it's obviously a scalar product over polynomials (the solutions of $(f|f)=0$ are not polynomials, and all other defining properties of a scalar product are trivial).
We can therefore find an orthonormal basis of polynomials, using Gram-Schmidt orthogonalization on the basis $[1,t,t^2,...]$. Let's call $P_n$ this orthonormal basis.
Now, after doing this numerically with a CAS, and after some plots, it looks like $\frac{P_n}{P_n(0)}$ tends uniformly to the function $t \to e^t$ (on $[0,1]$), a solution of $y'-y=0$, which is closely related to the definition of the scalar product above.
But I don't know how to prove this. Any idea?

Also, if I do the same with $(u|v)=\int_0^1 (u''(t)-u(t))(v''(t)-v(t)) \,\mathrm{d}t$, I'm not absolutely sure, but numerically it looks like, up to a constant factor $P_{2n} \to \cos t$ and $P_{2n+1} \to \sin t$, uniformly on $[0,1]$, so the orthogonal basis of polynomials leads to a basis of solutions of $y''-y=0$.
I'd like to know if the process can be generalized to other linear differential equations (provided they don't have polynomial solutions). Does it work with equations with nonconstant coefficients?
This is not homework. It's something I discovered while playing with a Ti89 long ago (initially, to test a program doing Gram-Schmidt orthogonalization), and I have never been able to prove or disprove the numerical evidence.

$\color{red}{\textbf{Now, why this question?}}$ Apart from mathematical curiosity, there may be some motivation: the method seems quite robust, there is convergence, probably under certain conditions, for several weight functions and differential equations. And we get uniform approximation of solutions of linear differential equations, by polynomials, which may be quite useful in practice. Since many linear ODE give rise to special functions, that would give a cheap way to compute uniform approximations of their solutions, on arbitrary (but bounded) intervals.

To help visualize, here are some sample graphics done with Maxima, computed with a straightforward Gram-Schmidt algorithm.
The two following are $P_4(x)/P_4(0)-\exp(x)$ and $P_9(x)/P_9(0)-\exp(x)$, in the case $(u|v) = \int_0^1 (u(t)'-u(t))(v(t)'-v(t)) \,\mathrm{d}t$. The division by $P_n(0)$ is only done in order to "normalize" polynomials, since $\exp 0=1$.


Now example with "trigonometric" differential equation $y''+y=0$. And $(u|v) = \int_0^{2\pi} (u(t)''+u(t))(v(t)''+v(t)) \,\mathrm{d}t$. Notice the integration bound differ from the preceding, and now uniform convergence of $P_n$ takes place on $[0,2\pi]$. The normalization is here $P_{2n}(x)/P_{2n}(0)$ as $\cos 0 = 1$ ans $P_{2n+1}(x)/P_{2n+1}(\pi/2)$, as $\sin \pi/2=1$, and precisely, $P_{8}(x)/P_{8}(0)-\cos x$ and $P_{9}(x)/P_{9}(\pi/2)-\sin x$ are respectively plotted.



Notice, it's also possible to add a weighting function to our scalar product, e.g. $(u|v) = \int_0^1 (u(t)'-u(t))(v(t)'-v(t)) \mathcal{W}(t)\,\mathrm{d}t$. For example, with$\mathcal{W}(t)=\frac{1}{\sqrt{t(1-t)}}$ on $[0,1]$, there is still convergence, and for example $P_9(x)/P_9(0) -\exp x$ yields the following plot.

 A: We prove the conjecture using the following more general formula.
Suppose $Y_n$ are orthogonal polynomials with respect to the
inner product on $[t_0,t_1]$ weighted by $w$:
$$
(f,g) = \int_{t=t_0}^{t_1} f(t) \, g(t) \, w(t)
$$
Then (up to a constant factor)
$$
\Upsilon_n(t) = e^t \int_t^\infty e^{-x} Y_n(x) dx
$$
are polynomials of degree $n$ orthogonal with respect to the inner product
$$
(f|g) = \int_{t=t_0}^{t_1} (f'(t)-f(t)) \, (g'(t)-g(t)) \, w(t).
$$
Note that the formula for $\Upsilon_n$ has the equivalent form
$$
\Upsilon_n(t) = e^t \Bigl( C_n - \int_0^t e^{-x} Y_n(x) dx \Bigr)
$$
where $C_n = \Upsilon_n(0) = \int_0^\infty e^{-x} Y_n(x) dx$.
In our setting, $t_0=0$, $t_1=1$, $w(t)=1$, and the $Y_n$ are
shifted
Legendre polynomials.
It is known that $Y_n$ has leading coefficient $2n\choose n$
and satisfies $|Y_n(t)| \leq 1$ for $0 \leq t \leq 1$.
Thus $\left|\int_0^t e^{-x} Y_n(x) dx \right| < 1$ for $0 \leq t \leq 1$.
But $C_n$ grows quickly with $n$: by general properties of orthogonal
polynomials, all the roots of $Y_n$ are real and in $(0,1)$;
therefore $Y_n(x) > {2n\choose n}(x-1)^n$ for $x>1$, whence
$\int_1^\infty e^{-x} Y_n(x) dx > {2n\choose n} n!/e$,
and $C_n$ is within $1$ of this.  Therefore $C_n \rightarrow\infty$
as $n \rightarrow \infty$ (numerically $C_4 = 1001$ and
$C_9 = 10622799089$; this agrees with OEIS sequence
A002119 up to a factor $(-1)^n$).  Hence
$$
\sup_{0 \leq t \leq 1} \left| \frac{\Upsilon_n(t)}{\Upsilon_n(0)} - e^t \right|
< 1/C_n
$$
decays rapidly to zero.
To prove the formula for $\Upsilon_n$, we observe that
$(u, \Upsilon'_n - \Upsilon_n) = 0$ for every $u$ of the form $f'-f$
with $f \in {\cal P}_n$, the vector space of polynomials of degree $<n$.
But the space of such $f'-f$ is just ${\cal P}_n$ itself (because the
linear map $f \mapsto f' - f$ of ${\cal P}_n$ has triangular matrix
with nonzero diagonal entries; or alternatively because it has kernel zero).
Therefore $\Upsilon'_n - \Upsilon_n$ is a scalar multiple of $Y_n$.
The integral formula for $\Upsilon_n$ is then obtained by solving
the differential equation $\Upsilon'_n - \Upsilon_n = -Y_n$:
multiply by the integrating factor $e^{-t}$ to get
$\frac{d}{dt}(e^{-t}\Upsilon_n(t)) = -e^{-t} Y_n(t)$, and note that
$e^{-t}\Upsilon_n(t) \rightarrow 0$ as $t \rightarrow \infty$
because $\Upsilon_n$ is a polynomial.
