Unexpected numerical damping often occurs when an implicit time-integration is employed for stiff simulation model with large time step. The energy decreases when the equation is established using the information of the next time step when estimating the integral value. This paper presents a novel method for stable numerical time-integration in simulation of stiff dynamic systems with large time step. The proposed method, first, linearly interpolates the position solutions computed by the forward and backward Euler integrations. Derivative of the solution is, then, numerically obtained to compute the velocity. Blending ratio between the two solutions is determined by employing energy conservation equations and constraints. This makes the resulting blended solution closer to the solution obtained from the backward Euler integration. The constraints allow the blended solution to have physic-based behavior and low error. Two ratios are independently computed for the position and the velocity. The conventional blending and the proposed blending methods are compared in terms of energy exchanges and errors. A solution from an off-line accurate FEM is used to compare the errors in the energy norms during simulation of free vibration. The proposed time-integration method resolves the unplanned numerical damping with small additional computation time, and results in low errors during the stiff system simulation with large time step.