Third-order density-functional perturbation theory: a practical   implementation with applications to anharmonic couplings in Si
  We present a formulation of third-order density-functional perturbation theory which is manifestly invariant with respect to unitary transfomations within the occupied-states manifold and is particularly suitable for a practical implementation of the so called `2n+1' theorem. Our implementation is demonstrated with the calculation of the third-order anharmonic coupling coefficients for some high-simmetry phonons in Silicon. 
Alberto Debernardi and Stefano Baroni