In this work, an improved domain decomposition method is developed to address workload imbalance when implementing the parallel computing of a four‐dimensional lattice spring model (4D‐LSM) to solve problems in rock engineering on a large scale. A cubic domain decomposition scheme is adopted and optimized by a simulated annealing algorithm (SAA) to minimize the workload imbalance among subdomains. The improved domain decomposition method is implemented in the parallel computing of the 4D‐LSM. Numerical results indicate that the proposed domain decomposition method can further improve the workload balance among processors, which is helpful to supersede the limit of computational scale when solving large‐scale geotechnical problems and decrease the runtime of the parallel 4D‐LSM by at most 40% compared to the original cubic decomposition method. This shows the practicability of the proposed method in parallel computing. Two types of target functions of SAA are tested, and their influence on the performance of the parallel 4D‐LSM is investigated. Finally, a computational model with one billion particles for one actual engineering application of using 4D‐LSM is realized, and the result shows the advantages of parallel computing.