In view of stability requirements of compliant mechanisms, based on the Heaviside three-field mapping scheme and the compliance change rate constraint measure, the topology optimization problem of the compliance minimization of hinge-free multi-input multi-output compliant mechanisms with buckling constraints is investigated. Firstly, in order to solve the problem of pseudo-buckling modes and calculation efficiency caused by low-density elements, smooth penalty functions of the element stiffness matrix and stress matrix are constructed, and comprehensive measures to solve the corresponding problem are proposed. Furthermore, cooperating with the Pian mixed stress element formula, the aggregation function method, the Heaviside three-field mapping scheme and the compliance change rate constraint measure, the topology optimization model of hinge-free multi-input multi-output compliant mechanisms considering buckling is constructed. Then, the sensitivities of the objective function and buckling constraints are derived, and the MMA algorithm is adopted to optimize the topology optimization model. Finally, the numerical examples are given to illustrate the feasibility and effectiveness of this method.