Resonant Machine Learning Based on Complex Growth Transform Dynamical Systems

TL;DR

This paper introduces a novel energy-efficient machine learning framework inspired by electrical resonance, utilizing complex growth transform dynamical systems to optimize learning while maintaining zero reactive power.

Contribution

It proposes a new energy-based learning model combining active and reactive energies, driven by complex-domain growth transforms, with applications to resonant SVMs.

Findings

Supports the design of self-sustained resonant SVMs

Demonstrates active-power regularization improves learning

Shows convergence controlled by an annealing process

Abstract

Traditional energy-based learning models associate a single energy metric to each configuration of variables involved in the underlying optimization process. Such models associate the lowest energy state to the optimal configuration of variables under consideration, and are thus inherently dissipative. In this paper we propose an energy-efficient learning framework that exploits structural and functional similarities between a machine learning network and a general electrical network satisfying the Tellegen's theorem. In contrast to the standard energy-based models, the proposed formulation associates two energy components, namely, active and reactive energy to the network. This ensures that the network's active-power is dissipated only during the process of learning, whereas the reactive-power is maintained to be zero at all times. As a result, in steady-state, the learned parameters are stored and self-sustained by electrical resonance determined by the network's nodal inductances and capacitances. Based on this approach, this paper introduces three novel concepts: (a) A learning framework where the network's active-power dissipation is used as a regularization for a learning objective function that is subjected to zero total reactive-power constraint; (b) A dynamical system based on complex-domain, continuous-time growth transforms which optimizes the learning objective function and drives the network towards electrical resonance under steady-state operation; and (c) An annealing procedure that controls the trade-off between active-power dissipation and the speed of convergence. As a representative example, we show how the proposed framework can be used for designing resonant support vector machines (SVMs), where we show that the support-vectors correspond to an LC network with self-sustained oscillations.

Tables

Table 1. TABLE 1: Notations

Variable	Definition
${ℝ}_{+}$	One-dimensional positive real vector space
${ℝ}^N$	N-dimensional real vector space
$ℂ$	One-dimensional complex vector space
$\left|z\right|$	magnitude of a complex variable $z$
$Re\left(z\right)$	real part of a complex variable $z$
$Im\left(z\right)$	imaginary part of a complex variable $z$
$z^{\ast }$	complex conjugate of a complex variable $z$
$z\left(t\right)$	a continuous-time complex variable at time $t$
$z_n$	a discrete-time complex variable at the $n^{\text{th}}$ time step

Table 2. TABLE 2: Performance on real benchmark datasets for ν = 0.1 𝜈 0.1 \nu=0.1

Dataset	Majority Class Size	Correct	Outliers	SVs
Iris	(50,4)	48	2	8
Heart	(139,13)	1210	18	56
Diabetes	(500,4)	466	34	90
Ecoli	(327,7)	309	18	42
Adulta3a	(2412,123)	2111	301	356
Mammography	(10923,6)	9608	1315	5969

Equations

$$\text{Total Network Power}{S}_N=\text{Active Power}{P}_N$$

$$+j\times \text{Reactive Power}{Q}_N$$

$$\sum _i{V}_{ij}{I}_{ij}^{\ast }=0,$$

$$\sum _{i\ne j,0;j\ne 0}(V_{ij})I_{ij}^{\ast }+\sum _i(V_{i0})I_{i0}^{\ast }=0$$

$$⟹\sum _{i=1}^N{V}_{i0}{I}_{i0}^{\ast }=-\sum _{i,j=1}^N{V}_{ij}{I}_{ij}^{\ast }$$

$$⟹S_T=-S_N$$

$$⟹P_T+j{Q}_T=-P_N-j{Q}_N$$

$$⟹|P_T|=|P_N|$$

$$⟹|Q_T|=|Q_N|$$

$$\underset{\{V_i,I_i\in \mathbb{C}\}}{\min }\mathcal{D}=\sum _{i=1}^N|\Re \{V_i{I}_i^{\ast }\}{|}^2$$

$$\text{s.t.}\sum _i\text{imaginary}\{V_i{I}_i^{\ast }\}=0$$

$$\sum\limits_{i=1}^{N}V_{i}\Big{(}C_{i}\dfrac{dV_{i}}{dt}\Big{)}^{*}+\Big{(}L_{i}\dfrac{dI_{i}}{dt}\Big{)}I_{i}^{*}=0,$$

$$\displaystyle\sum_{i=1}^{N}\Big{(}\dfrac{1}{2}C_{i}\lvert V_{i}\rvert^{2}+\dfrac{1}{2}L_{i}\lvert I_{i}\rvert^{2}\Big{)}=E_{0},$$

$$\underset{\{|V_i|,|I_i|\in {\mathbb{R}}_{+},{\varphi }_i\in \mathbb{R}\}}{\min }\mathcal{D}=\sum _{i=1}^N|V_i{|}^2|I_i{|}^2{\cos }^2{\varphi }_i$$

$$\displaystyle\textit{s.t.}\quad\sum_{i=1}^{N}\Big{(}\dfrac{1}{2}C_{i}\lvert V_{i}\rvert^{2}+\dfrac{1}{2}L_{i}\lvert I_{i}\rvert^{2}\Big{)}=E_{0},\enskip\lvert\phi_{i}\rvert\leq\pi\enskip\forall i$$

$$\sum\limits_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}=1,$$

$$\underset{\{|V_i|,|I_i|,{\varphi }_i\}}{\min }\mathcal{L}(\{|V_i|,|I_i|,{\varphi }_i\})=\mathcal{H}(\{|V_i|,|I_i|\})+\beta \mathcal{D}$$

$$\displaystyle\textit{s.t.}\quad\sum_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}=1,\enskip\lvert\phi_{i}\rvert\leq\pi\enskip\forall i,$$

$$\underset{\{|V_i|,|I_i|,{\varphi }_i\}}{\min }\mathcal{L}(\{|V_i|,|I_i|,{\varphi }_i\})=\underset{\{|V_i|,|I_i|\}}{\min }\mathcal{H}(\{|V_i|,|I_i|\})$$

$$\underset{\{|V|,|I|,\varphi \}}{\min }{\mathcal{L}}_1=(|V{|}^2-|I{|}^2{)}^2+\beta |V{|}^2|I{|}^2{\cos }^2\varphi$$

$$\text{s.t.}|V{|}^2+|I{|}^2=1,|\varphi |\le \pi$$

$$\underset{\{|V_i|,|I_i|,{\varphi }_i\}}{\min }\mathcal{L}(\{|V_i|,|I_i|,{\varphi }_i\})$$

$$\displaystyle\textit{s.t.}\quad\sum_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}$$

$$\frac{\partial V_i(t)}{\partial t}$$

$$\frac{\partial I_i(t)}{\partial t}$$

$$\text{and}{\tau }_i{\omega }_{{\varphi }_i}+{\varphi }_i(t)$$

$$\text{ensures that}\frac{\partial \mathcal{L}}{\partial t}\le 0,$$

