A fully coupled land surface hydrologic model, Flux-PIHM, is developed by incorporating a land surface scheme into the Penn State Integrated Hydrologic Model (PIHM). The land surface scheme is adapted from the Noah land surface model. Because PIHM is capable of simulating lateral water flow and deep groundwater at spatial resolutions sufficient to resolve upland stream networks, Flux-PIHM is able to represent heterogeneities due to topography and soils at high resolution, including spatial structure in the link between groundwater and the surface energy balance (SEB). Flux-PIHM has been implemented at the Shale Hills watershed (0.08 km2) in central Pennsylvania. Multistate observations of discharge, water table depth, soil moisture, soil temperature, and sensible and latent heat fluxes in June and July 2009 are used to manually calibrate Flux-PIHM at hourly temporal resolution. Model predictions from 1 March to 1 December 2009 are evaluated. Both hydrologic predictions and SEB predictions show good agreement with observations. Comparisons of model predictions between Flux-PIHM and the original PIHM show that the inclusion of the complex SEB simulation only brings slight improvement in hourly model discharge predictions. Flux-PIHM adds the ability of simulating SEB to PIHM and does improve the prediction of hourly evapotranspiration, the prediction of total runoff (discharge), and the predictions of some peak discharge events, especially after extended dry periods. Model results reveal that annual average sensible and latent heat fluxes are strongly correlated with water table depth, and the correlation is especially strong for the model grids near the stream.