Decades of intensive research effort has been dedicated towards developing more accurate, low computational cost turbulence models. However, due to a current plateau in progress, turbulent flows are still often modelled using Reynolds-averaged Navier-Stokes (RANS) approaches that give unsatisfactory predictions for complex flows. More recently, advancement in machine learning (ML) has propelled the use of data-driven methods for modelling Reynolds stress. However, many of these methods do not predict more than one region of flow physics accurately. In this work, a zonal ML approach for RANS turbulence modelling based on the divide-andconquer technique was deployed. This approach involved partitioning the flow domain into regions of flow physics called zones, training one ML model in each zone, then validating and testing them on their respective zones. The approach was demonstrated with the tensor basis neural network (TBNN) and another neural net called the turbulent kinetic energy neural network (TKENN) [1]. These were used to predict Reynolds stress anisotropy and turbulent kinetic energy respectively in test cases of flow over a solid block, which contain regions of different flow physics including separated flows. The results show that the combined predictions given by the zonal TBNNs and TKENNs were significantly more accurate than their corresponding standard non-zonal models. Most notably, shear anisotropy component in the test cases was predicted at least 20% and 55% more accurately on average by the zonal TBNNs compared to the non-zonal TBNN and RANS, respectively. The Reynolds stress constructed with the zonal predictions was also found to be at least 23% more accurate than those obtained with the non-zonal approach and 30% more accurate than the Reynolds stress predicted by RANS on average. These improvements were attributed to the shape of the zones enabling the zonal models to become highly locally optimized at predicting the output.