Physics Successfully Implements Lagrange Multiplier Optimization

The principle of Least Action is the most fundamental principle in physics. Newton’s Laws of Mechanics, Maxwell’s Equations of Electromagnetism, Schrödinger’s equation in Quantum Mechanics, and Quantum Field Theory can all be interpreted as minimizing a quantity called action. For the special case of light propagation, this reduces to the principle of Least Time, as shown in Fig. 1.

A conservative system without friction or losses evolves according to the principle of Least Action. The fundamental equations of physics are reversible. A consequence of this reversibility is the Liouville Theorem which states that volumes in phase space are left unchanged as the system evolves.

Contrary-wise, in both a computer and an optimization solver, the goal is to have a specific solution with less uncertainty or a smaller zone in phase space than the initial state, an entropy cost first specified by Landauer and Bennett. Thus, some degree of irreversibility, or energy cost, is needed, specified by the number of digits in the answer in the Landauer/Bennett analysis. An algorithm has to be designed and programmed into the reversible system (either via software or via hardcoding of hardware) to effect the reduction in entropy needed to solve the optimization problem. Coming up with a fast algorithm for NP-hard problems is still an open problem in the field of reversible computing, which includes quantum computing.

This would require an energy cost, but not necessarily a requirement for continuous power dissipation. We look forward to computer science breakthroughs that would allow the principle of Least Action to address unsolved problems. An alternative approach to computing would involve physical systems that continuously dissipate power, aiding in the contraction of phase space toward a final solution. While exceeding the Landauer limit, such systems might have the advantage of speed and simplicity. This brings us to the principle of Least Power dissipation.

If we consider systems that continuously dissipate power, we are led to a second optimization principle in physics, the principle of Least Entropy Production or Least Power Dissipation. This principle states that any physical system will evolve into a steady-state configuration that minimizes the rate of entropy production given the constraints (such as fixed thermodynamic forces, voltage sources, or input power) that are imposed on the system. An early version of this statement is provided by Onsager in his celebrated papers on the reciprocal relations [2]. This was followed by further foundational work on this principle by Prigogine, [21], and de Groot, [22]. This principle is readily seen in action in electrical circuits and is illustrated in Fig. 2. We shall frequently use this principle, as formulated by Onsager, in the rest of the paper.

This technique is widely used in materials science and metallurgy and involves the slow cooling of a system starting from a high-temperature. As the cooling proceeds, the system tries to maintain thermodynamic equilibrium by reorganizing itself into the lowest energy minimum in its phase space. Energy fluctuations due to finite temperatures help the system escape from local optima. This procedure leads to global optima when the temperature reaches zero in theory but the temperature has to be lowered prohibitively slowly for this to happen.

The Adiabatic Method involves the slow transformation of a system from initial conditions that are easily constructed to final conditions that capture the difficult problem at hand.

More specifically, to solve the Ising problem using this algorithm, one initializes the system of spins in the ground state of a simple Hamiltonian and then slowly varies the parameters of the Hamiltonian to end up with the final Ising Hamiltonian of interest. If the parameters are varied slowly enough, the system ends up in the ground state of the final Hamiltonian and the problem gets solved. In a quantum mechanical system, this is sometimes called ‘quantum annealing’. Several proposals and demonstrations, including the well-known D-Wave machine [23], utilize this algorithm.

The slow rate of variation of the Hamiltonian parameters that is required for this method to work is determined by the minimum energy spacing between the instantaneous ground state and the instantaneous first excited state that occurs as we move from the initial Hamiltonian to the final Hamiltonian. The smaller the gap is for a particular variation schedule, the slower the rate at which we need to perform the variation to successfully solve the problem. It has been shown that the gap can get exponentially small in the worst case, implying that this algorithm can take exponential time in the worst case for NP-hard problems.