$$\underset{\{|V_i|,|I_i|,{\varphi }_i\}}{\min }{\mathcal{L}}_N=\sum _{i=1}^N(|V_i{|}^2-|I_i{|}^2{)}^2+\beta \sum _{i=1}^N|V_i{|}^2|I_i{|}^2{\cos }^2{\varphi }_i$$

$$\displaystyle\textit{s.t.}\quad\sum\limits_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}=1,\quad\lvert\phi_{i}\rvert\leq\pi\enskip\forall i$$

$$\frac{\partial V_{ik}(t)}{\partial t}=j{\omega }_k{\sigma }_{V_{ik}}(t)V_{ik}(t)-\Delta {\sigma }_{V_{ik}}(t)V_{ik}(t)$$

$$\frac{\partial I_{ik}(t)}{\partial t}=j({\omega }_k+{\omega }_{{\varphi }_{ik}}){\sigma }_{I_{ik}}(t)I_{ik}(t)-\Delta {\sigma }_{I_{ik}}(t)I_{ik}(t)$$

$${\tau }_{ik}{\omega }_{{\varphi }_{ik}}+{\varphi }_{ik}(t)=g_{{\varphi }_{ik}}(t)\forall i,k,$$

$$\frac{\partial V_i(t)}{\partial t}$$

$$\frac{\partial I_i(t)}{\partial t}$$

$$\frac{\partial V_i(t)}{\partial t}$$

$$\frac{\partial I_i(t)}{\partial t}$$

$$\frac{{\partial }^2{V}_i}{\partial t^2}=j\omega {\sigma }_{V_i}\dot{V_i}+j\omega V_i\dot{{\sigma }_{V_i}}-\Delta {\sigma }_{V_i}\dot{V_i}+\dot{{\sigma }_{V_i}}{V}_i$$

$$=\underset{\text{limit cycles}}{\underbrace{-{\omega }^2{\sigma }_{V_i}^2{V}_i}}+$$

$$\underset{\text{dissipation}}{\underbrace{[1+j\omega ]{\dot{\sigma }}_{V_i}{V}_i-2j\omega {\sigma }_{V_i}\Delta {\sigma }_{V_i}{V}_i+(\Delta {\sigma }_{V_i}{)}^2{V}_i}},$$

$$\frac{{\partial }^2{I}_i}{\partial t^2}=j{\omega }_{{\varphi }_i}^{\prime }{\sigma }_{I_i}\dot{I_i}+j{\omega }_{{\varphi }_i}^{\prime }{I}_i\dot{{\sigma }_{I_i}}-\Delta {\sigma }_{I_i}\dot{I_i}+\dot{{\sigma }_{I_i}}{I}_i$$

$$=\underset{\text{limit cycles}}{\underbrace{-({\omega }_{{\varphi }_i}^{\prime }{)}^2{\sigma }_{I_i}^2{I}_i}}+$$

$$\underset{\text{dissipation}}{\underbrace{[1+j{\omega }_{{\varphi }_i}^{\prime }]{\dot{\sigma }}_{I_i}{I}_i-2j{\omega }_{{\varphi }_i}^{\prime }{\sigma }_{I_i}\Delta {\sigma }_{I_i}{I}_i+(\Delta {\sigma }_{I_i}{)}^2{I}_i}},$$

Full text

I Introduction

From an energy point of view, the dynamics of a machine learning framework is analogous to an electrical network since both the networks evolve over a conservation manifold to attain a low-energy state. In literature, this analogy has served as the basis for energy-based learning models, where the learning objective function is mapped to an equivalent network energy. [1, 2]. The network variables then evolve according to some physical principles subject to network constraints to seek out an energy optimal state.

However, even if the solution reached by such energy-based learning models may be optimal from the objective function point of view, it may not be the most energy-efficient solution when considering an actual physical implementation. This is because most of these formulations assume that the energy being minimized in the network is dissipative in nature. Whereas, in a physical electrical network, the total power $S_{N}$ (also known as the apparent power) comprises not only of the dissipative component (also referred to as the active-power) but also a latent or stored non-dissipative component (also referred to as the reactive-power) [3, 4]. This is illustrated in Figure 1(a) and can be mathematically expressed as:

[TABLE]

where $j=\sqrt{-1}$ denotes the imaginary component. While the active-power $P_{N}$ represents the rate-of-energy loss in the network, the reactive-power $Q_{N}$ represents the rate of change in stored energy in the network’s electric and magnetic fields (typically modeled as lumped capacitive and inductive elements). In the design of electrical networks, reactive-power is generally considered to be a nuisance since it represents the latent power that does not perform any useful work [5, 6]. However, from the point of view of learning, the reactive-power could be useful not only for storing the learned parameters of the network, but could also improve the dynamics of the network during learning. In this paper, we propose such a framework which exploits both active and reactive network power for learning and memory. The objective will be to achieve a network power profile as illustrated in Figures 1 (b) and (c). During learning, the network will optimize its active-power ($P_{N}$) and under steady-state condition or post-learning, ensure $P_{N}=0$. The reactive-power ($Q_{N}$), on the other hand, will always be maintained at zero. This implies that the stored energy (mathematically - the time-integral of the reactive-power) is conserved across the learning and post-learning phases respectively. Thus, during the post-learning phase or in steady-state, the network will not dissipate any power, and the reactive energy is used to maintain its current network state or memory.

This steady-state condition corresponds to a state of electrical resonance, and in Section III we generalize this concept to a framework of resonant machine learning. To reach this steady-state, in Section IV, we present a dynamical system based on complex-domain continuous-time growth-transforms, which extends our previous work on growth transform networks using real-variables[7]. The complex-domain formulation allows manipulation of the relative phase between the voltage and current variables associated with the network nodes and thus is used to optimize the active-power dissipation during learning. While the approach could be applied to different learning networks, in Section V, we use this framework for designing resonant one-class support vector machines (SVMs) [8]. In this context, we also compare the performance of the resonant SVM model with its non-resonant variant on both synthetic and real datasets. Finally, Section VI concludes the paper with a brief discussion on the implication of the proposed model when applied to other network-based models that do not involve learning, for instance, the coupled oscillator networks[9].

The key contributions of the paper can be summarized as follows:

• •

Contrary to standard energy-based models which optimize for a single energy metric, we map a learning problem to an electrical network with two energy components: dissipative(active) and non-dissipative(reactive) energy.

• •

Active power dissipation is introduced as a regularizer in the original learning problem while enforcing zero total reactive power by driving the network into electrical resonance. This ensures zero active power dissipation in the post-learning stage.

• •

We present a novel growth transform-based complex dynamical system for optimizing the cost function that ensures a resonant condition in steady-state.

• •

We propose an annealing schedule that can trade-off speed of convergence with active-power dissipation for different applications.

I-A Notations

For the rest of the paper, we will follow the mathematical notations summarized in Table 1 below:

II Background and Related works

II-A Energy-based Learning Models

