## 1 Introduction  
The human heart is a complex electro-pump that, through cycles of depolarization and repolarization, makes the propagation of the action potential, and consequently the contraction of the cardiac muscle tissue, possible. The electrical activity of the heart can be indirectly measured by the electrocardiogram (ECG), a non-invasive method based on 12 electrodes positioned on the body surface [1, 2]. Considering the crucial role that the ECG signal plays in the clinical practice, different techniques were presented in the literature in order to model the dynamics of the heartbeat. The electrical and muscular activities of the heart were analyzed through both mathematical modeling and time series analysis [3]. By using mathematical modeling, it is possible to understand the electrical activation of the heart obtaining, by simulation, both normal and abnormal rhythms [4, 5, 6]. The generation of synthetic ECG signals with a wide range of waveform shapes and heart rates allows the modeling of the characteristics of each subpart composing the cardiac conduction system, and the understanding of its behaviour [7, 8, 9, 10, 11]. The nonlinearity and nonstationarity of the cardiovascular system make the use of nonlinear techniques for the modelling of heart activity useful [12, 13, 14]. [PERSON] (VdP) and [PERSON] described and modeled the behavior of heart using nonlinear relaxation oscillators [15]. [PERSON] [16] proposed an extended version of the Van der Pol equations for the generation of the action potential. In his work, he proposed a model for the emulation of the signal observed in an excitable cell of a living organism. The model is composed of two coupled, nonlinear ordinary differential equations. The first equation describes the evolution of the neuronal membrane voltage (fast action), while the second one the recovery action of the sodium and potassium channels (slow action) [17]. [PERSON] and [PERSON] modified the Van der Pol model to allow the reproduction of the time series of the action potential generated by a natural pacemaker [18; 19]. In their model, the intrinsic frequency of the two pacemakers can be changed without the need to change the length of the refractory period. Although the model allows the manipulation of both the diastolic and the refractory period, the frequency obtained is too low if compared with the physiological values. [PERSON] and [PERSON] [20] proposed a mathematical model always based on the Van der Pol equations, but composed of three modified oscillators. This model is able to emulate the cardiac conduction system composed of a Sino-Atrial node (\(SA\)), Atrio-Ventricular node (\(AV\)) and His-Purkinje system (\(HS\)). The oscillators are connected to each other through time-delayed couplings. This model allows for generating the electrical response of the main cardiac pacemakers and obtaining the ECG wave as a composition of these signals. Although this model reproduces normal and pathological ECG signals, it does not model the activity of atrial and ventricular muscles. The characterization of muscle electrical responses proposed by [PERSON] [21] uses a quiescent excitable FitzHugh-Nagumo-type (FHN) oscillator. Successively, in [22; 23], the authors improved this model by including the depolarization and repolarization waves of the atrial and ventricles generated by modified FHN systems for each ECG wave. The advanced model was able to generate normal ECG signals and several well-known rhythm disorders. In particular, it reproduces sinus tachycardia, sinus bradycardia, complete \(SA\)-\(AV\) block and complete \(AV\)-\(HP\) block [24]. In the models proposed in [20; 21; 22], the cardiac conduction system was treated as a network of three modified Van der Pol oscillators representing \(SA\), \(AV\) nodes and the \(HP\) system.