Multi-Oscillator Arrays, subject to Parametric Gain were introduced in [4] and [7] for solving Ising problems. This can be regarded as a subset of the Principle of Minimum Entropy Generation, which is always subject to a non-zero input power constraint. In this case, gain acts as a boundary condition or constraint for the principle of minimum entropy generation, and the oscillator array must arrange itself to dissipate the least power subject to that constraint. In the context of a multi-coupled-oscillator arrays with gain, a certain oscillator mode will have the least loss. That mode will grow in amplitude most rapidly. This least loss mode can be individually selected by increasing the gain until it matches the intrinsic loss of that mode. Then further nonlinear evolution amongst all the modes will occur. If the oscillator array is bistable, as is the case for parametric gain which drives oscillation along the $\pm$real axis, this becomes the analog of magnetic bistability in an Ising problem. Then, we seek a solution for the lowest magnetic energy state with some oscillators locked at zero phase shift and others locked at $\pi$-phase shift. This mechanism will be the main point of Section 3.

An early optical machine for solving the Ising problem was presented by Yamamoto et al. [4] and [8]. Their system consists of several pulses of light circulating in an optical fiber loop, with the phase of each light pulse representing an Ising spin. In parametric oscillators, gain occurs at half the pump frequency. If the gain overcomes the intrinsic losses of the fiber, the optical pulse builds up. Parametric amplification provides phase dependent gain. It restricts the oscillatory phase to the Real Axis of the complex plane. This leads to bistability along the positive or negative real axis, allowing the optical pulses to mimic the bistability of magnets.

In the Ising problem, there is magnetic coupling between spins. The corresponding coupling between optical pulses is achieved by controlled optical interactions between the pulses. In Yamamoto’s approach, one pulse $i$ is first plucked out by an optical gate, amplitude modulated by the proper connection weight specified in the $J_{ij}$ Ising Hamiltonian, and then reinjected and superposed onto the other optical pulse $j$, producing constructive or destructive interference, representing ferromagnetic or anti-ferromagnetic coupling.

By providing saturation to the pulse amplitudes, the optical pulses will finally settle down, each to one of the two bistable states. We will find that the pulse amplitude configuration evolves exactly according to the Principle of Minimum Entropy Generation. If the magnetic dipole solutions in the Ising problem are constrained to $\pm 1$ then each constraint is associated with a Lagrange Multiplier. We find that each Lagrange Multiplier turns out to be equal to the gain or loss associated with the corresponding oscillator.

As an example, Yamamoto et al. [7] analyze their parametric oscillator system using coupled wave equations for the slowly varying amplitudes of the circulating optical modes. We now show that the coupled wave equation approach reduces to an extremum of their system ‘Entropy Generation’ or ‘power dissipation’. The coupled-wave equation for parametric gain of the slowly varying amplitude $c_{i}$ of the in-phase cosine component of the $i$-th optical pulse (representing magnetic spin in an Ising system) is as follows:

where the weights, $J_{ij}$, are the magnetic ferromagnetic or anti-ferromagnetic cross-couplings, and $\gamma_{i}$ represents the phase dependent parametric gain given by the pump to the $i$-th circulating pulse, whereas $\alpha_{i}$ is the corresponding dissipative loss. The amplitudes $c_{i}$ represent the optical cosine wave electric field amplitudes. In the case of circuits, the voltage amplitudes can be expressed as $V_{i}=c_{i}+\text{j}s_{i}$ where $c_{i}$ and $s_{i}$ represent the cosine and sine quadrature components of voltage, and j is the unit imaginary. For clarity of discussion, we dropped the cubic terms in (1) that Yamamoto et al originally had. A discussion of these terms in given in Appendix C.

Owing to the nature of parametric amplification, the quadrature sine wave components $s_{i}$ of the electric field amplitude die out rapidly. The rate of entropy generation, or net power dissipation $h$, including the negative dissipation associated with gain can be written:

If we minimize the entropy generation $h(\bm{c})$ without invoking any constraints, that is, with $\gamma_{i}=0$, the amplitudes $c_{i}$ simply go to zero, which generates the minimum entropy.