The principle of minimum energy states that all self-organizing physical systems evolve dynamically over an intrinsic conservation manifold to attain their minimum energy configuration in the steady state[10]. Energy-based machine learning models follow a similar cue, where the goal is to find the optimal configuration of a predetermined energy landscape determined by the learning problem at hand[1, 11]. Some of the most notable examples of energy-based models are those based on the Ising model of statistical mechanics, like the Hopfield network[2], and its stochastic counterpart the Boltzmann machine [12] and its variants[13]. In these models, the network converges to the local minimum of a Lyapunov function in steady state. Another class of energy-based models (EBMs)[1] propose a unified learning framework by mapping the learning problem to a global scalar “energy function”. The algorithm minimizes a loss functional to seek out an optimal energy landscape that associates lowest energy to the observed configuration of variables involved in the underlying optimization process. These methods are essentially non-probabilistic, and often involve intractable integrals which require MCMC methods for normalization when applied to probabilistic frameworks. Different classes of optimization algorithms like contrastive divergence, contrastive hebbian learning, equilibrium propagation etc.[14] have been proposed for a variety of supervised, unsupervised and reinforcement learning applications[15, 16, 17, 18] in the context of energy-based models. However, all of these approaches consider an energy metric that solely depends on the learning task under consideration. This implies that in an analog implementation of the model, the optimal solution may not be the most energy-efficient one.

II-B Complex-domain Machine Learning Models

A variety of complex domain learning models have been proposed in literature for different learning algorithms, e.g., complex neural networks, complex support vector machines, complex deep networks etc[19, 20, 21, 22, 23]. In addition to providing the phase variables which allow for another degree of freedom, other advantages of the complex domain operation have been shown in these works. It has been demonstrated that complex learning models lead to a richer set of dynamics, noise robustness and better convergence properties in the case of classification problems[20, 24, 25]. Moreover, phase information might provide additional insights in the context of many complex-valued physical signals (or complex domain transforms of physical signals, e.g., Fourier or wavelet transforms)[20, 26]. However, most of these complex learning models either treat the real and imaginary components of the signal separately, or do not utilize the entire phase space when operating in the complex domain.

III Optimization and Electrical Resonance

Consider an electrical network as shown in Figure 2(a), comprising of $N$ internal nodes. The voltage difference between the $i^{\text{th}}$ and $j^{\text{th}}$ nodes is denoted by $V_{ij}$, with $V_{i0}$ being the $i^{\text{th}}$ nodal voltage with respect to ground terminal (referred to as the $0^{\text{th}}$ node). Similarly, the current flowing between the $i^{\text{th}}$ and $j^{\text{th}}$ nodes is given by $I_{ij}$, and $I_{i0}$ denotes the current flowing out of the $i^{\text{th}}$ node into the ground terminal. Then, according to the Tellegen’s theorem [27]:

[TABLE]

which states that the total complex electrical power or apparent power is zero. Isolating the apparent power flowing from the nodes to the ground terminal from that flowing between the internal nodes, we have:

[TABLE]

where $S_{T}=\sum_{i}V_{i0}I_{i0}^{*}$ is the nodal apparent power, and $P_{T}=\sum_{i}\real\{V_{i0}I_{i0}^{*}\}$ and $Q_{T}=\sum_{i}\imaginary\{V_{i0}I_{i0}^{*}\}$ are the total active and reactive power consumed at the nodes. Similarly, $S_{N}$, $P_{N}$ and $Q_{N}$ represent the apparent, active and reactive power consumed due to current flow between the network nodes (other than the ground terminal). Note that this result holds even if active-power sources are embedded in the network, as shown in Figure 2(a). Thus, Equations (3) imply that if we minimize the active-power at the nodes of the network $P_{T}$ subject to the constraint that the nodal reactive power $Q_{T}=0$, then the formulation would be equivalent to minimizing the network active-power $P_{N}$ while ensuring that the network reactive power $Q_{N}=0$. This result can be expressed as:

[TABLE]

where we have used the notations $V_{i0}=V_{i}$ and $I_{i0}=I_{i}$ for the sake of brevity. If we assume that the $i^{\text{th}}$ node is associated with a lumped capacitance $C_{i}$, and a lumped inductance $L_{i}$, ensuring zero nodal reactive-power implies

[TABLE]

where \Big{(}C_{i}\dfrac{dV_{i}}{dt}\Big{)} and \Big{(}L_{i}\dfrac{dI_{i}}{dt}\Big{)} represent the current flowing through $C_{i}$ and the voltage across $L_{i}$ respectively. Equation (5) is equivalent to:

[TABLE]

where $\lvert V_{i}\rvert^{2}=V_{i}V_{i}^{*}$ and $\lvert I_{i}\rvert^{2}=I_{i}I_{i}^{*}$, implying that the total network reactive energy is conserved to be equal to some constant value $E_{0}$. Satisfying this constraint is equivalent to sustaining a condition of electrical resonance. The optimization problem (4) can be transformed as:

[TABLE]

where $\phi_{i}$ denotes the phase-angle between the voltage and current phasors at the $i^{\text{th}}$ node. Note that the optimization in Equation (7) admits only three types of solutions in steady-state: (a) ($\lvert I_{i}\rvert\neq 0,\lvert V_{i}\rvert\neq 0,\lvert\phi_{i}\rvert=\pi/2$) which corresponds to a resonant LC tank; (b) ($\lvert I_{i}\rvert=0,\lvert V_{i}\rvert\neq 0$) which corresponds to an electrically isolated or floating node; and (c) ($\lvert I_{i}\rvert\neq 0,\lvert V_{i}\rvert=0$) which corresponds to a short-circuit. In Appendix A we illustrate steady-state resonance conditions using a simple LC tank. Note that in all cases, active power is dissipated only during the learning phase, where $C_{i}$ and $L_{i}$ adapt to change the relative magnitude of the voltage and current variables. Figure 2(b) illustrates this resonant condition of the network in the post-learning stage, whereby active power sources in the network get electrically isolated from the rest of the network. The lumped capacitance and inductance associated with the LC network at the terminal nodes adapt such that the resonant frequency condition in maintained in steady-state (see Appendix A). This implies that in steady-state, the learned network parameters are stored and sustained by the resonant electric and magnetic fields of the LC tanks.

The constraint in Equation (8) can be simplified by normalizing with respect to $E_{0}$ such that:

[TABLE]

where $V_{i}\leftarrow\sqrt{\dfrac{C_{i}}{2E_{0}}}V_{i}$ and $I_{i}\leftarrow\sqrt{\dfrac{L_{i}}{2E_{0}}}I_{i}$ represent the dimension-less voltages and currents. Henceforth, unless stated otherwise, we will use dimension-less quantities in our derivations.

We now extend the optimization framework in Equation (7) to include a general optimization function $\mathcal{H}$ as:

[TABLE]

