  For the moment, there is no exact description of van der Waals (vdW)
interactions. ACFD-RPA \cite{Gould1} is expected to better describe vdW
bonding, but it is not exact. The PBE/DFT-D2 method is less satisfactory,
however, its results are in good agreement with experimental data. Although our
fitting technique may weaken (not neglect) the vdW interactions and produce
interlayer potentials with weakened vdW, the obtained interlayer potentials
reproduce energetics of graphite near the equilibrium interlayer distance very
well, as shown in Ref. \cite{Chen}. If having inputs which fully include vdW
interactions and having better fitting functions, we believe that interlayer
potentials can also fully include vdW interactions in graphite system.
