For $f(x)=|x|^2$, you get $\varepsilon\sim\tfrac C {k^{2/n}}$, where $k$ is the number of pieces fro for your PL-function. In particular, you will not get linear bound for $n\ge 3$.
I think for any convex 1-Lipschitz function you should get $\varepsilon=O(k^{-2/n})$.
An easy construction gives $\varepsilon=O(k^{-1/n})$, simply take the maximum of supporting linear functions with gradients $\{v_1,v_2,\dots,v_k\}$ which form a $C k^{-1/n}$ dense set in the unit ball.
(Note that the worse case for this approximation is a linear function with gradient sufficiently far from $v_i$.)
