We used 2D coupled sine-Gordon equations combined with 3D heat diffusion equations to numerically investigate the thermal and electromagnetic properties of a 250 × 70 µm 2 intrinsic Josephson junction stack. The 700 junctions are grouped to 20 segments; we assume that in a segment all junctions behave identically. At large input power a hot spot forms in the stack. Resonant electromagnetic modes, oscillating either along the length ((0, n) modes) or the width ((m, 0) modes) of the stack or having a more complex structure, can be excited both with and without a hot spot. At fixed bath temperature and bias current several cavity modes can coexist in the absence of a magnetic field. The (1, 0) mode, considered to be the most favorable mode for THz emission, can be stabilized by applying a small magnetic field along the length of the stack. A strong field-induced enhancement of the emission power is also found in experiment, for an applied field around 5.9 mT.