Amoeba-inspired analog electronic computing system integrating resistance crossbar for solving the travelling salesman problem

Combinatorial optimization to search for the best solution across a vast number of legal candidates requires the development of a domain-specific computing architecture that can exploit the computational power of physical processes, as conventional general-purpose computers are not powerful enough. Recently, Ising machines that execute quantum annealing or related mechanisms for rapid search have attracted attention. These machines, however, are hard to map application problems into their architecture, and often converge even at an illegal candidate. Here, we demonstrate an analogue electronic computing system for solving the travelling salesman problem, which mimics efficient foraging behaviour of an amoeboid organism by the spontaneous dynamics of an electric current in its core and enables a high problem-mapping flexibility and resilience using a resistance crossbar circuit. The system has high application potential, as it can determine a high-quality legal solution in a time that grows proportionally to the problem size without suffering from the weaknesses of Ising machines.


Amoeba-based computing system
(a) depicts the basic scheme of an amoeba-based computing system, called the "bounceback control" 38 . The amoeba-based computing system is constructed with two components. The first is a single-celled amoeboid organism, a plasmodium of the true slime mould Physarum polycephalum (hereafter "amoeba" for simplicity), placed inside a stellate container that has several lanes (grooves), each of which enables the amoeba to expand and shrink its pseudopod-like branch to indicate a value of each corresponding state variable. The second is a controller that monitors the shape of the amoeba, governs the interactions between the variables in accordance with certain rules and irradiates visible light on a number of branches of the amoeba to induce their shrinking behaviour.
The amoeba on an agar plate optimizes its body shape to maximize its body area for maximal nutrient absorption while minimizing the risk of being illuminated by light, where the volume of the amoeba in the container is almost conserved. When the amoeba expands its branch in lane so that the proportion of the area covered by the branch exceeds a threshold, state variable takes a value of 1, otherwise 0. The controller always monitors the amoeba branches and determines which branches are to be flipped from 1 to 0 in accordance with certain rules. Once the controller illuminates a branch to be flipped, the branch gradually shrinks due to its photoavoidance response. This system reaches a legal solution when the light illumination pattern and the shape of the amoeba become stable. Note that the expanding and shrinking behaviour of each amoeba branch to the light fluctuates stochastically; with a small probability, the branch expands even under illumination and hesitates to expand even under the nonilluminated condition. The fluctuation is found to play an important role in exploring a wider state space 33,36,37 .
In the amoeba-based computing system for solving -city TSP, the variable ∈ The bounceback signal is defined as follows: where , ,

⁄
. Equations (S1) and (S2) determine whether the lane is illuminated or not. For example, when 0.5 , the branch is illuminated and the pseudopod shrinks, otherwise the branch is not illuminated, and the pseudopod expands.
Parameters and in equation (S2) correspond to the inverse temperature in statistical mechanics and threshold value, respectively. In the amoeba-based computing system 38 , the parameters were set as follows: 1000, 0.5, 35, and 0.6. It is noted that any other parameter optimization is not necessary 38 . In the case of the electronic amoeba, we use a step function instead of the sigmoid function , in the crossbar IMC, since the sigmoid function is the almost same as a step function when is large enough.
The single-celled organism in the amoeba-based computing system attempts to find a shape where the pseudopods in branches fully elongate along the non-illuminated lanes by accumulating the history of the illumination. In accordance with the policy of the bounceback control, the amoeba shape representing the shortest route, that is the optimal solution, corresponds to the most stable state. At this state the amoeba body area is maximized, where nutrient absorption is maximized and the risk of illumination is minimized 35 .

