Many applications require diode lasers that are able to generate high output powers, for example in LiDAR (Light Detection and Ranging) systems or in future inertial confinement fusion power plants. In order to achieve high powers while maintaining low injected current densities, it is beneficial to stack multiple active regions and separate them using tunnel junctions. As the tunnel junctions are highly absorbing due to the small band gap and high doping levels, their accurate positioning in the nodes of the vertical mode is crucial. In order to answer important design questions such as the optimal distance between active regions and the tunnel junctions, a simulation tool is required. For our tool we solve the standard drift-diffusion equations augmented with a non-local band-to-band tunneling model and additional terms to describe effects such as stimulated emission. In this contribution we present the simulation model and demonstrate its application to a multi-junction laser diode.