Late-stage seismic slip in geothermal reservoirs has been shown as a potential mechanism for inducing seismic events of magnitudes to ~2.6 as late as two decades into production. We investigate the propagation of fluid pressures and thermal stresses in a prototypical geothermal reservoir containing a centrally located critically stressed fault from a doublet injector and withdrawal well to define the likelihood, timing, and magnitude of events triggered by both fluid pressures and thermal stresses. We define two bounding modes of fluid production from the reservoir. For injection at a given temperature, these bounding modes relate to either low- or high-relative flow rates. At low relative dimensionless flow rates the pressure pulse travels slowly, the pressure-driven changes in effective stress are muted, but thermal drawdown propagates through the reservoir as a distinct front. This results in the lowest likelihood of pressure-triggered events but the largest likelihood of late-stage thermally triggered events. Conversely, at high relative non-dimensional flow rates the propagating pressure pulse is larger and migrates more quickly through the reservoir but the thermal drawdown is uniform across the reservoir and without the presence of a distinct thermal front, and less capable of triggering late-stage seismicity. We evaluate the uniformity of thermal drawdown as a function of a dimensionless flow rate QD that scales with fracture spacing s (m), injection rate q (kg/s), and the distance between the injector and the target point L∗ (QD ∝ qs2/L∗). This parameter enables the reservoir characteristics to be connected with the thermal drawdown response around the fault and from that the corresponding magnitude and timing of seismicity to be determined. These results illustrate that the dimensionless temperature gradient adjacent to the fault dTD/dxD is exclusively controlled by the factor QD. More significantly, this temperature gradient correlates directly with both the likelihood and severity of triggered events, enabling the direct scaling of likely magnitudes and timing to be determined a priori and directly related to the characteristics of the reservoir. This dimensionless scaling facilitates design for an optimum QD value to yield both significant heat recovery and longevity of geothermal reservoirs while minimizing associated induced seismicity.