A mathematical model which incorporates the processes that influence water flow and heat transfer in subfreezing snow was developed. Among the aspects of snow included are density and grain‐size heterogeneities, capillary‐pressure gradients, meltwater refreezing, time dependent hydraulic and thermal parameters, and heat conduction. From this conceptual mathematical model a numerical model of two‐dimensional meltwater infiltration was developed. Results from various test cases show which data are most important to measure accurately in the field, in order to determine how the snowpack will respond to an introduction of meltwater. These simulations also show the importance of the orientation of the various layers which make up the snowpack and how randomly distributed heterogeneities can produce two‐dimensional flow of meltwater under unsaturated conditions. Finally, it is demonstrated that various assumptions related to density and porosity variations, dimensionality of flow, capillary effects, etc., which have been made by past investigators for ideal situations may not be valid under many circumstances, and several suggestions are made for improving predictions of meltwater behavior. Sensitivity analysis showed that the model is most sensitive to changes in bulk density, residual saturation of wet snow and meltwater supply rates, whereas changes in snow temperature and mean grain size had less marked effect. Copyright 1990 by the American Geophysical Union.