## 6 Conclusions  
In this study, we assumed a definition of \"multiplet\" specifically tuned for the application to the seismicity of our study area. A computer code was developed for counting the number of multiplets detected applying such a definition to any earthquake catalog. For the CPT15 1650-2020 catalog, the code detected eight multiplets, which is a number significantly higher than the average number of multiplets detected in the same way on a set of 500 data sets obtained randomizing the occurrence times of the original catalog (i.e., 4.63 \(\pm\) 1.88 in Table 5). The result of this comparison supports the hypothesis that the occurrence of sequences containing multiple mainshocks is not just a casual circumstance. In this study, we also developed a new earthquake simulation code, paying particular attention to the enhancement of stress interaction among rupturing fault elements, and increased the number of multiplets in the simulated catalog. In this way, the number of multiplets detected by the above mentioned computer code on a 100,000-year simulated catalog (176) is 2.13 times higher than the average number of multiplets detected in the same way on a set of 500 randomized catalogs (Table 5). Besides the production of a significant number of multiplets, the simulated catalog exhibits long- and short-term spatiotemporal features that can be considered realistic imitations of those commonly observed in the real seismicity. We use a stacking procedure in which we compute the number of \(M_{w}\geq 4.2\) events that preceded and followed an \(M_{w}\geq 5.2\) earthquake in bins of five years considering the origin time at the time of every strong event in the 100 keyears simulated catalog. Our results related to northern and central Apennines show an acceleration of seismic activity some centuries before a mainshock, a modest quiescence starting 50 years before the mainshock and a strong aftershock occurrence in the following five years (Figure 7). In the same way, we analyzed the short term patterns in periods of about one month before and after every \(M_{w}\geq 5.2\) earthquake. Our results confirm the capacity of the simulator  
Figure 8: **(a) Stacked number of \(M_{w}\geq 4.2\) earthquakes that preceded and followed an \(M_{w}\geq 5.2\) earthquake within an epicentral distance of 50 km in 100,000 years simulated catalog, zooming on a time scale spanning only 0.1 years (36.5 days); (b) average b-value in time bins of 0.365 days before and after an earthquake of \(M_{w}\geq 5.2\) containing at least 10 events. **code to reproduce typical foreshocks-aftershocks sequences (Figure 8a). Additionally, in the simulated catalog for northern and central Apennines, the average b-values show a decrease lasting a few weeks before the strong earthquakes, followed by an instantaneous increase at the time of the earthquakes (Figure 8b). This pattern was observed by [PERSON] et al. (2015), [PERSON] et al. (2016), [PERSON] and [PERSON] (2017) in real earthquake sequences. The following items are available online at [[https://www.mdpi.com/article/10.3390/app12042062/s1](https://www.mdpi.com/article/10.3390/app12042062/s1)]([https://www.mdpi.com/article/10.3390/app12042062/s1](https://www.mdpi.com/article/10.3390/app12042062/s1)), Table S1: Geometric and kinematic parameters of the 198 quadrilateral fault segments derived from the 43 Composite Seismogenic Sources (CSS) of DISS v. 3.3.0. Figure S1: Cumulative and density magnitude-frequency distributions of \(M_{w}\geq 4.2\) earthquake simulated catalogs, for the nine combinations of free parameters of Table 3 considered in this study. The straight dotted lines show the best-fit Gutenberg-Richter distributions. Figure S2: Stacked number of \(M_{w}\geq 4.2\) earthquakes that preceded and followed an \(M_{w}\geq 5.2\) earthquake within an epicentral distance of 20 km in the 100,000 years simulated catalog, for nine combinations of free parameters of Table 3 considered in this study.