Amoeba-inspired electronic computing system
Figures 2(a) and S1(b) shows a schematic of the electronic computing system for solving the travelling salesman problem (TSP), called an electronic amoeba 32,33 , which we designed on the basis of a mathematical model of the solution-searching dynamics of the amoeba-based computing system, AmoebaTSP, formulated by Aono et al 35,38 . In the amoeba core, the current flowing into unit represents the expanding behaviour of the branch, and state variable take In this study, although we did not implement any external fluctuation source on the system intentionally, we found that the system intrinsically generates a fluctuated behaviour in its dynamic operation.
Figures S2(a) and S2(b) show the circuit of the pseudopod-like branch and the experimentally measured response, respectively. When the MOSFET was turned off, corresponding to 0, the capacitor charged the current from the current source and the output of the sigmoid-like function circuit gradually approached low level. On the other hand, when it was turned on, corresponding to 1, the capacitor discharged and the output immediately became high level.
The configuration in the IMC for solving the TSP appears to be similar to that of the Hopfield's recurrent neural network 12,13 . However, the interactions between the variables in the electronic amoeba are derived from a different formulation; in the bounceback control, the interaction is defined to prohibit state variables from violating the constraints of the TSP to exclude revisiting a once-visited city and to exclude simultaneous visits to multiple cities 38 .
Aono et al. formulated another algorithm for solving the Boolean satisfiability problem (SAT), named AmoebaeSAT 36,37 , and implemented it on analogue and digital electronic circuits to maximize their rapid search performance 38,39 . These electronic amoebae for solving the TSP and SAT can easily handle arbitrary problem instances without any costly pre-processing, and their dynamics stabilize only at legal solutions that satisfy all given constraints. As electronic amoebae can be composed of conventional complementary metal-oxide semiconductor (CMOS) devices, they are highly scalable, energy-efficient, and suitable for both cloud-and edge-computing applications.

Figure S1
(a) Basic scheme of the amoeba-based computing system. The amoeboid organism expands its pseudopod-like branches. However, the branch shrinks when illuminated by light. The solution-searching of the amoeba-based computing system progresses with the changes in the amoeba shape and illumination pattern. (b) Circuit diagram of the electronic amoeba for solving the TSP.  S3(d). This is because the variation effect of the resistance increases relatively when the current is increased and/or the capacitance is decreased. Therefore, we can avoid the degradation of the solution quality by reducing the resistance variation.

Figure S3
Dependence

Amoeba-GA hybrid approach
As a possible approach to improve the quality of the solution found by the electronic amoeba, we have investigated the dynamic adjustment of the resistance variation in the amoeba core by introducing a genetic algorithm (GA). A method that incorporates a GA for dynamic adjustment of the parameters to the base algorithm is called a "hybrid-GA" or "memetic algorithm S1-S3 ." This approach has been reported to be effective for Ising machines S4 and is expected to be beneficial for our purpose as well, although there is an overhead for the communication between the electronic amoeba and a conventional computer. Figure S4 shows a flowchart of the electronic amoeba-GA hybrid approach. The system has many electronic amoebae and assigns random resistance values to all or some of the units. The randomly assigned resistance values are regarded as individuals in the GA. Then, the electronic amoebae start the solution search simultaneously. When we obtain several legal solutions, their qualities are evaluated. The combination of the resistance values in the amoeba core that finds a higher-quality solution is selected to take over to the next generation. To promote the diversity in each generation, we apply the "crossover" to swap the resistance value between the selected combinations and introduce the "mutation" to change the values randomly. The combinations of the updated resistance values are fed into the units of the electronic amoebae, and the parallel solution search start again. Repeating these operations, it is expected that the qualities of the solutions found by the electronic amoebae would be improved collectively.

Figure S4
Flowchart of Amoeba-GA hybrid approach.

Fluctuation-imposing approach
To improve the solution quality, the mathematical model of the amoeba-based computing system, AmoebaTSP 38 , introduced stochastic noise to fluctuate the state variable temporally, which resulted in effective trial-and-error dynamics for exploring a wider state space. In the case of the electronic amoeba, if the time scale of the dynamics of the amoeba core could be set to be larger than that of the temporal fluctuation, similar trial-and-error dynamics would be feasible and effective. We can verify this hypothesis by numerical simulation. We formulate the solution-searching dynamics of the electronic amoeba as the following time evolution equations: where ℎ 1 is the number of the units without applying the bounceback signals (i.e., an analogy of the number of amoeba branches at which the lights are turned off) at an iteration time step 1, Δ and Δ correspond to the expanding and shrinking rates (velocities) of the branch, respectively, and is random noise to represent stochastic fluctuation. The numerical simulation is conducted by iterating the calculations of equations (S3) and (S4) successively.
Examples of the time series obtained from the simulation for solving a 10-city TSP instance is shown in Figure S5, where a representative dynamics of state variable , and the outputs of the sigmoid function  respectively. The sigmoid curve fluctuates largely when large noise is imposed. We stress here that when averaging , the dynamic range of the nonlinear output is extended by adding large noise, which is sometimes explained in the context of stochastic resonance S5,S6 . We consider that the fluctuation-imposing approach presented above enables the electronic amoeba to explore a larger variety of candidate solutions, which has to be verified experimentally in the future.