In this formulation, the active-power dissipation $\mathcal{D}$ in Equation (III) acts as a regularization function with $\beta\geq 0$ being a hyper-parameter. Note that the objective function $\mathcal{H}(\{\lvert V_{i}\rvert,\lvert I_{i}\rvert\})$ is only a function of the magnitudes of the voltages and currents and is independent of the phase-angle $\phi_{i}$. This ensures independent control of the magnitudes and the phase to achieve the desired objective of optimizing the active-power dissipation. This is illustrated in Figure 3 where controlling the phase allows different paths from the initial to the final state, whereas evolution over the real-domain allows only one possible path. The complex-domain approach thus results in steady-state limit cycle oscillations that encode the final solution. Compared to other complex-domain machine learning frameworks [20, 23], the proposed formulation avoids non-linear processing of phase/frequency which produces unwanted higher-order harmonics. This would have made it difficult to maintain the network in the state of resonance (at specific set of frequencies) under steady-state.

The two important properties of the optimization framework in Equation (III) are as follows:

$\bullet$ For a convex cost function $H$, we have

[TABLE]

This result follows from the three possible solutions of optimization problem (7).

$\bullet$ If $\beta$ is slowly annealed to a sufficiently high value, then $\phi_{i}\rightarrow\pi/2$ under steady state for $i:\lvert V_{i}\rvert\lvert I_{i}\rvert\neq 0$. This implies that network achieves zero active power dissipation in steady state. Note that the method also holds for non-convex objective functions as well. In this case however, the network might show resonant behavior at a locally optimal solution.

Example 1: Consider a single-variable quadratic optimization problem of the form $\mathcal{H}_{1}(x)=x^{2}$, subject to the constraint $\lvert x\rvert\leq 1,\enskip x\in\mathbb{R}$. Substituting $x=\lvert V\rvert^{2}-\lvert I\rvert^{2}$, the problem can be mapped (please see Appendix B for more details) into the form equivalent to Equation (III) as:

[TABLE]

Figures 4(a)-(c) plots $\mathcal{L}_{1}$ for different values of $\beta$. As shown in 4(a), and as expected for $\beta=0$, the cost function has several minima (or attractors), whereas for $\beta>0$, the minima corresponds to $\phi=\pm\pi/2$, for which the active-power dissipation is zero. The Figures 4(b)-(c) show that controlling $\beta$ will control the optimization landscape (without changing the location of the attractors) and will determine the attractor trajectory. This feature has been exploited in later sections to optimize the active-power dissipation profile during the learning phase.

IV Complex Growth Transforms

The problem in Equation (III) involves complex phasors and hence entails the use of learning models operating in the complex domain for reaching the optimal solution. To this end, in this section, we propose a dynamical system that can be used to solve this optimization problem. The main result is summarized in Table 1 and the details of the proof is given in Appendix C.

Theorem 1: The system of nonlinear dynamical equations given by Equations (15)-(18) in Table 1 converge to the optimal point of Equation (1) in the steady state, with zero energy dissipation i.e., $\sum\limits_{n=1}^{N}\lvert V_{i}\rvert\lvert I_{i}\rvert\cos\phi_{i}=0$ (Proof given in Appendix C) .

Consider an optimization framework which is an multi-variable extension of Example 1, given by Equations  (12)-(13):

[TABLE]

The optimal solution is reached when $\phi_{i}=\pm\frac{\pi}{2}\enskip\forall i$, which implies $\mathcal{L}_{N}=0$. For the sake of comparison, we will consider two variants: (a) the non-resonant model ($\mathcal{M}_{nr}$) where $\beta=0$, and (b) the resonant model ($\mathcal{M}_{r}$), where $\beta\neq 0$. Figures 5(a) and (b) show a comparison of the cost $\mathcal{H}_{N}$ and the regularization metric $\mathcal{D}_{N}$ ($\sum\limits_{i=1}^{N}\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}$ for $\mathcal{M}_{nr}$ and $\sum\limits_{i=1}^{N}\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}\cos^{2}\phi_{i}$ for $\mathcal{M}_{r}$), for $N=5$ and $\omega=\pi/10$. In the case of $\mathcal{M}_{nr}$, $\mathcal{L}_{N}=\mathcal{H}_{N}$ and in the case of $\mathcal{M}_{r}$, $\mathcal{L}_{N}=\mathcal{H}_{N}+\beta\mathcal{D}_{N}$, with $\beta=1$. The solid line indicates the mean response across $10$ trials with random initializations of $V_{i},I_{i}$ and $\phi_{i}$. The shaded region shows the standard deviation across the trials in each case. Also shown in Figure 5(b) are the true active-power dissipation profiles ($\sum\limits_{n=1}^{N}\lvert V_{i}\rvert\lvert I_{i}\rvert$ for $\mathcal{M}_{nr}$ and $\sum\limits_{n=1}^{N}\lvert V_{i}\rvert\lvert I_{i}\rvert\cos\phi_{i}$ for $\mathcal{M}_{r}$) over the 10 trials.

It can be observed that under steady-state conditions, the model $\mathcal{M}_{nr}$ dissipates power. However for the model $\mathcal{M}_{r}$ the steady-state active-power goes to zero. This is illustrated in Figure 5(b). Figure 5(c) shows the initial phasor configuration for the currents and voltages at each node of the network for a single trial for $\mathcal{M}_{r}$. Here, we assume the voltage phasor $V_{i}$ to be aligned with the horizontal axis, and the current phasor $I_{i}$ to be at an angle $\phi_{i}$ with respect to the horizontal, $\forall i$. Figure 5(d) shows the steady-state phasor configuration for $\mathcal{M}_{r}$ for the same trial. The voltage and current phasors are orthogonal to each other for all the nodes, thus verifying the zero active-power dissipation.

In the next set of experiments, we annealed the hyperparameter $\beta$ and evaluated its impact on the active-power dissipation metric $\mathcal{D}_{N}$ and the convergence of the object function $\mathcal{H}_{N}$ for the model $\mathcal{M}_{r}$. Figure 6 presents a comparison of the performance of $\mathcal{M}_{r}$ for different choices of annealing schedule for $\beta$, with the angular frequency $\omega=\pi/10$ as before. Figure 6(a) shows the time evolution of the objective function $\mathcal{H}_{N}$, 6(b) shows the time evolution of the dissipation metric $\mathcal{D}_{N}$ and 6(c) shows the annealing schedules adopted for $\beta$. In all the cases, the optimization process starts after time $t=0.1\enskip a.u.$ from the onset of the simulation. The curves corresponding to $\beta=1$ denote the case when $\beta$ takes a constant value from $t=0.1\enskip a.u.$; $\beta=\text{logistic}$ corresponds to the case when $\beta$ is slowly increased from $\beta_{\text{min}}=0$ following a logistic curve of the form $\beta(t)=\beta_{\text{min}}+\dfrac{\beta_{\text{max}}-\beta_{\text{min}}}{1+\exp(-k(t+t_{0}))}$ from $t=0.1\enskip a.u.$, and takes on a maximum value of $\beta_{\text{max}}=1$ ($k$ and $t_{0}$ are hyperparameters determining the steepness and mid-point of the sigmoid respectively); $\beta=\text{switching}$ corresponds to the case when $\beta$ switches from a minimum value ($\beta_{\text{min}}=0$) to a maximum value ($\beta_{\text{max}}=1$) at $t=0.3\enskip a.u.$, after the system has converged to the optimal solution. We can observe that in all of the cases, the model converges to the optimal solution, irrespective of the choice of $\beta$. However, different annealing schedules for $\beta$ lead to different active-power dissipation profiles. For example, a constant value of $\beta$ throughout the duration of the experiment would lead to faster minimization of the active-power dissipation metric, but at the cost of slower convergence. The opposite trend can be seen when $\beta$ is slowly annealed to a sufficiently high value throughout the course of the optimization. The annealing schedule thus acts as a trade-off between the speed of convergence, and the rate of minimization of the active power.

