The van der Pol physical reservoir computer

The van der Pol oscillator has historical and practical significance to spiking neural networks. It was proposed as one of the first models for heart oscillations, and it has been used as the building block for spiking neural networks. Furthermore, the van der Pol oscillator is also readily implemented as an electronic circuit. For these reasons, we chose to implement the van der Pol oscillator as a physical reservoir computer (PRC) to highlight its computational ability, even when it is not in an array. The van der Pol PRC is explored using various logical tasks with numerical simulations, and a field-programmable analog array circuit for the van der Pol system is constructed to verify its use as a reservoir computer. As the van der Pol oscillator can be easily constructed with commercial-off-the-shelf circuit components, this PRC could be a viable option for computing on edge devices. We believe this is the first time that the van der Pol oscillator has been demonstrated as a PRC.


Introduction
The van der Pol oscillator was first described approximately a century ago [1]. It was later proposed as a model for heart oscillations [2]. Since then, it has become a prolific model in a wide range of systems. Other than the heart, the van der Pol model has been used to describe the driven dissipative nonequilibrium quantum systems [3][4][5], vocal cords [6], neurons [7], and single reed instruments [8]. As the van der Pol oscillator is easily implemented as a circuit, there have been many different types of experimental circuit models as well. These including circuits using vacuum tubes [9], operational amplifiers [10], resistor-inductor-capacitor [11], memristors [12], field programmable gate arrays [13], and field programmable analog arrays (FPAAs) [14]. In part because of its ability to be constructed as a circuit, the chaotic behavior of the van der Pol oscillator has been studied extensively [15][16][17][18]. Moreover, a comprehensive history leading to the discovery of this relaxation oscillator is provided by [19].
In relationship to machine learning, the van der Pol oscillator is often used as a benchmarking waveform for neural networks [20][21][22]. Neural networks have been used to solve a modified van der Pol heart model [23]. Van der Pol time series data has also been used to benchmark reservoir computers [24,25]. Conversely, networks of van der Pol oscillators have been constructed to produce central pattern generators for walking robots [26]. Van der Pol oscillator arrays have been used as spiking neural networks, where the nodes are van der Pol oscillators [27,28].
Historically, non-physical reservoirs are a set of hidden layers within a deep neural network, which has untrained weights. The more traditional technique of non-physical reservoir computing is used for many applications, including dynamic spectrum management [29], network traffic prediction [30], and unmanned aerial vehicle (UAV) trajectory design [31]. In general, a non-physical reservoir computing architecture has several advantages that include robustness to overfitting and lower training costs.
As a natural extension of non-physical reservoirs, physical reservoir computers (PRCs) are a form of artificial intelligence that repurposes the nonlinearity of a physical dynamic system to perform computation [32][33][34][35]. PRC is the product of merging two machine learning techniques, 'echo state networks' [36] and 'liquid state machines' [37]. PRCs are usually composed of the input (an external force that perturbs the nonlinear system), the reservoir (the nonlinear physical system), and the readout function. The computation is performed by mapping a lower dimensional input space to a higher dimensionality via a nonlinear system, which is known as the reservoir. In addition to the advantages already mentioned for non-PRCs, PRCs have several other notable advantages. As it is difficult to implement advanced machine learning techniques on edge devices (such as smart phones) [38], PRCs provide an attractive alternative to conventional neural networks because of their faster learning, a simpler architecture, and cheaper training costs [39]. A single PRC can replace a sizable neural network because of the underlying nonlinear dynamics of the reservoir. Furthermore, many complex systems have integral components that can be repurposed for computation (e.g. the microphone in a smart phone can be repurposed as a compute-on-the-sensor audio recognition device [40]). As one additional advantage, analog systems can simultaneously provide both memory storage and computation.
Rather than use an array of van der Pol oscillators whose connections must be trained, we propose a simple van der Pol PRC. The problem of interest is to quantify the computational ability of this simple PRC and provide limits to its working range. This knowledge can be used to inform the construction of PRCs from existing physical systems, which can be modeled as van der Pol oscillators. Alternatively, this knowledge also provides design considerations for engineered van der Pol oscillators that can be implemented as electronic circuits [10] or quantum systems [5], which can be used as PRCs. As such, this van der Pol PRC approach could have applications to powerful reconfigurable computation for edge devices.
Instead of using a van der Pol array as a neural network, this paper will explore a single van der Pol PRC by discarding the regularly used delay line [35]. Although quite simple when compared to a full neural network, this PRC is still capable of computational tasks and memory storage. Moreover, this PRC is created from one of the most widely used nonlinear oscillators, the van der Pol oscillator, which can readily be implemented as an electronic circuit. The results presented here provide a framework for further exploration of van der Pol PRCs. Specifically, this work highlights the potential to leverage the van der Pol PRC as a single node in a spiking neural network, which could provide substantially different network architectures and training procedures. On the other hand, the van der Pol oscillator is a prolific model for a wide variety of physical systems, as mentioned previously. Because of this, it is likely that many of these van der Pol-like systems could be repurposed as PRCs (e.g. the heart as a computer).
The rest of the paper is organized as follows. In section 2, the equations of motion for the van der Pol PRC are provided; a description of how information is encoded and sent to the PRC are also described. In section 3, numerical simulations are used to gain an understanding of how system parameters affect the computational abilities of the van der Pol PRC. A logical task is used to benchmark this computer, so Shannon's information rate is used as the performance metric. In section 4, a FPAA circuit is described for the van der Pol PRC, and an experimental validation of the efficacy of the van der Pol PRC is provided. In section 5, concluding remarks are provided.

