A non-linear diffusion model of Fujita is adapted to obtain an analytic solution describing the temperature distribution and position of any number of phase boundaries as a material cools on an effectively semi-infinite base material. Each material is initially homogeneous and at a uniform temperature. The solution method may incorporate any materials with temperature-dependent thermal properties undergoing any number of phase changes. As an example, we incorporate transitions through five phases of iron with non-linear heat conduction, as the iron cools on a copper base.