when you put tag open, you should not compute the next v_perp_ssdna1 because it might change in the wrong way the direction of the first bead when ind=0 (ind is c+1)
R = utils.get_rotation_matrix(dist[ind, :], gamma)
v_perp_ssdna1[ind, :] = np.dot(R , p[ind, :])
v_perp_ssdna2[ind, :] = -v_perp_ssdna1[ind, :]
when you put tag open, you should not compute the next v_perp_ssdna1 because it might change in the wrong way the direction of the first bead when ind=0 (ind is c+1) R = utils.get_rotation_matrix(dist[ind, :], gamma) v_perp_ssdna1[ind, :] = np.dot(R , p[ind, :]) v_perp_ssdna2[ind, :] = -v_perp_ssdna1[ind, :]