Van der Pol oscillator
The van der Pol oscillator has a nonlinear damping term, which causes limit cycle motion. The forced van der Pol oscillator is described by the following equation [54]: In this form of equation (1), the van der Pol oscillator has only one system parameter, µ. Here, f (t) is an external force on the van der Pol oscillator. The full external force, f (t), is constructed such that it provides both a periodic component to entrain the van der Pol oscillator, f s (t), and an encoded signal, f u (t). These two portions of the full external force are discussed in the next subsection. f (t) is written as: In this paper, the velocity state of equation (1) is used generate virtual nodes. The velocity state, y, may be seen more easily when equation (1) is written in state space as:

External forcing
In this subsection the two portions of the full external force, f u and f s , are described. The f u (t) force is used to encode 'True' and 'False' bits of information, which are discrete, in the continuous domain to send this information to the oscillator. Whereas, the f s (t) force is used to synchronize the oscillator to the external force. In some sense, this sinusoidal force acts similar to a clock frequency. The f u (t) force is defined as: The function u(n) is a discrete function that is a random series of +1 (which corresponds to 'True') and −1 (which corresponds to 'False') values. In the next portions of the paper, these random 'True' and 'False' values will be used to calculate various logical tasks. This discrete function is used to define a continuous function, f u (t). This encoded signal function takes the value of u(n − 1) for a specific, arbitrary length of time, which is defined as the pseudo-period, T p . This pseudo-period, T p , and the encoding amplitude, a u , are constants. The pseudo-period is also used to define a pseudo-frequency, such that ω p = 2π Tp . During a pseudo-period, f u (t) holds the value of u(n − 1).
The sinusoidal portion of the external force, which is used to synchronize the oscillator with the external force, is given as: Here, a s is the sinusoid's amplitude, and Ω is the forcing frequency. For this van der Pol PRC, there are five parameters: µ, a u , a s , Ω, and ω p (which corresponds to T p ). Although an additional parameter is sometimes included before the x term in equation (1), it was set equal to unity here for simplicity.

Training procedure and information rate
A conceptual schematic for the van der Pol PRC is shown in figure 1. On the left, a non-physical reservoir is shown, which has random connections of the nodes. On the right, the van der Pol oscillator PRC is shown. The van der Pol oscillator is perturbed with the external force, f (t). The oscillator's dynamics are then used to collect virtual nodes over each pseudo-period at the readout, Σ (here, there are five virtual nodes for visualization purposes). These virtual nodes are then trained to calculate different logical tasks.
For the following simulations and experiment, 100 virtual nodes were collected over each pseudo-period from the y state. A matrix of these virtual nodes was then trained with a simple ridge regression with a regularization penalty equal to 0.001. The target for the training is an exclusive OR (XOR) task. This parity task is a common benchmark for PRCs. Eighty percent of the data was used for training, and the remaining 20% was used for testing. After the reservoir was trained in the continuous domain, the output was averaged over each pseudo-period to translate the continuous data back to a discrete domain.
For logical task, the predicted value must be 'True' or 'False' . Because of this, the root mean square error is an ill-posed metric, as it inherently measures the distance between numerical values. For this reason, Shannon's information rate [55] is modified to measure the performance of the van der Pol PRC at computing logical tasks. This metric is well-posed for logical operations, but this metric also rather severely penalizes incorrect calculations. The information rate, IR, is defined as IR = H(x) − H y (x), where the Shannon entropy is H(x) and the conditional entropy is H y (x). The Shannon entropy quantifies the amount of information in a sent signal; for this modified calculation, the Shannon entropy is the amount of information from the output of the logical task. The conditional entropy quantifies the probability of an incorrect bit in the received signal; for this modified calculation, the conditional entropy is the probability of an incorrect bit in the predicted output of the logical task. In the following equations, p i is the probability of getting a particular bit, i; this provides the following equations: To modify this for the logical task here, a bit i is from the target and a bit j is from the prediction. Here, p i (j) = p(j|i) = p(i,j) ∑ j p(i,j) , and p(i, j) is the joint probability distribution of the two variables. Since the target was an XOR task and u(n) was a random binary input, the Shannon entropy is equal to 1 for this task; thus, the maximum value of the information rate for this task is also equal to 1. In all figures, only the resulting information metric for the testing data was plotted.

