Compute $\lim\limits_{n \to \infty} n\sum\limits_{k=1}^n(f(k/n) - f((k-1)/n))\int_{(k-1)/n}^{k/n}f(t)dt$ 
Let $f : [0,1] \to \mathbb{R}$ be a continuous function and let $(a_n)_{n>0}$ and $(b_n)_{n>0}$ be two sequences such that $$\displaystyle{ a_n = \sum_{k=1}^n{f \left(\frac{k-1}{n}\right) \cdot \int_\frac{k-1}{n}^\frac{k}{n}}{f(t)dt}},$$
$$\displaystyle  b_n = \sum_{k=1}^n{f \left(\frac{k}{n}\right) \cdot \int_\frac{k-1}{n}^\frac{k}{n}}{f(t)dt}, $$ $\forall n \in \mathbb{N}^*$.
a) Prove that $\displaystyle{\lim_{n \to \infty}{(b_n - a_n)} = 0}$.
b) Compute $\displaystyle{\lim_{n \to \infty}{n(b_n - a_n)}}$.

I have managed to solve a).
Proof for a) :
From the mean value theorem, we know that $\displaystyle{\exists c_k \in \left(\frac{k-1}{n}, \frac{k}{n}\right)}$ such that $\displaystyle{\int_{\frac{k-1}{n}}^{\frac{k}{n}}{f(t)dt} = \frac{1}{n}f(c_k)}$
So, $b_n = \displaystyle{\frac{1}{n}\sum_{k=1}^n {f(c_k) \cdot f\left(\frac{k}{n}\right)} = \frac{1}{2} \cdot \frac{1}{n}\sum_{k=1}^n{\left(\frac{f(c_k) + f\left(\frac{k}{n}\right)}{2}\right)^2 \cdot 4 - (f(c_k))^2 - \left(f\left(\frac{k}{n}\right)\right)^2}}$.
Since $f$ has the intermediate value property, then $\exists x_k \in \displaystyle{\left(c_k, \frac{k}{n}\right) \subset \left(\frac{k-1}{n}, \frac{k}{n} \right)}$ such that $\displaystyle{\frac{f(c_k) + f\left(\frac{k}{n}\right)}{2} = f(x_k)}$.
Therefore, $\displaystyle{\lim_{n \to \infty}b_n = \int_0^1{(f(t))^2dt}}$.
Using the same method, $\displaystyle{\lim_{n \to \infty}a_n = \int_0^1{(f(t))^2dt}}$, so $\displaystyle{\lim_{n \to \infty}{(b_n - a_n)} = 0}$.
I have trouble solving b). I tried using the same method but it doesn't work, not unless $f$ is differentiable (so I can use Lagrange's theorem).
$\displaystyle{n(b_n - a_n) = \sum_{k=1}^n{f(c_k) \left( f\left(\frac{k}{n}\right) - f\left(\frac{k-1}{n} \right) \right)}}$ and if $f$ is differentiable, then $\displaystyle{f\left(\frac{k}{n}\right) - f\left(\frac{k-1}{n} \right) = \frac{1}{n}f'(c_k)}$.
 A: I think this is less easy than it sounds. So far, I only have a solution when $f$ is ${\mathcal C}^2$. The density argument suggested in the comments unfortunately seems to break down along the way (see Martin R's comment below). 
We will show that the limit of $d_n=n(b_n-a_n)$ is $\frac{f(1)^2-f(0)^2}{2}$. With my additional hypothesis on $f$, $|f|+|f''|$ is bounded by some constant $M$ on $[0,1]$. For any subinterval $(a,b) \subseteq [0,1]$, we have 
the following error estimate for the trapezoidal rule (see for example,
here for a proof) 
$$
\Bigg| \int_a^b f(t)dt - (b-a)\bigg(\frac{f(a)+f(b)}{2}\bigg) \Bigg| \leq \frac{M(b-a)^3}{12} \tag{1}
$$
Using (1) with $a=\frac{k-1}{n}$ and $b=\frac{k}{n}$, we deduce :
$$
\Bigg| \int_{\frac{k-1}{n}}^{\frac{k}{n}} f(t)dt - \frac{1}{n}\bigg(\frac{f(\frac{k-1}{n})+f(\frac{k}{n})}{2}\bigg) \Bigg| \leq \frac{M}{12n^3} \tag{2}
$$
Multiplying by $f(\frac{k}{n})-f(\frac{k-1}{n})$ :
$$
\Bigg| \bigg(f(\frac{k}{n})-f(\frac{k-1}{n}\bigg)\int_{\frac{k-1}{n}}^{\frac{k}{n}} f(t)dt - \frac{1}{n}\bigg(\frac{f(\frac{k}{n})^2-f(\frac{k-1}{n})^2}{2}\bigg) \Bigg| \leq 
\frac{M}{12n^3}(|f(\frac{k}{n})|+|f(\frac{k-1}{n})|) \leq
\frac{M^2}{6n^3} \tag{3}
$$
and summing and multiplying by $n$, we deduce
$$
\Big|d_n -\frac{f(1)^2-f(0)^2}{2}\Big| \leq \frac{M^2}{6n^2}
$$
which finishes the proof.
A: It seems, if I didn't made any mistake, in general $n(b_n-a_n)$ may not converge!  I'll construct an explict counterexample to it. Consider periodic function $g:\mathbb{R}\to \mathbb{R}$ defined as:
$$g(x) = \begin{cases}0, \{x\} \le \frac{1}{4}\\ 4x - 1, \frac{1}{4} \le \{x\} \le \frac{1}{2}\\ 1, \frac{1}{2} \le \{x\} \le \frac{3}{4}\\ 4 - 4x, \frac{3}{4} \le \{x\} \le 1\end{cases}$$
and put $g_n(x) = g(2^{9n}x)$. We will construct our function $f$ as series $f(x) = \sum\limits_{n = 1}^\infty h_ng_n(x)$ with $\sum\limits_{n = 1}^\infty |h_n| < \infty$. We will prove that $2*2^{9n}(b_{2*2^{9n}} - a_{2*2^{9n}})$ is unbounded for suitable choise of $h_n$ thus sequence from the post can not converge. We are interested in the following sum:
\begin{equation}
\sum\limits_{k = 1}^{2*2^{9n}} \left(f\left(\frac{k}{2*2^{9n}}\right) - f\left(\frac{k-1}{2*2^{9n}}\right)\right)\int_{\frac{k-1}{2*2^{9n}}}^{\frac{k}{2*2^{9n}}}f(x)dx.
\end{equation}
Since for fixed $n$ everything is absolutely convergent we are basically interested in computing
\begin{equation}
\sum\limits_{k = 1}^{2*2^{9n}} \left(g_r\left(\frac{k}{2*2^{9n}}\right) - g_r\left(\frac{k-1}{2*2^{9n}}\right)\right)\int_{\frac{k-1}{2*2^{9n}}}^{\frac{k}{2*2^{9n}}}g_m(x)dx\qquad (1)
\end{equation}
for any $n, r, m > 0$. If $r > n$ then $ \left(g_r\left(\frac{k}{2*2^{9n}}\right) - g_r\left(\frac{k-1}{2*2^{9n}}\right)\right) = 0$ and so we are not interested in this cases. If $r < n$ then $ \left(g_r\left(\frac{k}{2*2^{9n}}\right) - g_r\left(\frac{k-1}{2*2^{9n}}\right)\right)$ has intervals of $k$ on which it is constant(those intervals of the length $2^{-9r}$ up to some factor of $2$ or $4$) and up to constant on them it equals to $0, 1, 0, -1, 0, 1, 0, -1, \ldots$ Now we have three cases:
1) $m > r$. In this case every copy of $g$ contained in $f_m$ contains in some interval of $k$ and moreover every interval is divided into some number(fixed for fixed $m, r$) of copyes of $g$ so we have sum of the kind $0 + I + 0 - I + ... = 0$ since in every tuple of $4$ intervals we get $0$(from now on lets call those tuples "quads").
2) $m < r$. In this case every quad is contained in some copy of $g$ in $f_m$ and moreover this copy(or, more preciesly, interval, to which it correspondes)  is divided into some number of this quads. Even more since $9$ is a huge number we can say that every quater of copy of $g$ is divided into some number of our quads. On the first and third quaters of $g$ we clearly has zero sum. With some computations one can show that second and fourth quaters also adds up to $0$.
3) $m = r$. In this case we up to some constant interested in $0*\int_0^{\frac{1}{4}}g(x)dx + 1*\int_\frac{1}{4}^\frac{1}{2}g(x)dx + 0*\int_\frac{1}{2}^\frac{3}{4}g(x)dx + (-1)*\int_\frac{3}{4}^1g(x)dx = 0$.
So the remainig case is $r = n$. In this case we would have $ \left(g_r\left(\frac{k}{2*2^{9n}}\right) - g_r\left(\frac{k-1}{2*2^{9n}}\right)\right) = (-1)^{k-1}$. First and second cases could be treated the same as before while in the third case we would have:
$$1*\int_0^\frac{1}{2}g(x)dx + (-1)*\int_\frac{1}{2}^1g(x)dx = \frac{-1}{4}.$$
So (1) is nonzero iff $n = r = m$ and equals to some fixed constant $c < 0$. Thus we have:
$$2*2^{9n}(b_{2*2^{9n}} - a_{2*2^{9n}}) = 2ch_n^2*2^{9n}.$$
Choosing $h_n = 2^{-n}$ we get what we want. Moreover since over some subsequence $h_n$ can decay arbitary slow we can not say anything better than $b_n - a_n \to 0$
In fact, when I was thiking about this problem I treated (1) as some kind of bilinear form on the space of continious functions so I wanted to find sequence of functions $(g_m)$ such that $<g_r, g_m>_n = c\delta_{n, m}\delta_{n, r}$ so I was searching for some kind of continious Rademacher functions.
A: It seems we can compute limit in the case when $f$ is of bounded variation which means that even if counterexample exists it must look very ugly. $f$ is continious so it is uniformly continious so we can assume that  $|\int_{\frac{k-1}{n}}^{\frac{k}{n}} f(t)dt - \frac{1}{2n}(f(\frac{k}{n}) + f(\frac{k-1}{n}))| < \frac{\epsilon}{n}$ with $\epsilon \to 0$ as $n\to \infty$. So we have:
$|n(b_n - a_n) - \frac{1}{2}\sum\limits_{k=1}^n (f(\frac{k}{n}) - f(\frac{k-1}{n}))(f(\frac{k}{n}) + f(\frac{k-1}{n}))|\le \epsilon\sum\limits_{k=1}^n |f(\frac{k}{n}) - f(\frac{k-1}{n})| \le \epsilon Var(f)$
On the other hand $\frac{1}{2}\sum\limits_{k=1}^n (f(\frac{k}{n}) - f(\frac{k-1}{n}))(f(\frac{k}{n}) + f(\frac{k-1}{n})) = \frac{1}{2}(f(1)^2 - f(0)^2)$ so the limit is equal to $\frac{1}{2}(f(1)^2 - f(0)^2)$.