IV-A Model Properties and Extensions

The dynamical system represented by Equations (15)-(18) and the resonant optimization framework also exhibits the following properties and extensions:

Energy constraints can be imposed over subgroups of nodes in the network: We can have the reactive energy conservation constraint between subgroups of nodes, instead of on the network as a whole, i.e., \sum\limits_{i=1}^{N_{k}}\Big{(}\lvert V_{ik}\rvert^{2}+\lvert I_{ik}\rvert^{2}\Big{)}=1\enskip\forall k=1,\ldots,M. Here $M=$number of subgroups and $N_{k}=$ number of nodes in the $k^{\text{th}}$ subgroup. The update equations in this case are given by:

[TABLE]

where $\omega_{k}$ is the constant angular frequency of the $k^{\text{th}}$ subgroup of nodes, and $\omega_{\phi_{ik}}=\dfrac{d\phi_{ik}(t)}{dt}$. 2. 2.

System dynamics and reactive-energy constraints remain invariant under the imposition of a global phase: The network dynamics remain invariant to the imposition of a global phase component on all the network variables, and the conservation constraint is also satisfied in this case. The governing equations are given by:

[TABLE]

where $\phi_{g}$ is the global phase and $\omega_{\phi_{g}}=\dfrac{d\phi_{g}}{dt}$. 3. 3.

Reactive-energy constraints remain invariant with varying network dynamics under the imposition of a relative phase: The conservation constraints are satisfied on the introduction of a relative phase component between the voltage and current phasors of each node, even though the overall network dynamics change. The governing equations are given by:

[TABLE]

where $\phi_{i}=\phi_{I_{i}}-\phi_{V_{i}}$ is the relative external phase shift applied between the voltage and current phasors of the $i^{\text{th}}$ node, $\omega_{\phi_{V_{i}}}=\dfrac{d\phi_{V_{i}}}{dt}$ and $\omega_{\phi_{I_{i}}}=\dfrac{d\phi_{I_{i}}}{dt}$. 4. 4.

The model is dissipative and converges to limit cycle oscillations in steady state: The second order time derivatives of Equations (15) and (16) lead to the following:

[TABLE]

where $\omega^{\prime}_{\phi_{i}}=\omega+\omega_{\phi_{i}}$. The first terms in the RHS of Equations (24) and (25) correspond to stable limit cycle oscillations of all the phasors, while the other terms correspond to the dissipative effects in the network. This demonstrates that the network as a whole is essentially a coupled dissipative system that is capable of self-sustained oscillations under steady-state. Each individual state variable describing the network thus returns to the same position in its respective limit cycle at regular intervals of time, even when subjected to small perturbations.

V Resonant Machine Learning Framework

In this section, we show how the framework introduced in Section III can be applied for constructing resonant machine learning networks. In general, the framework can be applied to any learning network that optimizes a cost function defined over a set of learning variables $\alpha_{i}$ as

[TABLE]

Here $\mathcal{H}(\{\alpha_{i}\})$ represents a loss-function [28] which depends on the learning problem (e.g. supervised, unsupervised or semi-supervised) and the dataset under consideration (e.g., training data). The second term $\Psi(\cdot)$ in the objective function is any linear or nonlinear function which represents either (a) a regularization function, or (b) a penalty function used to satisfy optimization constraints. $h$ is a hyperparameter which acts as a trade-off between $\mathcal{H}(\cdot)$ and $\Psi(\cdot)$. Because $\alpha_{i}$ could be viewed as probability measure, the optimization framework in Equation (26) naturally lends itself to probabilistic learning models [29, 30, 31].

The above problem can be mapped to the resonant learning framework in Section III by substituting $\alpha_{i}=\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}$, to arrive at the following problem:

[TABLE]

Note that non-probabilistic learning problems can also be mapped to the probabilistic framework by imposing an additional constraint, as discussed in Appendix A.

V-A One-class resonant SVM

We now show how the framework in Equation (27) can be used to design a resonant one-class SVM.

The solution of a generic one-class SVM is obtained by solving the following optimization problem [8, 32, 33]:

[TABLE]

where $\mathbf{X}=[\mathbf{x_{1}},\ldots,\mathbf{x_{i}},\ldots,\mathbf{x_{N}}]\in\mathbb{R}^{N\times D}$ is the $D-$dimensional input dataset of size $N$, $\nu\in\{0,1\}$ is a parameter which controls the size of the decision surface, $K(\cdot,\cdot)$ is a positive definite Kernel function satisfying Mercer’s conditions, and $\alpha_{i}$’s are the learning parameters. The optimization problem above can be reformulated by replacing the inequality constraint with a smooth penalty or loss function $\Psi(\cdot)$ like the logarithmic barrier, e.g., $\Psi(\alpha_{i},\nu)=-\log(\dfrac{1}{\nu N}-\alpha_{i})$.

[TABLE]

The parameter $h$ determines the steepness of the penalty function, where a lower value of $h$ implies an almost-accurate inequality constraint.

An equivalent complex-domain representation in terms of voltages and currents in an LC network can be arrived at if we consider $\alpha_{i}=\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\enskip\forall i$. In this case, we consider that the network is globally energy constrained, and all the individual units in the network have the same frequency $\omega$. The redefined learning problem is as follows:

[TABLE]

Introducing the active-power dissipation regularization we arrive at the following problem:

[TABLE]