If the gain $\gamma_{i}$ is large enough, some of the amplitudes might go to infinity. To avoid this, we may employ the $n$ constraint functions $g_{i}(c_{i})=\left(1-c_{i}^{2}\right)=0$, which enforce a digital $c_{i}=\pm 1$ outcome. (Actually, a constraint of the form $g_{i}(c_{i})=\left(1-c_{i}^{2}\right)=0$ is quite general in the sense that $\left(1-c_{i}^{2}\right)$ can represent the first two terms of the Taylor Series of an arbitrary constraint.) Adding the constraint function to the entropy generation, yields the Lagrange function including the constraint functions times the respective Lagrange Multipliers:

Comparing the unconstrained (2) to the constrained (3), they only differ in the final (-1) term which effectively constrains the amplitudes and prevents them from diverging to $\infty$. Expression (3) is the Lagrange function given at the end of Section 4. Surprisingly, the gains $\gamma_{i}$ emerge to play the role of Lagrange Multipliers. This means that each mode, represented by the subscripts in $c_{i}$, must adjust to a particular gain $\gamma_{i}$ which minimizes the overall entropy generation, and the respective gains $\gamma_{i}$ represent the Lagrange Multipliers. Minimization of the Lagrange function (3) provides the final steady state of the system dynamics.

If the circuit or optical system is designed to dissipate power, or equivalently generate entropy, in a mathematical form that matches the magnetic energy in the Ising problem, then the dissipative system will seek out a corresponding local optimum configuration of the magnetic Ising energy.

Such a physical system, constrained to $c_{i}=\pm 1$, is digital in the same sense as a flip-flop circuit, but unlike the von Neumann computer, the inputs are resistor weights for power dissipation. Nonetheless a physical system can evolve in a direct manner, without the need for shuttling information back and forth as in a von Neumann computer, providing faster answers. Without the communications overhead but with the higher operation speed, the energy dissipated to arrive at the final answer will be less, in spite of the circuit being required to generate entropy during its evolution toward the final state.

To achieve minimum entropy production, the amplitudes $c_{i}$, and the Lagrange Multipliers $\gamma_{i}$, must all be simultaneously optimized. While a circuit will evolve toward optimal amplitudes $c_{i}$, the gains $\gamma_{i}$ must arise from a separate active circuit. Ideally, the active circuit which controls the Lagrange Multiplier gains $\gamma_{i}$, would have its entropy production included with the main circuit. A more common method is to provide gain that follows a heuristic rule using an external feedback circuit. For example, Yamamoto et al. follow the heuristic rule $\gamma_{i}=a+bt$. It is not yet clear whether the heuristic-based approach toward gain evolution will be equally effective as simply lumping together all main circuit and feedback components and simply minimizing the total power dissipation.

We conclude this subsection by noting that the Lagrangian, Eq (3), corresponds to Lagrange multiplier optimization using the following merit function and constraints:

1. 1.

   Physical systems minimize the entropy production rate, or power dissipation, subject to input constraints of voltage, amplitude, gain, etc.

2. 2.

   These systems actually perform Lagrange Multiplier optimization.

3. 3.

   Indeed, it is the gain $\gamma_{i}$ in each oscillator $i$ that plays the role of the corresponding Lagrange Multiplier.

4. 4.

   If the Lagrange function is split in such a way that only the Ising function, $\sum_{i,j}J_{ij}c_{i}c_{j}$, is treated as the merit function $f(\bm{c})$, then the Lagrange multipliers corresponding to the constraints $g_{i}$ are the net gains, $\gamma_{i}-\alpha_{i}$.

5. 5.

   Under the digital constraint, amplitudes $c_{i}=\pm 1$, entropy generation minimization schemes are actually binary, similar to a flip-flop.

A coupled LC oscillator system with parametric amplification was analyzed in the circuit simulator, SPICE, by Xiao et al., [12]. This is analogous to the optical Yamamoto system but this system consists of a network of radio frequency LC oscillators coupled to one another through resistive connections. The LC oscillators contain linear inductors but nonlinear capacitors which provide the parametric gain. The parallel or cross-connect resistive connections between the oscillators are designed to implement the ferromagnetic or anti-ferromagnetic couplings $J_{ij}$ between magnetic dipole moments $\mu_{i}$ as shown in Fig. 7. The corresponding phase of the voltage amplitude $V_{i}$, 0 or $\pi$ determines the sign of magnetic dipole moment $\mu_{i}$.