Parametric studies using numerical simulations
In this section, numerical simulations are used to probe the relationship between parameter values and the computational ability of the van der Pol PRC. An example of one of the state space plots is shown in figure 2, along with the XOR (second order parity task) prediction. For the example given, the van der Pol PRC predicts the XOR task with an IR equal to 1. In all these tasks, the van der Pol oscillator is synchronized with the f s force and perturbed by the encoded information force, f u . The van der Pol oscillator stores some of the previous information in its dynamic states, and its response acts as a type of calculation. Interestingly, the response can be reprogrammed to act as different tasks.
In the remainder of this section, parametric studies are conducted to provide a practical basis for creating van der Pol PRCs. For instance, the effects of the limit cycle radius, µ, the relationship of frequencies, Ω and ω p , and bounds on amplitudes, a u and a s , are explored in figures 3-5, respectively. For each of these numerical experiments, other values of interest were held constant, and the same logical task (XOR) was used to test the computational ability of the PRC. For all cases, the modified Shannon's information rate was used to measure the performance.
The effects of µ on the computation are considered. For figure 3, a u = 1, a s = 1, ω p = π 2 , and Ω = ωp 1/4 . In this figure, a plateau of high IR may be seen, which extends from approximately 0.3 to 4. For the other chosen parameter values, the µ term has a very strong effect on the computational ability of the van der Pol PRC. When µ is near zero, it is likely that the nonlinearity of the system is too small to compute nonlinear tasks.
In figure 4, the effects of the frequency combinations, ωp Ω , are explored. It is found that the van der Pol PRC only has good performance when this frequency ratio is equal to a Farey sequence ratio. The Farey sequence is often found to be a necessary synchronization condition for nonlinear oscillators, so this suggests  that choosing the pseudo-frequency to forcing frequency ratio equal to a Farey sequence ratio is most likely to result in good performance of the PRC. It should also be noted, however, that many of the Farey sequence ratios do not correspond to high performance, so the frequency ratio being equal to a Farey sequence ratio is not a sufficient condition for a high IR. In figure 5, the amplitudes of the sinusoidal force, a s , and the encoded information, a u , are explored. The other parameter values are set as µ = 1, ω p = π 2 , and Ω = ωp 1/4 . It can be seen that there is a rather marked delineation between poor performance to maximal IR when a u is larger than approximately 1. Also, it  This figure suggests that the sinusoidal forcing, as, is necessary to properly synchronize the van der Pol PRC. It also shows that au needs to have a certain magnitude to convey information to the PRC. When au is bigger than approximately 5, the IR is unreliable.
Here, µ = 1, ωp = π 2 , and Ω = appears that a small, non-zero value of a s is necessary to aid in synchronizing the van der Pol system to the encoded information in f u , which is shown in figure 4. Larger values of a u give unreliable computation, but a broad region is found in which the IR is maximal for the PRC. In figure 6, the effects of the frequency combinations, ωp Ω , are explored by considering different orders of parity tasks. It is observed that the van der Pol PRC can store more information in its dynamics states when Figure 6. The effect of frequency on the memory capacity of the van der Pol oscillator is studied. By setting ωp equal to a constant, the relationship between Ω and ωp can be found. The tick marks in this plot are the 20th order Farey sequence ratios. From this plot, it can be seen that the van der Pol PRC mostly has higher memory storage capability when ωp Ω is equal to a Farey sequence ratio. Here, µ = 1, au = 1, as = 1, and ωp = π 2 .
this frequency ratio is equal to a Farey sequence ratio. For a δ-delayed, nth order parity function given by equation (7), the memory storage capacity can be calculated using equation (8) [49]: where δ ∈ Z + is the delay and Par n,δ is the nth order parity. The memory capacity, MC n , for the nth order parity task can be calculated as: p n,δ is the probability of predicting the correct bit for the δ-delayed nth order parity. In figure 6, parity tasks to the fifth order with a maximum delay of δ = 6 are considered to study the memory capacity performance of the reservoir.
So far, the van der Pol oscillator has been studied with the XOR gate (⊻) to gauge the effects of various system parameters. However, the van der Pol PRC is reprogrammable, which can be seen in figure 7. Here, a single van der Pol oscillator simulation is used, but its response is trained to predict the three basic logic gates: NOT (¬), AND (∧), and OR (∨). With these three gates, any other logic gates can be constructed. For each of these tasks, the PRC provided perfect performance after training (e.g. IR = H(x) for each task).

