The successful detection of gravitational waves is the product of a combination of experimental technological advances and theoretical research breakthroughs. The waveform theoretical model of the wave source plays a key role in the experimental data processing. In order to build the waveform theoretical model, one needs to solve Einstein's equations. Due to the extreme strong field and strong dynamic spacetime behavior of gravitational wave sources, common approximation methods such as post-Newtonian approximation and perturbation methods encounter various limitations in solving the Einstein equations corresponding to gravitational wave sources. In contrast, the analytical solution of Einstein's equation is exceptionally difficult. So far, the only solutions with typical astronomical significance are the Kerr-Newman black hole solution (Schwarzy black holes, etc. are its special cases) and the FRW cosmological solution. The use of large-scale scientific calculations for the solution of the corresponding Einstein equations was then born as a general tool for modeling gravitational wave sources. Unlike other topics in computational mathematics and computational physics, numerical relativity typically faces the problem of computational stability in addition to the problems of computational accuracy and efficiency. The problem of computational stability makes the study of numerical relativity difficult. The existing algorithms that can maintain the computational stability include BSSN algorithm, Z4c algorithm, CCZ4 algorithm and GHG algorithm. Typical numerical relativity software includes AMSS-NCKU, Einstein-Toolkit, SpEC, etc.