Figure S5
Examples of output waveforms when solving the 10-city TSP instance at Δ

Solution-searching performance of the Ising models
The Ising model searches for a spin assignment that minimizes the Hamiltonian of energy.
The Hamiltonian of the Ising model is defined with a cost function and penalty term S7 as follows: , S5 where and are the cost function and penalty term, respectively, and is a coefficient to tune the weight of the penalty term. The cost function is responsible for the quality of a solution and takes a lower value when the length of a route selected is shorter. The penalty term is to ensure the legality of the solution and is designed to increase the energy when the system goes into an illegal state violating the TSP constraints. Therefore, it is a very important to optimize to obtain a high-quality legal solution 29 . When is relatively large, the Ising model tends to find a legal solution that is close to that obtained by random sampling due to the minimal impact of the cost function on the Hamiltonian. On the other hand, when is too small, the system tends to reach an illegal state, which represents an infeasible route. Hence, to determine the optimal value, a user of the Ising model has to perform many solution-searching runs until it succeeds in finding a legal solution. For this reason, D-Wave's Ising machine 17 requires to perform the post-processing for repairing the converged spin alignment to satisfy the constraints and for improving the solution quality 2,S8 . The similar treatments are also necessary for the Hopfield's neural network to obtain a legal solution. In the case of the electronic amoeba, however, these costly pre-processing and post-processing are not necessary.
We evaluate the solution-searching performance of the Ising model-based algorithms for comparison with that of the electronic amoeba. The Ising model for solving the TSP is formulated on the basis of quadratic unconstrained binary optimization (QUBO) given as follows S9 : where , 1 denotes that a salesman visits city at the th order, is the distance between cities and and is a parameter to tune the weight of the first two penalty terms.
An assignment of the variables at the minimum value of represents an optimal solution of the TSP. Equation (S6) can be converted to the Hamiltonian of the Ising model as , , 1 2 ⁄ , , ∈ 1,1 , where , represents the spin at lattice site , . Figures S8(a) and S8(b) show flowcharts of the two algorithms for exploring a variable assignment that minimizes in equation (S6), which are the Ising-model-based algorithms without and with simulated annealing (SA), respectively (referred to as "Ising model" and "Ising model SA"). Figure S9 shows the simulation results of the Ising model. For each , when of the Ising model was set at 100 that resulted in finding a legal solution whose quality was comparable to that found by the electronic amoeba ( Fig. S9(a)), the rate of simulation trials in which the Ising model found a legal solution out of 100 trials was only around 20% ( Fig. S9(b)).
Note that the Ising model succeeds in converging at a legal solution by setting large enough, such as 0 max . For each , when was set at max 1 that ensured to find a legal solutions with a success rate of 100% ( Fig. S9(d)), the quality of a legal solution found by the Ising model degraded and could not be comparable to that found by the electronic amoeba ( Fig. S9(c)). On the other hand, a success rate of the electronic amoeba to find a legal solution was always 100%. Figure S10 shows results of the Ising model SA. When 100, the solution quality of the Ising model SA was further improved by reducing the temperature gradient and became higher than that of the electronic amoeba ( Fig. S10(a)). However, with the increase in the number of cities , the success rate of finding a legal solution decreased as shown in Fig. S10(b).
When max 1, although the Ising model SA found a legal solution with a success rate of 100% (Fig. S10(d)), the solution quality degraded considerably to a level that was even worse than that of the electronic amoeba (Fig. S10(c)). These results confirm that the electronic amoeba is more useful than the Ising machines for practical applications that require to obtain a legal solution reliably and swiftly.

Figure S8
Flowcharts of Ising model-based algorithms for TSP (a) without SA and (b) with SA.