The sphere tracing algorithm provides a fast and high-quality strategy for visualizing surfaces encoded by signed distance functions (SDFs), which have become a centerpiece in a wide range of visual computing algorithms. In this paper we introduce a sphere tracing algorithm for a completely different class of functions, harmonic functions, opening up a whole new set of possibilities. Harmonic functions are found throughout geometric and visual computing, where they arise naturally as the solution to interpolation problems, and in the physical sciences, where they appear as solutions to the Laplace equation. Our algorithm for harmonic functions is similar in spirit to the sphere tracing algorithm for SDFs: by using a conservative lower bound on the distance to the level set, we can take much larger steps than with naïve ray marching. Our key observation is that for harmonic functions such a bound is given by Harnack's inequality. Unlike Lipschitz bounds used in traditional sphere tracing, this Harnack bound requires only the value of the function at a point—we use this bound to develop a sphere tracing algorithm that can also handle jump discontinuities arising in angle-based harmonic functions. We show how this algorithm can be used to directly visualize smooth surfaces reconstructed from point clouds (via Poisson surface reconstruction) or polygon soup (via generalized winding numbers) without performing linear solves or mesh extraction. We also show how it can be used to render nonplanar polygons (including those with holes), and to visualize key objects from mathematics, including knots, links, spherical harmonics, and Riemann surfaces.