In phonon theory, what is the physical significance of the force constants tensors? I want to compute phase equilibria of crystalline solids using First-principle methods (DFT in particular). The methodology for this computation is described in this chapter.
In order to calculate the vibrational energy of the lattice one may use the harmonic approximation:

I dont understand where that tensor came from or what is its significance.
Can someone explain?
 A: The atoms in each cell interact with the atoms in other cells. That matrix is assigns a Hessian matrix to every pair of atomic interactions. The Hessian matrix can be thought of as the 2nd order analogue to the Jacobian matrix.
As to why it's called the force constant matrix, and why it is written in that way, refer to this question.
Basically, for small displacements we can treat each atom as a kind of harmonic oscillator, and the second derivative gives us the spring constant $k$. Intuitively, this tells you how 'squished' each inter-atomic 'spring' is in the steady-state situation, which says how much potential energy is stored in the 'springs'.
A: That tensor mathematically arises from Taylor expanding the crystal potential energy $U(\mathbf{R})$ about equilibrium lattice positions $\mathbf{R}$ (which includes positions of all the atoms). We want to do this expansion in all three Cartesian directions, so use the 3D Taylor expansion:
$U(\mathbf{R}+\mathbf{u})= U(\mathbf{R})+ \mathbf{u}\cdot\nabla U(\mathbf{R}) + \frac{1}{2}(\mathbf{u}\cdot\nabla)^2U(\mathbf{R})+\frac{1}{3!}(\mathbf{u}\cdot\nabla)^3 U(\mathbf{R}) + ...$
where $\mathbf{u} = \mathbf{(R+u)} - \mathbf{R}$ is a small displacement about equilibrium. Since $U(\mathbf{R})$ is an arbitrary energy scale, we can ignore it and $\mathbf{u}\cdot\nabla U(\mathbf{R})=0$ if we are expanding about equilibrium positions (i.e., at the potential function minimum). Ignoring these terms, people usually write this in summation notation:
$U(\mathbf{R}) = \sum_{ij}\sum_{\alpha\beta}\Phi_{ij}^{\alpha\beta}u_i^{\alpha}u_j^{\beta}+\sum_{ijk}\sum_{\alpha\beta\gamma}\Psi_{ijk}^{\alpha\beta\gamma}u_i^{\alpha}u_j^{\beta}u_k^{\gamma} + ...$
where we sum over atom indices $ijk$ and Cartesian directions $\alpha\beta\gamma$, so that $u_i^{\alpha}$ is the displacement of atom $i$ in the $\alpha$ direction.
The $\Phi_{ij}^{\alpha\beta}$ (2nd order or harmonic force constants) is a $3\times3$ tensor for every $ij$ atom pair, representing 2nd order derivatives in every Cartesian direction combination $\alpha\beta$. Similarly, $\Psi_{ijk}^{\alpha\beta\gamma}$ (third order force constants) is a $3\times3\times3$ tensor for every $ijk$ triplet, representing 3rd order derivatives in every Cartesian direction combination $\alpha\beta\gamma$.
For a more qualitative understanding, let's extend this to a 1D harmonic mass-spring system with equilibrium position $x^0$ and general position $x$.
The potential energy is $U=\frac{1}{2}ku^2$, where $u=x-x^0$ is the displacement from equilibrium. The force is $F=-\frac{dU}{dx}=-\frac{dU}{du}=-ku$, which is zero for the equilibrium position. 
The 2nd derivative is simply $k$, which is the spring harmonic "force constant" or stiffness. It has the same physical meaning as $\Phi_{ij}^{\alpha\beta}$ in the 3D case, except masses can move in all 3 directions there. 
The harmonic force constants give rise to phonon dispersion and density of states. Anharmonic force constants (mainly 3rd order) give rise to thermal conductivity, thermal expansion, etc.
