Discretization of Hamiltonian using finite difference always justified? I have this continuum version $$ H_{R}=\int dx\psi^{\dagger}(x)(\frac{p^{2}}{2}+V)\psi(x) $$ with $V$ as constant potential.
Is it always justified to go from this to $$ \sum_{i}c_{i}^{ \dagger }\left[c_{i+1}+c_{i-1}-2c_{i}\right] +V \sum_{i}c_{i}^{ \dagger }c_{i} $$
using the finite difference form of the one-dimensional second derivative? 
Ignore the factor $-1/2 $ and assume lattice constant $a= 1 $ and $\hbar =1 $. 
Actually I am thinking whether it is justified to use this even when the derivative of eigenfunction is discontinuous at some of the points in real space like for the delta-function barrier. Will that affect the second derivative of the field operators ? 
 A: [I answered this on the assumption you were talking about approximating a solid-state type system, so the answer is pretty particular that. The issue of actually defining a continuum limit of a quantum field theory is non-trivial.]
The usual way to go the continuum limit is to go to Fourier space. Since you are on a lattice you will end up in a Brillouin Zone. Taking the continuum limit is now equivalent to making the Brillouin zone infinite, instead of periodic. You can then take a regular Fourier transform to get a continuous real space theory. Going from continuum to lattice is the inverse process.
This, I think, makes it clear when you can make this transformation. You can do it if
1) You are concerned with wavevectors in a region of the the BZ which doesn't know about the size of the BZ. This includes not just the naive guess of small wave vectors, but also say fermions close a large-ish Fermi-surface. And of course all external perturbation need to couple to this small region.
2) You are not concerned about any processes that know about the periodicity of the lattice. So no Umklapp scattering, no Fermi-surface nesting, etc... This can both be very important. For example, your lattice model on the square lattice at half-filling has perfect nesting vectors.  These lead to antiferromagnetic states that would have no equivalent in the continuum model (I think).
Another way of thinking about pt. (2) is that your continuum and lattice models have completely different symmetry groups. Your continuum model has a continuous translations $\mathbb{R}^n$ whereas the lattice model has symmetries $\mathbb{Z}^n$. Any time you are changing the symmetries you should be cautious. 
They also have different rotational symmetries. I see in your particular model, you have also made a particular expansion of the the dispersion relation, by taking $E = p^2/2m - \mu$ in the continuum. This gives you a spherical Fermi-surface, whereas the lattice model may have an extremely aspherical shape (like the half-filling model). In a lot of life this is actually fine, and leads only to numerical discrepancies, but in extreme case you should be aware. You can always change the dispersion to be $E = v_F(p-p_F)$ where $v_F$ and $p_F$ depend on the angle. Also be careful of accidentally maintaing or destroying rotational symmetries that don't exist in the lattice model. These usually correspond to irrelevant operators, but again be aware, be cautious.