The nonlinear capacitors are pumped by voltage $V(2\omega_{0})$ at frequency $2\omega_{0}$, where the LC oscillator natural frequency is $\omega_{0}$. Second harmonic pumping leads to parametric amplification in the oscillators. As in the optical case, parametric amplification plays the dual role of generating/sustaining voltage oscillations as well as imposing phase bistability to the ac voltages in the oscillators. The bistability refers to gain along the $\pm$Real-Axis defined by time synchronization with the $2\omega_{0}$-pump.

The $2\omega_{0}$-pump induces gain $\gamma_{i}$ in the Real-Axis quadrature. As in the case of optical parametric amplifier machines, ideally, an active circuit would control the Lagrange Multiplier gains $\gamma_{i}$, and the gain control circuit would have its entropy production included with the main circuit. A more common approach is to provide gain that follows a heuristic rule.

As in the optical parametric amplifier case, a mechanism is needed to prevent the parametric gain from producing infinite amplitude signals. Zener diodes can be inserted into the circuit to restrict the amplitudes to finite saturation values or to digitally defined values. With the diodes in place, the circuit settles into a voltage phase configuration, 0 or $\pi$, that minimizes net power dissipation for a given pump gain.

The amplitudes can be defined by the voltages across the LC oscillator capacitors and derived from Kirchoff’s voltage and current equations as in Xiao, [12]. Performing the slowly-varying amplitude approximation on the cosine component of these voltages, $c_{i}$, Xiao obtains the following equation of motion:

where the $c_{i}$ are the peak voltage amplitudes in units of a reference voltage, $R_{c}$ is the resistance of the coupling resistors, the cross-couplings $J_{ij}$ are assigned binary values $J_{ij}=\pm 1$, $C_{0}$ is the linear part of the capacitance in each oscillator, $n$ is the number of oscillators, $\omega_{0}$ is the natural frequency of all the oscillators, and the parametric gain constant $\gamma=\omega_{0}|\Delta C|/4C_{0}$ where $|\Delta C|$ is the capacitance modulation at the second harmonic. In the decay constant $\alpha=(n-1)/(4R_{c}C_{0})$, there are $(n-1)$ resistors $R_{c}$ in parallel, since it is assumed that each oscillator can give up energy to all the other $(n-1)$ oscillators, with the coupling resistors acting in parallel. In this simplified model all decay constants $\alpha$ are taken as equal, and moreover each oscillator experiences exactly the same parametric gain $\gamma$; conditions that can be relaxed if needed.

The number 4 is present in the first denominator, since for two coupled LC circuits as shown in Fig. 7, the decay is controlled by the RC time of the capacitors in parallel and the $R_{c}$ resistors in series. Likewise, the number 4 appears for loss and parametric gain, but for different reasons in each case.

We note that equation (4) above performs gradient descent on the net power dissipation function:

which is very similar to Section 5.1. On the right-hand side, the reason for the 4 in the first denominator is to compensate for double counting in $i\neq j$. The first two terms on the right-hand side together represent the dissipative losses in the coupling resistors while the third term is the negative of the gain provided to the system of oscillators.

Next, we obtain the following Lagrange function through the same replacement of $\left(-c_{i}^{2}\right)$ with $\left(1-c_{i}^{2}\right)$ that we performed in Section 5.1:

The above Lagrangian corresponds to Lagrange multiplier optimization using the following merit function and constraints:

Again, we see that the gain coefficient $\gamma$ is the Lagrange multiplier of the constraint $g=0$.

Although the extremum of Eq (6) represents the final evolved state of a physical system representing an optimization outcome, it would be interesting to examine the time evolution toward the optimal state. The optimization occurs by iterative steps, where each iteration can be regarded to take place in a time interval $\Delta t$. At each successive iteration, the voltage amplitude $c_{i}$ takes a step whose magnitude is proportional to the gradient of the Lagrange function:

where the minus sign on the right hand side drives the system toward minimum power dissipation. As the Lagrange function comes closer to its minimum, the gradient $\frac{\partial}{\partial c_{i}}L(\bm{c},\gamma)$ diminishes and the amplitude steps become smaller and smaller. The adjustable proportionality constant $\kappa$, controls the size of each iterative step; it also calibrates the dimensional units between power dissipation and voltage amplitude. (Since $c_{i}$ is voltage amplitude, $\kappa$ has units of reciprocal capacitance.) Converting Eq (7) to continuous time,

where the $\gamma_{j}$ play the role of Lagrange multipliers, and the $g_{j}=0$ are the constraints. Taking $L(\bm{c},\gamma)$ from Eq (6), the gradient in the voltage amplitudes becomes

Substituting into Eq (8), the time derivative of the voltage amplitudes becomes

The constant $\kappa$ can be absorbed into the units of time to reproduce the dynamical equation (4), the slowly varying amplitude approximation for the coupled radio oscillators. Thus, it is interesting that the slowly-varying time dynamics can be reproduced from iterative optimization steps on the Lagrange function.

1. 1.

   As in other cases, the coupled LC oscillator system [12] minimizes the entropy production rate, or power dissipation, incorporating the power input from the pump gain.

2. 2.

   The coupled LC oscillator system actually performs Lagrange Multiplier optimization whose Merit Function is the Lagrange Function.

3. 3.

   The gain $\gamma_{i}$ in each oscillator $i$ plays the role of the corresponding Lagrange Multiplier.

4. 4.

   Under the amplitude constraint on the voltage cosine component, $c_{i}=\pm 1$, the entropy generation minimization scheme is actually binary, similar to a flip-flop.

5. 5.

   Successive iterations toward a power dissipation minimum, employing gradient descent in the amplitudes of the cosine components of the voltages $c_{i}$ , surprisingly, yield the same time-dependent Slowly Varying amplitude differential equation (4).

The Ising solver designed by Babaeian et al., [13], makes use of coupled laser modes in a multicore optical fiber. Polarized light in each core of the optical fiber corresponds to each magnetic moment in the Ising problem. This means that the number of cores needs to be equal to the number of magnets in the given Ising instance. The right-hand circular polarization and left-hand circular polarization of the laser light in each core represent the two polarities (up and down) of the corresponding magnet. The mutual coherence of the various cores is maintained by injecting seed light from a master laser.

The coupling between the fiber cores is achieved through amplitude mixing of the laser modes by Spatial Light Modulators at one end of the multicore fiber, [13]. These Spatial Light Modulators couple light amplitude from the $i$-th core to the $j$-th core according to the prescribed connection weight $J_{ij}$.

As in prior physical examples the electric field amplitudes can be expressed in the slowly-varying polarization modes of the $i$-th core, $E_{iL}$ and $E_{iR}$, where the two electric field amplitudes are in-phase temporally, are positive real, but have different polarization. They are,

where $\alpha_{i}$ is the decay rate in the optical core in 1/sec units, and $\gamma_{i}$ is the gain supplied to the $i$-th core. The first term on the right in both the equations represents optical fiber losses while the second term represents the gain provided. The third term represents the coupling between the $j$-th and $i$-th cores that is provided by the Spatial Light Modulators. We next define the degree of polarization as $\mu_{i}\equiv E_{iL}-E_{iR}$. (This differs from the usual definition of degree of polarization of an optical beam which contains intensity rather than electric field amplitude.) Subtracting the two equations above, we obtain the following evolution equation for $\mu_{i}$:

The power dissipation, or entropy generation, is proportional to two orthogonal components of electric field, each squared, $|E_{iL}|^{2}+|E_{iR}|^{2}$. But this can also be written $|E_{iL}-E_{iR}|^{2}+|E_{iL}+E_{iR}|^{2}=|\mu_{i}|^{2}+|E_{iL}+E_{iR}|^{2}$. But $|E_{iL}+E_{iR}|^{2}$ can be regarded as relatively constant as energy switches back and forth between right and left circular polarization. Then the changes in power dissipation, or entropy generation, $h(\bm{\mu})$ would be most influenced by quadratic terms in $\bm{\mu}$:

As before, we add the $n$ digital constraints $g_{i}(\mu_{i})=1-\mu_{i}^{2}=0$ where $\mu_{i}=\pm 1$ represents fully left or right circular optical polarization. With the constraints, the corresponding Lagrange function is:

Once again, the gains $\gamma_{i}$ play the role of Lagrange multipliers. Thus, a minimization of the power dissipation, subject to the optical gain $\gamma_{i}$, solves the Ising problem defined by the same $J_{ij}$ couplings.

The power dissipation function and constraint function in the Lagrange function above are:

1. 1.

   The coupled multicore fiber system [13] minimizes the entropy production rate, or power dissipation, which includes the power input from the pump gain, as in Sections 5.1 and 5.2.

2. 2.

   The coupled multicore fiber system actually performs Lagrange Multiplier optimization.

3. 3.

   The gain $\gamma_{i}$ in each fiber $i$ plays the role of the corresponding Lagrange Multiplier.

4. 4.

   Under the constraint of completely light polarization, $\mu_{i}=\pm 1$, the entropy generation minimization procedure is actually binary, similar to a flip-flop.

In this section, we consider a network of interconnected, nonlinear, but amplitude-stable electrical oscillators designed by Roychowdhury et al. [14] to represent a conventional Ising system for which we seek a digital solution with each magnetic dipole $\mu_{iz}=\pm 1$ along the z-axis in the magnetic dipole space. To solve this, Roychowdhury et al. provide a dissipative system of LC oscillators, somewhat similar to the Optical Parametric Oscillators in Section 5.1, but with oscillation amplitude clamped, and oscillation phase $\phi_{i}=0$ or $\pi$ revealing the preferred magnetic dipole orientation $\mu_{iz}=\pm 1$. It is noteworthy that Roychowdhury goes beyond Ising machines and constructs general digital logic gates using these amplitude-stable oscillators in [28].

In the construction of their oscillators, Roychowdhury et al. [14] use nonlinear elements that behave like negative resistors at low voltage amplitudes but with saturating resistance at high voltage amplitudes. This produces of amplitude-stable oscillators. In addition, Roychowdhury et al. [14] provide a second harmonic pump and use a form of parametric amplification (referred to as sub-harmonic injection locking in [14]) to obtain bistability with respect to phase.

The dynamics of the amplitude saturation is purposely very fast and the oscillators are essentially clamped to an amplitude limit dictated by the nonlinear resistor in each oscillator. The phase dynamics induced by the 2^nd harmonic pump are slower. It is the readout of these phase shifts, 0 or $\pi$ that provides the magnetic dipole orientation $\mu_{iz}=\pm 1$. One key difference between this system and Yamamoto’s Optical Parametric Oscillator (OPO) system is that the OPO system had fast phase dynamics and slow amplitude dynamics, while the injection locking system has the reverse.

Roychowdhury et al. [14] derived the dynamics of their amplitude stable oscillator network using perturbation concepts developed in [29]. While a circuit diagram is not shown, [14] invokes the following dynamical equation for the phases of their electrical oscillators:

where $R_{c}$ is a coupling resistance in their system, $\phi_{i}$ is the phase of the $i$-th electrical oscillator, and the $\lambda_{i}$ are decay parameters that dictate how fast the phase angles settle towards their steady state values.

We shall now show that Eq (11) can be reproduced by iteratively minimizing the power dissipation in their system. Power dissipation across a resistor is $(V_{1}-V_{2})^{2}/R_{c}$ where $(V_{1}-V_{2})$ is the voltage difference across resistor $R_{c}$. Since $V_{1}$ and $V_{2}$ are sinusoidal, the power dissipation consists of constant terms and a cross-term of the form:

where $f(\phi_{1},\phi_{2})$ is the power dissipated in the resistors. Magnetic dipole orientation parallel or anti-parallel is represented by whether electrical oscillator phase difference is $\phi_{1}-\phi_{2}=0$ or $\phi_{1}-\phi_{2}=\pi$ respectively. We are permitted to choose an origin for angle space at $\phi=0$ which implies $\phi_{i}=0$ or $\phi_{i}=\pi$. This can be implemented through a constraint of the form:

Combining the power dissipated in the resistors with the constraint function $g_{i}(\phi_{i})=0$ we obtain a Lagrange function:

where $\lambda_{i}$ is the Lagrange Multiplier corresponding to the phase angle constraint, and $J_{ij}$ is a digital multiplier $\pm 1$ to the conductance $\frac{1}{R_{c}}$.

The Lagrange function above is isomorphic with the general form in Section 4. The effective merit function $f$ and constraints $g_{i}$ in this correspondence are:

1. 1.

   As in other cases, the amplitude-stable oscillator system [14] minimizes the entropy production rate, or power dissipation, subject to constraints of amplitude and phase reference.

2. 2.

   The amplitude-stable oscillator system actually performs Lagrange Multiplier optimization.

3. 3.

   The phase decay time constant $\lambda_{i}$ in each oscillator $i$ plays the role of the corresponding Lagrange Multiplier.

4. 4.

   Under the phase reference constraint, $\phi_{i}=0$ or $\pi$, the entropy generation minimization scheme is actually binary, similar to a flip-flop.

Kalinin and Berloff [15] proposed a system consisting of coupled polariton condensates to minimize the XY Hamiltonian. The XY Hamiltonian is a 2 dimensional version of the Ising Hamiltonian and is given by:

where the $\bm{\mu_{i}}$ represents the magnetic moment vector of the $i$-th spin restricted to the spin-space XY plane.

To represent the spin system, Kalinin and Berloff pump a grid of coupled semiconductor microcavities with laser beams and observe the formation of strongly coupled exciton-photon states called polaritons. For our purposes, the polaritonic nomenclature is irrelevant. For us, these are simply coupled electromagnetic cavities similar to optical and LC resonators that we have already discussed. The electromagnetic system operates by the principle of minimum entropy generation similar to the previous cases. The complex electromagnetic amplitude in the $i$-th microcavity can be written $E_{i}=c_{i}+\text{j}s_{i}$, where $c_{i}$ and $s_{i}$ represent the cosine and sine quadrature components of $E_{i}$, and j is the unit imaginary. In this case we identify $\text{Re}(E)=c$ as representing the X-component of the magnetic dipole vector, and $\text{Im}(E)=s$ representing the Y-component of the magnetic dipole vector. The electromagnetic microcavity system settles into a state of minimum entropy generation as the laser pump and optical gain are ramped up to compensate for the intrinsic cavity losses. The phase angles in the complex plane of the final electromagnetic modes are then reported as the corresponding $\bm{\mu}$-magnetic moment angles in the XY plane.

An important point to note here is that the electromagnetic cavities experience normal phase-independent gain and not parametric gain which happens to be phase-dependent. As a consequence, this system does NOT seek phase bistability as appropriate for binary up/down spin orientations. This is because we are searching for the magnetic dipole vector angles in the XY plane which would minimize the corresponding XY magnetic energy.

Ref. [15] uses ‘Ginzburg-Landau’ equations to analyze their system resulting in equations for the complex amplitudes $\Psi_{i}$ of the polariton wavefunctions. But the $\Psi_{i}$ are actually the complex electric field amplitudes $E_{i}$ of the corresponding $i$-th cavity. The electric field amplitudes satisfy the following slowly-varying-amplitude equation:

where $\gamma_{i}$ represents the optical gain, $\beta$ represents nonlinear attenuation, $U$ represents nonlinear phase shift, and $J_{ij}$ is a dissipative cross-coupling-term representing linear loss. We note from the above that both the amplitudes and phases of the electromagnetic modes are coupled to each other and evolve on comparable timescales. This is in contrast to ref. [14] where the main dynamics were embedded in phase—amplitude was fast and almost fixed—or conversely [12] embedded in amplitude—phase was fast and almost fixed.

We shall now show that the method of ref. [15] is essentially the method of Lagrange multipliers with an added ‘rotation’. The entropy generation or power dissipation rate is:

If we add a saturation constraint, $g_{i}(E_{i})=\nobreak\left(1-|E_{i}|^{2}\right)=\nobreak 0$, then by analogy to the previous sections, $\gamma_{i}$ is reinterpreted as a Lagrange Multiplier:

