The paper demonstrates the development of a numerical algorithm to solve a set of nonlinear hyperbolic equations in time domain in Eulerian frame of reference using a characteristics based finite difference implicit scheme to analyze density wave oscillations in boiling water nuclear reactor. The present algorithm removes the uncertainties existing in literature over the treatment of boundary conditions while simulating parallel channel instability of a reactor core by providing requisite mathematical support. The model is used to simulate parallel channel instability of a boiling water reactor core undergoing in-phase and out-of-phase modes of oscillations for both forced and natural circulation systems and numerical investigation confirms the existence of type-I and type-II instabilities in appropriate conditions. Next, the numerical simulations are conducted to evaluate the relative dominance of in-phase and out-of-phase modes of oscillations under various operational regime.