An ice-formation algorithm is implemented in the three-dimensional Estuary and Lake Computer Model, to allow simulation of hydrodynamics and the thermal structure beneath the ice during winter. The one-dimensional governing equation of heat conduction among the three layers of white ice, blue ice, and snow is solved for the formation of ice cover considering the heat flux through air and water. This algorithm is applied independently in each grid cell within the simulation domain, allowing for spatially variable ice formation. The model was validated against observed data from both a large and a small Canadian mid-latitude lake (Lake Ontario and Harmon Lake, respectively). The lake surface temperature and the distribution and thickness of ice cover on Lake Ontario were predicted successfully during the 2006-2007 winter period. The model also accurately simulated spring 2007 temperature profiles, as typically used for the initial conditions for a summer simulation. The variation of ice and snow thickness, and vertical temperature profiles, were well-simulated for Harmon Lake during winter 1990-1991. These comparisons demonstrate the applicability of the model for year-round simulation of mid-latitude lakes of varying size.Coupled hydrodynamic and biogeochemical computational lake models have been developed for simulation of lake circulation and management of water quality (Chapra 1997; Hodges 2009). However, attempts at simulation of water quality during winter have had mixed results and remain unpublished (see discussion in Gosink 1987). More recently, Hamilton et al. (2002) applied the one-dimensional (vertical) hydrodynamic and water quality model Dynamic Reservoir Simulation Model-Water Quality (DYRESM-WQ), coupled with an ice model to simulate the changes in tropic status in a small lake, resulting from atmospheric nutrient deposition and climate change. To our knowledge, simulation of lake biogeochemistry beneath ice cover in three dimensions has not been attempted. Consequently, processes such as winter primary production, the winter mixing of phosphorus (for example, which released from sediments during late-autumn hypoxia or introduced during the spring freshet) remain comparatively unexamined. Ice cover significantly modifies lake hydrodynamics, inhibiting wind stress, so that vertical mixing is sustained only by natural convection (Farmer 1975). Therefore, correct modeling of lake dynamics during the ice-covered season requires accurate simulation of ice cover and its effects on heat and momentum transfer.Several thermodynamic models for ice formation have been developed over the past few decades. We find the existing models to be unsuitable for research on lake management, because they are either overly simplified three-dimensional (3D) ice-formation and transport models or more detailed one-dimensional (vertical) models that neglect water column dynamics and biogeochemistry. Most of these models are applying simplified assumptions for energy fluxes in the formation of ice (Maykut and Unterste...