where $L$ is the Lagrange function that represents power dissipation combined with the amplitude constraint $|E_{i}|^{2}=1$. Thus, the scheme of coupled polaritonic resonators operates to find the state of minimum entropy generation in steady-state, similar to the parametric oscillator case of Section 5.1. That is, the steady-state solution of dynamics (12) renders (13) stationary with respect to changes in $\bm{E}$. The difference is that the coupled polaritonic system solves the XY Ising problem for a magnetic moment restricted to the magnetic XY plane, while the parametric oscillator system, in Section 5.1, solves the $\pm$Z Ising problem. For any particular mathematical optimization that we would perform, we still retain the burden of specifying that dissipative system whose entropy generation matches the optimization that we are seeking.

The imaginary ‘rotation’ term, $iU$, could potentially be of use in developing more sophisticated algorithms than the method of Lagrange multipliers and we discuss this prospect in some detail in Section 6.2 where a system with a similar, but more general, ‘rotation’ term is discussed.

Next, we discuss how Kalinin and Berloff et al. [15], adjusted their Lagrange multipliers $\gamma_{i}$ during the course of their optimization. In the method of Lagrange multipliers, the merit-function, eq (13), is used to optimize not only the electric field amplitudes $E_{i}$ but also the Lagrange Multipliers $\gamma_{i}$. The authors in Sections 5.1 and 5.2 used simple heuristics to adjust their gains/decay constants which we have proven to be the Lagrange multipliers. Kalinin and Berloff of this section employ the Lagrange function optimization itself to adjust the gain/losses, as in the complete Lagrange method discussed next.

To iteratively adjust the Lagrange multipliers, we shall briefly shift back to the tutorial notation of Section 4. Afterward, we shall translate it back to the notation of this section and show that Kalinin, Berloff et al., indeed use this same iterative procedure to adjust their Lagrange multipliers.

The Lagrange function is given by $L(\bm{x},\bm{\lambda})=f(\bm{x})+\sum_{i=1}^{p}{\lambda_{i}g_{i}(\bm{x})}$, where $x_{i}$ are the field variables and $\lambda_{i}$ are the Lagrange multipliers. Then, the procedure to find the optimal $\bm{x^{*}}$ and $\bm{\lambda^{*}}$ is to perform gradient descent of $L$ in $\bm{x}$ and gradient ascent of $L$ in $\bm{\lambda}$. The reason for ascent in $\bm{\lambda}$ rather than descent is to more strictly penalize deviations from the constraint. In the language of iterations, this leads to the expressions:

where $\kappa$ and $\kappa^{\prime}$ are suitably chosen step sizes.

With our identification that the Lagrange multipliers $\lambda$ are actually the same as the gains $\gamma$, we plug the expression for the full Lagrange function (13) into the second iterative equation, projecting out the constraint function $g_{i}(E_{i})=\left(1-|E_{i}|^{2}\right)=0$, and take the limit $\Delta t\to 0$. We obtain the following dynamical equation for the gains $\gamma_{i}$:

The above equation for the iterative evolution of the Lagrange multipliers is indeed the very same evolution that Kalinin, Berloff et al. employ in their coupled polariton system.

To Eq (14), we must add the iterative evolution of the field variables $x_{i}$:

Equations (14) and (15) represent the full iterative evolution, but in some of the earlier sub-sections, $\gamma_{i}(t)$ was sometimes assigned a heuristic time dependence.

We conclude this sub-section by splitting the Lagrange function into the effective merit function $f$, and the constraint function $g_{i}$. The extra ‘phase rotation’ U is not captured by this interpretation.

1. 1.

   As in other cases, the coupled-cavity system [15] minimizes the steady state entropy production rate, or power dissipation, which includes the power input from the pump gain.

2. 2.

   The gain $\gamma_{i}$ in each oscillator $i$ plays the role of the corresponding Lagrange Multiplier.

3. 3.

   There is a phase rotation term $iU$ in the dynamics that is not captured by the Lagrange function framework.