Multiple phase transitions, based on experimental and theoretical results, have been incorporated into an extended Boussinesq model of mantle convection. A formulation is used which facilitates the numerical computations of convection with multiple phase transitions and a triple point near the 670 km discontinuity. The propensity for layered convection is enhanced by including the additional phase transitions. The location of the triple point of the spinel to perovskite phase transitions plays an important role in controlling the fate of the ascending hot plumes. A plume with a core temperature above this triple point temperature Tc tends to pass easily through the phase boundary at 670 km. Plumes with a core temperature lower than Tc are blocked by the multiple phase transitions. Both depth-dependent thermal expansivity and internal heating increase the tendency of the convective system to layering. The dynamical time-scales in the transition zone are much shorter than those of the top and bottom boundary layers. Simulations with increasing Rayleigh number (Ra) from 106 to 5 × 107 show unambiguously the increased tendency toward layering with larger Ra. In the early Earth, when Ra was higher and internal heating was several times stronger, layered convection might have been the preferred mode of mantle convection. With time, as both Ra and radioactive heating decreased, the style of mantle convection would have become less layered.