This dissertation analyzes the conventional haptic rendering methods for simulations involving interaction with deformable object, and propose a haptic rendering method which improves the conventional haptic rendering problems. This dissertation analyzes the stability and performance of the conventional haptic renderings using the stability condition and impedance transparency through a two-port network model. The virtual tool, the numerical integration method, and the contact stiffness are modeled with an impedance that can be interpreted in the frequency domain through z transform. The stability condition of a simulation involving the virtual tool, the numerical integration method, and the contact stiffness is verified with an experiment. The standard deviation of the stiffness value obtained from the stability condition obtained in this dissertation is 0.52 N/m less than the passivity condition 4.92 N/m and the conventional stability condition 4.82 N/m. The stability condition obtained in this dissertation better shows the tendency of the maximum stable stiffness value over time interval than the conventional conditions. This dissertation proposes a control architecture which does not have a trade-off between stability and performance based on the stability and performance analysis of the conventional haptic renderings. The proposed control architecture solves the stability problem caused by discretization by computing the feedback force in a separate haptic loop. The feedback force is computed by using the equivalent impedance representing the impedance of the interacting object and the user input. This dissertation develops a haptic rendering method to compute the feedback force using the local region of the contact surface of the deformable object in order to implement the equivalent impedance of the proposed control architecture. The local area is determined based on the contact point and a local stiffness matrix is constructed when a collision occurs between the virtual tool and the deformed object. The haptic feedback is quickly computed by using the local stiffness matrix instead of the entire stiffness matrix. This dissertation develops methods to reduce the error caused by computing the feedback force using the local stiffness matrix. A condensation method is developed to reduce the error without increasing the amount of computation in the haptic loop by further condensing the neighbor region to the local stiffness matrix. A method of condensing on a point-by-point basis is developed to reduce the computational load of the condensation in the simulation loop. An equivalent spring is added to the boundary of the local area to replace the excluded area in the local area. The coefficients of the equivalent springs are determined so that the stiffness energy of the equivalent spring equals the stiffness energy of the region excluded from the local region. Interaction simulation was performed using the Stanford bunny object consisting of 2087 points and 9997 elements and a box object consisting of 1183 points and 4320 elements. The proposed method reduces the x, y, and z axis errors by 79%, 79%, and 83%, respectively, compared to the virtual coupling in the simulation that involving interaction with the box object and by 52%, 80%, and 70%, respectively, in the simulation involving interaction with the Standard bunny object. Both simulation result shows that the measured energy flow satisfies the passivity condition. The proposed haptic rendering method provides more sophisticated haptic feedback at 1 kHz update rate than the conventional haptic rendering by 3ms and 5ms additional computation of the simulation loop in the simulation involving interaction with the Stanford bunny object and the box object, respectively.