The update equations in this case are of the form shown in Equations (15)-(18). Figure 7 shows a comparison of the active and reactive power metrics of each node of the non-resonant model $\mathcal{M}_{nr}$ and the resonant model $\mathcal{M}_{r}$ for a synthetic one-class SVM problem on a two-dimensional dataset (Dataset I). The dataset was generated by uniformly selecting $300$ random points within a circle having a fixed radius. Here $N=300,\nu=0.1,\omega=\pi/4$, with random initializations for $V_{i},I_{i},\phi_{i}\enskip\forall i=1,\ldots,N$. A constant value of the regularization hyperparameter $\beta=1$ was considered throughout the duration of the optimization process for $\mathcal{M}_{r}$. Figures 7(a) and (c) show the values of the active power metric ($=\lvert V_{i}\rvert\lvert I_{i}\rvert$) at each node in the initial and final stages respectively for $\mathcal{M}_{nr}$, while Figures 7(b) and (d) show the values of the reactive power metric ($=0$) at each node in the initial and final stages respectively for $\mathcal{M}_{nr}$. Similarly, Figures 7(e) and (g) show the values of the active power metric ($=\lvert V_{i}\rvert\lvert I_{i}\rvert\cos\phi_{i}$) at each node in the initial and final stages respectively for $\mathcal{M}_{r}$ while Figures 7(f) and (h), finally, show the values of the reactive power metric ($=\lvert V_{i}\rvert\lvert I_{i}\rvert\sin\phi_{i}$) at each node in the initial and final stages respectively for $\mathcal{M}_{r}$.

Figure 8 shows a comparison of the performance of the resonant model $\mathcal{M}_{r}$ for different choices of annealing schedules for $\beta$ on three different two-dimensional synthetic datasets for a simulation duration of $t=4\enskip a.u.$. In all cases, $N=300,\nu=0.1,\omega=\pi/20$, and $K(\cdot,\cdot)$ was chosen to be a Gaussian kernel with parameter values $\sigma=1,10$ and $20$ respectively for synthetic Dataset I, II and III. Additionally for each dataset, $V_{i},I_{i},\phi_{i}$ were randomly initialized $\forall i=1,\ldots,N$. Figure 8(a) shows the contour plot with the decision boundary around the datapoints, along with the support vectors (SVs) for Dataset I, while Figures 8(b), (c) and (d) show the time evolutions of the cost $\mathcal{H}$, the dissipation metric $\mathcal{D}$ and the hyperparameter $\beta$ for different annealing schedules. The curves corresponding to $\beta=1$ and $10$ denote the cases when $\beta$ takes a constant value throughout the simulation duration. $\beta=\text{logistic}_{1}$ and $\beta=\text{logistic}_{2}$ correspond to the cases when $\beta$ is slowly increased following a logistic curve of the form $\beta(t)=\beta_{\text{min}}+\dfrac{\beta_{\text{max}}-\beta_{\text{min}}}{1+\exp(-k(t+t_{0}))}$, and takes on maximum values of $\beta_{\text{max}}=1$ and $\beta_{\text{max}}=10$ respectively starting from $\beta_{\text{min}}=0$. Finally, $\beta=\text{switching}$ corresponds to the case when $\beta$ switches from a minimum value ($\beta_{\text{min}}=0$) to a maximum value ($\beta_{\text{max}}=10$) at $t=2\enskip a.u.$, after the system has converged to the optimal solution. Figures 8(e)-(h) show similar plots on Dataset II, while Figures 8(i)-(l) show the plots corresponding to Dataset III. Dataset I, as described before, was generated by selecting $300$ points uniformly within a two-dimensional circle having a fixed radius. Datasets II and III also consist of $300$ data points generated using a Gaussian mixture model consisting of four different clusters, with a fixed cluster mean and variance associated with each cluster. Dataset II and III differ only in terms of the cluster means associated with their respective constituent clusters, while the cluster variances are the same for all the clusters for both the datasets. It can be seen that since the optimization problem is convex, the model always converges to the optimal solution for every dataset, irrespective of the the annealing schedule for $\beta$ or the dataset complexity. However, the dissipation profiles corresponding to a particular annealing schedule strongly depend on the complexity of the dataset. In general, however, higher values of $\beta$ would lead to lower dissipation profiles during the learning process. Also, the model shows a much slower convergence in terms of the actual objective function for a constant non-zero value of $\beta$ throughout the optimization compared to the case when $\beta$ is slowly annealed to a sufficiently high value. The choice of a proper annealing schedule for $\beta$ would thus involve a trade-off between the speed of convergence and the rate of power dissipation.

Figures 9(a) and (b) finally, show the robustness of the proposed framework to different initial conditions by providing a comparison of the time evolutions of the cost $\mathcal{H}$ and the dissipation metric $\mathcal{D}$ ($\sum\limits_{n=1}^{N}\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}$ for the non-resonant model $\mathcal{M}_{nr}$ and $\sum\limits_{n=1}^{N}\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}\cos^{2}\phi_{i}$ for the resonant model $\mathcal{M}_{r}$ respectively), when applied to Dataset I. We show the results for 10 random initializations of $V_{i},I_{i}$ and $\phi_{i},\enskip\forall i$. In all cases, $\nu=0.1,\omega=\pi/8$, and $K(\cdot,\cdot)$ was chosen to be a Gaussian kernel with parameter value $\sigma=1$. Note that even though the 10 trails had different initializations of $V_{i},I_{i}$, they were chosen to have the same initial value of $\alpha_{i}^{\prime}$s in all cases, because of which there is no deviation between the cost evolution curves for both $\mathcal{M}_{nr}$ and $\mathcal{M}_{r}$. The dissipation evolution is however different across different trials for both the models. However, the dissipation attains a final value of zero for $\mathcal{M}_{r}$ for all the trials, while there is a finite dissipation for $\mathcal{M}_{nr}$ in all cases.

We also conducted experiments on real-life benchmark datasets of varying sizes and dimensionality, considering the majority class as inliers in all cases. In all the experiments, we used similar parameter settings and annealing schedules as used in the experiments shown in Figure 8, and $\nu=0.1$. Table 2 shows the dataset description along with the performance (number of inliers, number of outliers and number of support vectors) of the one-class SVM classifier for the following datasets: ‘Iris’ (inlier: class ‘setosa’), ‘Heart’ (inlier: ‘healthy’ heart), ‘Diabetes’ (inlier: non-diabetic patients), ‘Ecoli’ (inlier: classes ‘cp’, ‘im’, ‘pp’, ‘imU’ and ‘om’), ‘Adulta3a’ (inlier: individuals with income\leq\50K$) and ‘Mammography’ (inlier: non-calcifications)[34, 35]. Figures 10 (a)-(f) show the dissipation profiles corresponding to different annealing schedules for the different datasets. It can be seen that an optimal choice of the annealing schedule depends on both the dataset size and complexity, even though the dissipation decreases over time for all the cases, irrespective of the schedule chosen.

Interestingly, the solution of the one-class SVM can be interpreted in terms of an electrical network as follows, as shown in Figure 11: (a) the support vectors have voltage and current magnitudes with a $\pm\pi/2$ phase shift between them, hence can be interpreted as resonant LC tanks; and (b) the interior points well inside the boundary have both zero voltage and current magnitudes, and can be essentially treated as floating sinks.

VI Discussions and Conclusions

