The numerical model of fluid-elastic instability of square tube array is studied in this paper. A 3D fluid structure interaction vibration model is established based on rigid body motion function and Newmark integral method. Influences of turbulence model and axial dimension are analyzed. A proper analysis model is obtained. The model is adopted computing critical velocity and the result is compared with reference experiment data. The result shows that computed critical velocity is well coincidence with experiment data. Computed amplitudes are precise at low velocity, but have a lager error at high velocity. The work provides reference and foundation for numerical model application.