In the process of development and utilization of a large-scale borehole heat exchanger (BHE) array system, the phenomenon of load shifting within BHE array can be observed. In this paper, OpenGeoSys software coupled with TESPy toolkit is used to establish a comprehensive numerical model of BHE system (without depicting the heat pump part), and the behaviors of load shifting between BHEs with different design parameters are studied. The results show that the outlet temperature of single BHE and BHE array is generally rising, and the soil temperature around the BHE has accumulated unbalanced heat. The soil temperature near the BHEs array fluctuates more obviously than the single BHE system, and the distribution is uneven. At the end of the 15th year, the soil temperature near the center BHE increased by 2 °C compared with the initial soil temperature, which was more favorable in winter, but was not conducive to the performance improvement in summer. Further analysis by changing the inter-borehole spacing shows that with the increase of the inter-borehole spacing, the load shifting behaviors are gradually weakened, and the maximum shifted load of the central BHE is linear with the change of the inter-borehole spacing. After changing the layout methods, we observe that the more intensive the layout is, the more load shifting behavior is, and the unbalanced rate of soil temperature distribution around the linear layout is lower than other layouts. With the increase in the number of BHEs, the load shifting behaviors are further enhanced. By analyzing the proportion of shifted load amount relative to the average value, it is found that the system will take a longer time to reach heat balance with the increase of BHEs’ number. A shutdown of part of BHEs for a certain period of time will help to improve the long-term operational efficiency of the large-scale shallow ground source heat pump (GSHP) system.