In this paper we proposed a complex-domain formulation of a machine learning network that ensures that the network’s active power is dissipated only during the process of learning whereas the network’s reactive power is maintained to be zero at all times. We showed that the active power dissipation during learning can be controlled using a phase regularization parameter. Also, the framework is robust to variations in the initial conditions and to the choice of the input/driving frequency $\omega$. The proposed approach thus provides a physical interpretation of machine learning algorithms, where the energy required for storing learned parameters is sustained by electrical resonances due to nodal inductances and nodal capacitances. Using the one-class support vector machine problem as a case study, we have shown how the steady state solution of a learning problem can be interpreted in terms of these nodal circuit elements.

Future directions involve exploring the implications of incorporating network dynamics in the free frequency variable $\omega$, and utilizing the phase information associated with each node in the learning process. Also, the experimental results presented in this paper were based on an unsupervised learning setting. Exploring an energy-efficient framework for supervised learning problems would also be a part of our future work.

In this paper we also proposed a dynamical system model based on complex-domain growth transforms. The formulation is general enough to be applicable to other complex-domain learning models[19, 20, 21, 22, 23]. Our proposed framework also preserves both the magnitude and phase information, and provides additional flexibilty compared to other complex domain learning models in terms of phase manipulation/cancellation [20].

In addition to implementing classic machine learning algorithms, the complex growth transform dynamical system can also be used for designing synchronized networks of coupled oscillators. Such networks can be potentially used for solving different computing tasks like optimization, pattern matching etc. as is achievable using coupled oscillatory computing models[36, 37]. An oscillator network designed in this fashion is capable of demonstrating stable, self-sustaining oscillatory dynamics, whereby the network can return to its initial stable limit cycle configuration following small perturbations, while simultaneously minimizing some underlying system-level objective function. The framework could also be used to study connections between learning and synchronization, or the emergence of a rhythmic periodic pattern exhibited by a group of coupled oscillators, which provides the key to understanding periodic processes pervading complex networks of different biological, physical, social and quantum ensembles [9, 38]. In this regard, the existing mathematical models for such collective behavior are mostly phenomenological or bottom-up, and in general do not provide a network-level perspective of the underlying physical process. The proposed growth-transform formulation, thus, could provide new network-level insights into the emergence of phase synchronization, phase cohesiveness and frequency synchronization in coupled-oscillator networks.

Note here that since we are implicitly assuming an analog implementation, the learning network should converge to the steady-state solution of the optimization, where the time of convergence depends on the network’s time-constant [39]. Hence the notion of time complexity is not well-defined. However, in a digital implementation, the time complexity would depend on the learning algorithm under consideration, dataset size and dimensionality, angular frequency $\omega$, and time constants $\tau_{i}$ for the phase updates.

VII Appendix

VII-A Resonance in an LC tank

Consider the parallel LC tank circuit shown in Figure 12, with $V_{C}$ and $V_{L}$ being the voltages across the capacitor $C$ and inductor $L$ respectively. $I_{C}$ and $I_{L}$ denote the corresponding currents flowing through the elements. Thus, $V_{S}=V_{L}=V_{C}$ and $I_{S}=I_{L}+I_{C}$. Considering the LC tank to be driven by the voltage source $V_{S}$ at frequency $\omega$, we have the following condition in steady state:

[TABLE]

Resonant condition of the circuit is achieved when

[TABLE]

This implies that the apparent power $S_{N}=P_{N}+jQ_{N}=V_{S}I_{S}^{*}+V_{L}I_{L}^{*}+V_{C}I_{C}^{*}$, where the active power $P_{N}=0$. Additionally at resonance, the reactive power $Q_{N}=Q_{C}+Q_{L}=V_{L}I_{L}^{*}+V_{C}I_{C}^{*}=-\dfrac{j}{\omega L}\lvert V(\omega)\rvert^{2}+\dfrac{j}{\omega L}\lvert V(\omega)\rvert^{2}=0$. Here $Q_{C}$ and $Q_{L}$ are the reactive powers associated with the capacitance and inductance respectively.

VII-B Mapping a generic optimization problem to the equivalent network model

Let us consider an optimization problem defined over a probabilistic domain, given by the following generic form:

[TABLE]

We can map the above to the electrical network-based model introduced in Section III by replacing $x_{i}=\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}$, which leads to the following problem in the $\{\lvert V_{i}\rvert^{2},\lvert I_{i}\rvert^{2}\}$ domain:

[TABLE]

Note that the method also works for optimization problems defined over non-probabilistic domains, of the following form:

[TABLE]

This can be done by considering $x_{i}=x_{i}^{+}-x_{i}^{-}\enskip\forall i$, where both $x_{i}^{+},x_{i}^{-}\geq 0$. Since by triangle inequality, $\lvert x_{i}\rvert\leq\lvert x_{i}^{+}\rvert+\lvert x_{i}^{-}\rvert$, enforcing $x_{i}^{+}+x_{i}^{-}=1\enskip\forall i$ would automatically ensure $\lvert x_{i}\rvert\leq 1\enskip\forall i$. Thus we have,

[TABLE]

In this case, we replace $x_{i}^{+}=\lvert V_{i}\rvert^{2},\enskip x_{i}^{-}=\lvert I_{i}\rvert^{2}$, and the equivalent problem in the $\{\lvert V_{i}\rvert^{2},\lvert I_{i}\rvert^{2}\}$ domain is thus as follows:

[TABLE]

For example, the variables $\{x_{i}\}$ can represent the Lagrangian multipliers in the primal space of a support vector machine network, or the weights and biases of a generic neural network.

VII-C Complex growth transform dynamical system (Table 1)

Let us consider the optimization problem in Equation (34) again. We can use the Baum-Eagon inequality [40, 41] to converge to the optimal point of $\mathcal{H}$ in steady state, by using updates of the form

[TABLE]

Here, $\mathcal{H}$ is assumed to be Lipschitz continuous [7] on the domain $D=\{x_{1},\ldots,x_{N}:\sum\limits_{i=1}^{N}x_{i}=1,\enskip x_{i}\geq 0\enskip\forall i\}\subset\mathbb{R}^{N}_{+}$. The constant $\lambda\in\mathbb{R}_{+}$ is chosen such that \Big{\lvert}-\dfrac{\partial\mathcal{H}(\{x_{i}\})}{\partial x_{i}}+\lambda\Big{\rvert}>0,\forall i.

We can solve the optimization problem given by Equation (III) by using the growth transforms discussed above. The outline of the proof is as follows: (1) We will start with a generic magnitude domain optimization problem without any phase regularizer and derive the form for the growth transform dynamical system which would converge to the optimal point asymptotically; (2) Next, we derive a complex domain counterpart of the above, again without phase constraints; (3) Finally, we derive the complex domain dynamical system by incorporating a phase regularizer in the objective function.

Since the time evolutions of $V_{i}$ and $I_{i}$ are symmetric because of the conservation constraints, for the rest of the section we will consider only the update equations for the voltages $V_{i}$, and similar results would also apply to the updates for $I_{i}$.

VII-C1 Condition 1

