Minimal surfaces in $mathbb{R}^n$ can be locally approximated by graphs of harmonic functions, i.e., functions that are critical points of the Dirichlet energy, but no analogous theorem is known for $H$-minimal surfaces in the three-dimensional Heisenberg group $mathbb{H}$, which are known to have singularities. In this paper, we introduce a definition of intrinsic Dirichlet energy for surfaces in $mathbb{H}$ and study the critical points of this energy, which we call contact harmonic graphs. Nearly flat regions of $H$-minimal surfaces can often be approximated by such graphs. We give a calibration condition for an intrinsic Lipschitz graph to be energy-minimizing, construct energy-minimizing graphs with a variety of singularities, and prove a first variation formula for the energy of intrinsic Lipschitz graphs and piecewise smooth intrinsic graphs.