Van der Pol FPAA circuit and results
As they can be dynamically programmed, FPAAs are very useful in prototyping nonlinear analog circuits. FPAAs have been popularly used to implement many nonlinear oscillators. These include a pendulum adaptive oscillator [56], a four-state adaptive oscillator [57], a fractional order chaotic system [58], the Lorenz oscillator [59], a chaotic system with infinite attractors [60], the Sprott N chaotic system [61], the Nahrain chaotic system [62], and a hyperjerk oscillator [63]. This FPAA implementation of the van der Pol oscillator was chosen to provide a rapid prototype, which could be easily reproduced by other researchers. In addition, this FPAA implementation could be modified for specialty edge devices that make use of the van der Pol PRCs, such as ultrasound applications, audio recognition devices, and structural health monitoring. As the analog counterpart of a field programmable gate array, a FPAA is an analog device that uses switched-capacitor technology [64]. FPAAs contain configurable analog blocks for analog operations, while configurable analog modules (CAMs) can perform specific math operations. An FPAA circuit design of the van der Pol oscillator is shown in figure 8.
An Anadigm QuadApex v2.0 board was used, which has four AN231E04 chips. For the van der Pol oscillator, only a single AN231E04 chip was used. CAMs for summation, multiplication, DC voltage (denoted by a circle with a ± symbol inside), integration, and sample-and-hold (denoted by a Z −1 symbol) were used to make the circuit. The external force, f (t), was first generated in MATLAB. Then, a National Instruments 9263 module was used to send f (t) to the FPAA. A National Instruments 9201 module was used to collect the x and y states of the van der Pol oscillator to MATLAB.
The phase portrait of the experiment is shown in figure 9. The same machine learning method that was used for the simulations in section 2 was also used for the experimental results. For the experimental results shown in figure 9, the IR was equal to 0.994. This provides an experimental validation of the FPAA circuit for the van der Pol PRC. However, it should also be noted that the phase portrait shown in figure 9 has somewhat 'sharp' corners. This is caused by the voltage limitations of the circuit (e.g. the circuit rails at

±3.3 V)
. Even with this physical limitation, the experimental van der Pol PRC still works with almost perfect accuracy for this XOR benchmark task.

Conclusions
In this manuscript, a van der Pol PRC is proposed. An XOR task, different logical tasks, and memory capacity calculation were used to provide a benchmark for its computational ability. These tasks require that the PRC dynamically store the information of multiple bits, while simultaneously calculating the different logical tasks. The method of encoding 'True' and 'False' bits was defined, and the implementation was described.
Although the van der Pol oscillator is only a second order ordinary differential equation, it is capable of having complex behavior due to its nonlinearity. The van der Pol PRC that is presented here has no delay lines, but it is still capable of storing memory and completing a XOR task. The van der Pol oscillator's response is influenced by the parameters, µ, a u , a s , ω p , and Ω. Numerical simulations were used for parametric studies to shed some light on how the parameter values affect the computational ability of the van der Pol PRC. µ was found to provide a broad range of good computational ability, for values from approximately 0.3 to 4.0. When µ was very small or larger than 4, the computational ability was quite poor. By considering the ratio of ωp Ω , the relationship between these two frequencies to computational ability could be seen. Farey sequence ratios, which often appear in the context of synchronization in nonlinear oscillators, play an important role in memory and computation. Lastly, the computational ability was explored with different combinations of a u and a s . Through these parametric studies, parameter combinations can be found that provide robust performance of the van der Pol PRC.
The implementation of this PRC can be realized through many different types of circuits. This manuscript explored an FPAA circuit implementation, since FPAAs are dynamically programmable analog devices for circuit prototyping. The circuit diagram for an Anadigm FPAA was provided. This circuit implementation provides a hardware validation of the van der Pol PRC.

Data availability statement
All data that support the findings of this study are included within the article (and any supplementary files).