Reverse osmosis (RO) membrane filtration systems have numerous applications in seawater desalination and advanced water treatment. However, fluid flow and mass transport in RO membrane filtration systems are difficult to simulate numerically due to the coupling between the transmembrane filtrate flow, pressure, and osmotic pressure. This coupling is challenging for real-world systems in which the local filtrate flow varies both spatially and temporally. Consequently, previous simulations tend to focus on two-dimensional steady cases. We present a new, high-accuracy numerical method specifically tailored to simulating RO. Using the method, we simulate the transient accumulation of rejected solutes at membrane surfaces and the resulting formation of concentration polarization layers. Specifically, we investigate polarization in plate and frame RO systems as a function of inlet Reynolds number, inlet solute concentration, and outlet transmembrane pressure. This talk will cover the numerical techniques used in our model and the unsteady results as they compare to current analytical and numerical models.