  We study the band structure of self-adjoint elliptic operators $\mathbb{A}_g=
-\nabla \cdot \sigma_{g} \nabla$, where $\sigma_g$ has the symmetries of a
honeycomb tiling of $\mathbb{R}^2$. We focus on the case where $\sigma_{g}$ is
a real-valued scalar: $\sigma_{g}=1$ within identical, disjoint "inclusions",
centered at vertices of a honeycomb lattice, and $\sigma_{g}=g \gg1 $ (high
contrast) in the complement of the inclusion set (bulk). Such operators govern,
e.g. transverse electric (TE) modes in photonic crystal media consisting of
high dielectric constant inclusions (semi-conductor pillars) within a
homogeneous lower contrast bulk (air), a configuration used in many physical
studies. Our approach, which is based on monotonicity properties of the
associated energy form, extends to a class of high contrast elliptic operators
that model heterogeneous and anisotropic honeycomb media.
  Our results concern the global behavior of dispersion surfaces, and the
existence of conical crossings (Dirac points) occurring in the lowest two
energy bands as well as in bands arbitrarily high in the spectrum. Dirac points
are the source of important phenomena in fundamental and applied physics, e.g.
graphene and its artificial analogues, and topological insulators. The key
hypotheses are the non-vanishing of the Dirac (Fermi) velocity $v_D(g)$,
verified numerically, and a spectral isolation condition, verified analytically
in many configurations. Asymptotic expansions, to any order in $g^{-1}$, of
Dirac point eigenpairs and $v_D(g)$ are derived with error bounds.
  Our study illuminates differences between the high contrast behavior of
$\mathbb{A}_g$ and the corresponding strong binding regime for Schroedinger
operators.