Considering $\beta=0$ in Equation (III) and $\mathcal{H}(\{\lvert V_{i}\rvert,\lvert I_{i}\rvert\})$ to be Lipschitz continuous over the domain D=\{\lvert V_{i}\rvert,\lvert I_{i}\rvert:\sum\limits_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}=1\}, we can use the growth transforms to arrive at the following discrete-time update equations in terms of the voltage and current magnitudes :

[TABLE]

where

[TABLE]

and $\lambda\in\mathbb{R}_{+}$ is chosen to ensure that \Big{(}-\dfrac{\partial\mathcal{H}}{\partial\lvert V_{i}\rvert^{2}}+\lambda\Big{)}>0 and \Big{(}-\dfrac{\partial\mathcal{H}}{\partial\lvert I_{i}\rvert^{2}}+\lambda\Big{)}>0,\forall i.

Writing $g_{V_{i},n-1}=g_{V_{i},n-1}(\{\lvert V_{i,n-1}\rvert^{2},\lvert I_{i,n-1}\rvert^{2}\})$ for notational convenience and taking $g_{V_{i},n-1}=\lvert V_{i,n-1}\rvert^{2}\sigma_{V_{i},n-1}^{2}$, we get:

[TABLE]

VII-C2 Condition 2

Considering $\beta=0$ in Equation (III) and V_{i},I_{i}\in D^{C}=\{V_{i},I_{i}\in\mathbb{C}:\sum\limits_{i=1}^{N}\Big{(}\lvert V_{i}\rvert^{2}+\lvert I_{i}\rvert^{2}\Big{)}=1\}, a time evolution of the form given below converges to limit cycles corresponding to the optimal point of a Lipschitz continuous objective function $\mathcal{H}(\{\lvert V_{i}\rvert,\lvert I_{i}\rvert\})$:

[TABLE]

where $\sigma_{V_{i}},\sigma_{I_{i}}\rightarrow 1\quad\forall i=1,\ldots,N$, in steady state, $\theta$ is the instantaneous global phase difference of the system of phasors with respect to an absolute stationary reference frame, and $\phi_{i}$ is the instantaneous phase difference between $V_{i}$ and $I_{i}$.

Proof:

Since $\lvert V_{i,n}\rvert^{2}=V_{i,n}V_{i,n}^{*}$, $\lvert I_{i,n}\rvert^{2}=I_{i,n}I_{i,n}^{*}$, the update equations can be rewritten as:

[TABLE]

where $\sigma_{{V_{i}},n-1}$ is used to represent $\sigma_{V_{i},n-1}(\{V_{i,n}V_{i,n}^{*},I_{i,n}I_{i,n}^{*}\})$ for ease of notation, and similarly for $\sigma_{{I_{i}},n-1}$. Considering $\mathcal{H}(\{V_{i}V_{i}^{*},I_{i}I_{i}^{*}\})$ to be analytic in $D^{C}$ and applying Wirtinger’s calculus [42], since

[TABLE]

we have

[TABLE]

where

[TABLE]

The discrete time update equations for $V_{i,n}$ is thus given by:

[TABLE]

Similar expressions can also be derived for the current phasors.

$\blacksquare$

Lemma: A continuous-time variant of the model given by Equation (VII-C2) is given by:

[TABLE]

Proof: The difference equations for the voltage phasors are given by:

[TABLE]

which can be rewritten as:

[TABLE]

where $\Delta\sigma_{V_{i},n-1}=1-\sigma_{V_{i},n-1}$, and $\theta=\omega\Delta t$, where $\omega$ is the common oscillation frequency of all the phasors, and $\Delta t$ is the time increment. In the limiting case, when $\Delta t\rightarrow 0$, this reduces to the following continuous time differential equation for the complex variable $V_{i}(t)$:

[TABLE]

$\blacksquare$

In the steady state, since $H$ is Lipschitz continuous in $D$, $\sigma_{V_{i}}(t)\overset{t\rightarrow\infty}{\longrightarrow}1$, the dynamical system above thus reduces to:

[TABLE]

which implies that the steady state response of the complex variable $V_{i}(t)$ corresponds to steady oscillations with a constant angular frequency of $\omega$.

The difference equations in terms of the nodal currents can be similarly written as:

[TABLE]

where $\omega_{\phi_{i}}=\dfrac{d\phi_{i}}{dt}$. The equivalent continuous domain differential equation is then given by:

[TABLE]

$\blacksquare$

VII-C3 Condition 3

Considering $\beta\neq 0$ in Equation (III), additional phase constraints can be incorporated in the dynamical system represented by using the update rules in Equations (15)-(18). In steady state, for $\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}\neq 0$, the system settles to $\phi_{i}=\pm\pi/2\enskip\forall i$. Additionally, for sufficiently small values of $\beta$ (or if $\beta$ is slowly increased during the optimization process), the system converges to the optimal point of the modified objective function $\mathcal{H}(\{\lvert V_{i}\rvert,\lvert I_{i}\rvert\})$.

Proof: Since $\mathcal{L}(\{\lvert V_{i}\rvert,\lvert I_{i}\rvert,\phi_{i}\})$ is Lipschitz continuous in both $\lvert V_{i}\rvert^{2}$ and $\lvert I_{i}\rvert^{2}$, the same form of the update equations proved in Lemma 2 can be applied. For arriving at the updates for the phase angles $\phi_{i}$, we will use a similar approach as shown in Equation (37). We can split $\phi_{i}$ as $\phi_{i}=\phi_{i}^{+}-\phi_{i}^{-},\enskip\phi_{i}^{+},\phi_{i}^{-}>0$, which implies that $\phi_{i}^{+}+\phi_{i}^{-}=\pi$. We can then apply the growth transform dynamical system [7] to get

[TABLE]

where $\omega_{\phi_{i}}=\dfrac{d\phi_{i}(t)}{dt}$ and

[TABLE]

Since $\dfrac{\partial\mathcal{L}}{\partial\phi_{i}}=\dfrac{\partial\mathcal{L}}{\partial\phi_{i}^{+}}=-\dfrac{\partial\mathcal{L}}{\partial\phi_{i}^{-}}$, the above can be simplified to $g_{\phi_{i}}=\pi\dfrac{\lambda\phi_{i}-\pi\dfrac{\partial\mathcal{L}}{\partial\phi_{i}}}{-\phi_{i}\dfrac{\partial\mathcal{L}}{\partial\phi_{i}}+\lambda\pi}$ This implies that the voltage and current phasors corresponding to the $i-$th node in the network may be phase shifted by an additional amount $\phi_{i}$ with respect to the absolute reference.

Since for optimality, $\phi_{i}=\pm\pi/2$ for $\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}\neq 0$ in the steady state, the net energy dissipation in steady state is zero, i.e., $\sum\limits_{n=1}^{N}\lvert V_{i}\rvert^{2}\lvert I_{i}\rvert^{2}\cos^{2}\phi_{i}=0$. Also, in the steady state, for sufficiently small values of the hyperparameter $\beta$,

[TABLE]

which implies that the system reaches the optimal solution with zero active power in the post-learning stage.

$\blacksquare$