Abstract. This paper proposes an improved algorithm for simulating the surface flow dynamics based on the flow-path network model. This algorithm utilizes the parallel-multi-point method to extract the critical points and the D8 algorithm to retrieve the drainage networks from the regular-grid digital elevation model (DEM) for constructing a drainage-constrained triangulated irregular network (TIN). Then, it combines the flow directions of triangular facets over TIN with resampled flow source points to track flow lines to generate the flow path network (FPN) based on the flow-path network model. On this basis, the proposed algorithm employs three terrain parameters (slope length factor, topographic wetness index and flow path curvature) to improve the classical Manning equation based on the analytic hierarchy process (AHP) to enhance the accuracy of the flow velocity calculation. The topographic wetness index and flow path curvature are derived by the flow-path-network-triangular-facet-network (FPN_TFN) algorithm, a new flow-path-network-topographic-wetness-index (FPN_TWI) algorithm and the flow-path-network-flow-path-curvature (FPN_C) algorithm, respectively. Finally, the velocity estimation function and surface flow discharge simulation function are parallelized by the Compute Unified Device Architecture (CUDA) to enhance its computational efficiency. The outcomes are compared with the algorithm before improvement (TIN_based algorithm) and the SWAT model. The results demonstrate that the speedup ratio reaches 15.7 compared to the TIN_based algorithm. The Nash coefficient increases by 6.49 %, the correlation coefficient decreases slightly, and the balance coefficient increases by 19.08 %. Compared with the SWAT model, the Nash coefficient and correlation coefficient increase by 97.56 % and 4.60 %, respectively. The balance coefficient is close to 1 and outperforms the compared algorithms.