Methodological aspects of using the driven Liouville-von Neumann (DLvN) approach for simulating dynamical properties of molecular junctions are discussed. As a model system we consider a non-interacting resonant level uniformly coupled to a single Fermionic bath. We demonstrate how a finite system can mimic the depopulation dynamics of the dot into an infinite band bath of continuous and uniform density of states. We further show how the effects of spurious energy resolved currents, appearing due to the approximate nature of the equilibrium state obtained in DLvN calculations, can be avoided. Several ways to approach the wide band limit, which is often adopted in analytical treatments, using a finite numerical model system are discussed including brute-force increase of the lead model bandwidth as well as efficient cancellation or direct subtraction of finite-bandwidth effect. These methodological considerations may be relevant also for other numerical schemes that aim to study non-equilibrium thermodynamics via simulations of open quantum systems.