A Circular Formation Method for Biomimetic Robotic Fish Inspired by Fish Milling

Circular motion phenomena, akin to fish milling, are prevalent within the animal kingdom. This paper delineates two fundamental mechanisms underlying such occurrences: forward following and circular topological communication. Leveraging these pivotal concepts, we present a multi-agent formation circular model based on a second-order integrator. This model engenders the attainment of homogeneous intelligence convergence along the circumferential trajectory. The convergence characteristics are intricately linked to the number of agents and the model parameters. Consequently, we propose positive and negative solutions for ascertaining the convergent circle property and model parameters. Furthermore, by integrating our proposed formation control methodology with a robotic fish dynamics model, we have successfully implemented simulations and experiments, demonstrating the circular formation of multiple biomimetic robotic fish. This study provides a mathematical explication for the circular motion observed in animal groups and introduces a novel approach to achieving circular formation in multiple robots inspired by biological phenomena.


Introduction
An intertwined and causal progression relationship exists between biological collectives, collective models, and multi-robot formation.Researchers draw inspiration from biological collective phenomena, establish corresponding group models, and apply them to multi-robot formation tasks.Among these, the step from natural collective phenomena to artificial group models is the most crucial.During this process, group models based on simple rules serve as a bridge to explore the mechanisms of animal collective motion.It is generally believed that the closer the simulation results of group models align with the observations of reality, the closer these simple rules that constitute the models come to the essence of biological collective motion.
Since Rynalods introduced the Boids rules (cohesion, separation, and alignment) in 1987, many studies in the field of animal collective behavior have considered "alignment", which refers to the consistency of movement direction, as a critical factor in the formation of group behavior [1].As a result, numerous related models have emerged [2][3][4][5].Among them, the most influential is the group model proposed by Vicsek in 1995, whose model considers self-propelled particles moving on a plane and coordinates neighbors' interactions within a certain distance based on the alignment rule to achieve overall consistency motion [2].
In the research field of control theory, starting from the Vicsek model, Morse et al. used non-negative matrices and stability theory to discuss the consensus problem of discrete systems with time-varying topology in 2003 [6].Subsequently, in 2004, Murray et al. researched the consensus problem of switching topology and time-delay systems [7].The studies by Morse and Murray ignited a research trend in the field of control theory for multi-agent systems, which has been widely applied to multi-robot formation tasks.This has led to the development of centralized leader-follower approaches, virtual structure methods, and various distributed methods.With the continuous research on consensus theory, various orderly formation methods for multi-robot formation have been developed [8][9][10][11].
It is essential to note that the Vicsek model based on the Boids rules only reproduces the polarized phenomenon in fish schools.In addition, typical ordered collective motion in fish schools also includes phenomena such as milling and bait balls.Some studies have explored models that do not adhere to the Boids rules and have achieved various group patterns.Romanczuk et al. developed a biologically motivated model based on only pursuit and escape interactions, which achieved spatial migration and vortex-like structures [12].Similarly, Strömbom et al. proposed a model based solely on mutual attraction, which formed structures resembling milling and chain rotations by limiting the individual's field of view [13].Barberis proposed a group model based on position attraction, introducing a conical field of view, which resulted in linear structures, rod-like structures, milling-like structures, and gas-like structures by changing the opening angle of the cone [14].Robert et al. presented a model with short-range velocity matching and long-range anti-alignment rules.This model generated directed ordered states, periodic vortex patterns, and medium-scale turbulence, resembling observations of swimming bacteria in dense suspensions [15].Bastein et al. proposed a general collective model based on visual projection, which accounted for visual occlusion without requiring explicit velocity matching.By varying model parameters, this model exhibited behaviors such as aggregation, polarization, milling, and swarming under different conditions [16].In our previous research, we proposed a fellow-following principle, established a collective model, and quantitatively compared it with the real fish school, demonstrating the similarity between the model and the real fish school [17].
However, most studies mentioned above have yet to extend the models into multirobot formation applications.This paper aims to apply the models established in previous works to multi-robot formation, not only to achieve bio-inspired multi-robot formation but also to delve deeper into the mechanism of fish school's milling behavior.In fact, as depicted in Figure 1, milling, or the so-called circular group motion, is not exclusive to fish.It is prevalent among various animal groups, including ants [18], chickens, ducks [19], cows, and sheep, among other animals [20].Therefore, further investigation of the general mechanism behind such milling-like motion holds significance both in enhancing our understanding of animal collective behavior and inspiring robotic swarm formation.On the one hand, in terms of animal behavior, while the previously proposed model accomplishes the group milling, it does not explicitly elaborate on its formation mechanism or explain why milling structures take on a circular shape [17].This is due to the self-organizing nature of the model's process, resulting in a certain level of randomness in the outcomes.Establishing a direct causal relationship between the proposed rules and the milling motion is challenging.Thus, in pursuit of a deeper understanding of animal collective mechanisms, proposing a group model with a concise mathematical representation becomes necessary [21].On the other hand, within the field of biomimetics, entirely self-organized group models can not be directly applicable to robot formation tasks, especially when dealing with a small robot group.
Considering these two points, the proposed model was analyzed and abstracted, leading to the development of the multi-robot circular formation method presented in this paper.The main contributions of this paper can be summarized as follows: 1.
Inspired by fish milling and our previous model, it is assumed that the formation of milling is the cause of the circular communication topology and forward-following rule.Based on the proposed hypothesis, the first-order and second-order kinematics models are established, and the convergence characteristics are analyzed.

2.
Based on the second-order kinematics model, combined with the robotic fish dynamics model, a robotic fish circular formation controller is proposed, and the simulation and multi-robotic fish circular formation experiment are carried out.
The remainder of this paper is organized as follows.The kinematics model and its analysis are presented in Section 2. In Section 3, the circular formation method based on the dynamics model of robotic fish is established and verified by simulation.Several groups of circular formation experiments were carried out based on robotic fish in Section 4, followed by a discussion in Section 5.

Kinematics Models
Inspired by the phenomenon of fish school and the group model proposed in our previous research, this section presents a group control protocol designed for kinematics models in a circular communication topology context.In the communication topology, where agents form a circular chain, each agent communicates with its neighboring agents in the front (or front and back) and adheres to the following rule.Through simulation, this protocol achieves a circular formation in which all agents are evenly distributed along the circumference.Protocol parameters determine the convergence radius of the circle, and the center position is related to the system's initial state.Agents only communicate with their adjacent neighbors, rendering the proposed control protocol highly scalable.It can form a circular formation with as few as three agents to an infinite number of agents.The formation approach holds significant inspirational value for circular formation tasks involving robots.Interestingly, upon further literature review, we discovered a rule-based circular formation approach that, while slightly different in its mathematical formulation compared to our approach, shares a similar underlying concept.These studies will be mentioned in the subsequent stability analysis.

First-Order Kinematics Model
Considering the scenario in which agents follow a first-order kinematics model, define N all as a set of n agents.This group of agents exhibits a circular communication topology.In a system composed of five agents, as shown in Figure 2, the agents are sequentially connected, forming a closed loop and creating a circular undirected graph.Information propagates between adjacent agents.The communication topology of the agents is represented by a connectivity matrix A: In the kinematics model, the agents are treated as particles without collisions and orientations.Their state is represented by x, and, in this paper, the agents move within a two-dimensional plane, with the state x being a point in that plane.The state update equation for the agents is depicted in Formula (2).
where x i t represents the position of agent i at time t, u i t indicates the control input for agent i at time t, and ζ is the control protocol function that depends on the positions of agents i − 1, i, and i + 1 at time t.For the first-order model, the control protocol ζ(x i−1 t , x i t , x i+1 t ) for the i-th agent is defined as follows: where δ i f is defined as the forward-following vector of agent i, δ i b represents the backwardrepulsion vector of agent i, and δ i is the weighted sum of these two vectors.k f indicates the weight for the forward-following vector, and δ i * is the normalized unit direction vector of δ i .The motion step size for each update of the agent is a constant v c .

Second-Order Kinematics Model
Considering the case in which agents follow a second-order motion model, the set of agents N all exhibits a circular communication topology.In a system composed of five agents, as illustrated in Figure 3, the agents are connected sequentially, forming a closed loop and creating a circular directed graph.Information is transmitted from the forward agents to the adjacent agents behind them.The communication topology of the agent set can be represented by a connectivity matrix A: In terms of the second-order case, consider a second-order unicycle model for the agents' motion.The agents move within a two-dimensional plane, and their state is represented by x and v, where x is the position and v is the velocity.The state update equation for the agents is as follows: where x i t and v i t are the position and velocity direction of agent i at time t.v i * represents the unit vector of v i .u i t denotes the control input for agent i at time t.v c indicates the motion step size for the agent.ζ is the control protocol for the second-order motion model, which depends on the position and velocity of agent i at time t and the position of agent i+1.
For the second-order model, the control protocol ζ for the i-th agent is defined as: where δ i f represents the forward-following vector of agent i. δ i indicates the weighted sum of δ i f and the normalized direction vector v i * .k f is the weight for the forward-following vector.δ i * represents the normalized direction vector of δ i .

Stability Analysis of Kinematics Models
The presented first-order and second-order models both adhere to a circular communication topology and the forward-following rule.In the context of the second-order model, it is evident from Equation ( 7) that, at each step, it weights the heading vector of the forward agent and its own to derive the subsequent updating direction, which can be expressed as follows: where φ i t+1 represents the direction angle of the i-th agent at time t+1, φ i t and φ i+1 t denotes the direction angle of the i-th and i+1-th agent at time t, respectively.
The first-order model can be seen as a specific case of the second-order model, where the difference between the positions of the i-th agent and i−1-th agent is defined as the direction for agent i.As a result, the first-order and second-order models are unified, consistent with the concept of nonlinear pursuit equations proposed in the related works [22][23][24].Therefore, the stability of the introduced first-order and second-order models can be justified using the stability analysis provided in their paper.Referring to Lemma 4 in Marshall's paper, it can be inferred that agents asymptotically converge to a regular n-sided polygon, implying that both the first-order and second-order kinematics models proposed here achieve uniform convergence on the circular trajectory [22].

Convergent Circle Analysis for First-Order Model
Simulation experiments and stability analysis indicate that for the first-order kinematics model, when the forward-following weight k f is less than or equal to 0.5, the model fails to converge.In fact, as can be seen from Figure 4, even if the system eventually converges to a circle, when k f is less than or equal to 0.5, the position of the i-th agent at the next moment will definitely fall outside the convergent circle instead of on the circle and will eventually be far away from the convergent circle.While k f is more significant than 0.5, the agents converge to a circle with a fixed center, forming a regular n-sided polygon evenly distributed along the convergent circle [22].When the model's step size v c and the forward-following weight k f are determined, the radius of the convergent circle is also determined.Based on the geometric relationship between agents during convergence, the following derivation is conducted.Based on the stability analysis, the agents converge to a circle with a fixed center, denoted as circle O, and the agents are positioned to form a regular n-sided polygon on this circle.The schematic diagram of the convergent circle is shown in Figure 4.In the diagram, the solid arc represents the convergent circle with center O. x i t , x i−1 t , and x i+1 t are located on the convergent circle, representing the positions of agent i, its neighbor i−1, and its neighbor i+1 at time t, respectively.Dashed lines represent the connections between them.x i t+1 represents the position of agent i at time t+1, and the connection between x i t and x i t+1 is represented by a solid line.The angle ∠ x i t Ox i t+1 corresponds to the angle that the agent rotates along the convergent circle in one time step, denoted as θ.The distance covered by the agent in each time step is constant and denoted as v c , so x i t x i t+1 = v c .Indeed, due to the movement of agents along circle O, x i t x i t+1 is a chord on the circle, and θ represents the corresponding angle along the circumference.Based on this, the radius of the convergent circle can be determined as Defining the tangent vector of point x i t on the convergent circle as γ i t , its direction aligns with δ i f + δ i b .Based on the geometric relationship between u i t and the vector γ i t , the angle between them is θ /2, which can be calculated using the following formula: where Then, expanding this formula, we have where α is the angle ∠ x i+1 t x i t x i−1 t+1 formed by δ i f and δ i b , and, by using the polygon interior angle sum formula, we obtain Substituting Formula (11) into Formula (10), it can further simplified by using basic properties of trigonometric functions as ).
Hence, the radius of the convergent circle and the convergent adjacent distance are calculated by the following formulas: By choosing a step size of v c =1, from Formula (13), the relationship between k f and the parameters of the convergent circle is derived, as shown in Figure 5.In this study, the distance between neighboring agents during convergence is defined as the convergence adjacent distance.The radius and the convergence adjacent distance can be used to characterize the properties of the convergent circle.For the convenience of indication, the vertical axis in Figure 5a represents the reciprocal of the radius of convergent circle r ∞ .Figure 5b represents the reciprocal of the convergent adjacent distance ∆d ∞ .It can be observed that, with a constant number of agents, as the forward-following weight k f increases, the radius of the convergent circle monotonically decreases, and the convergent adjacent distance also decreases.When k f =1, the convergent adjacent distance is independent of the number of agents and remains at one.Under the same k f , with an increasing number of agents, the convergent circle becomes larger, and the convergent adjacent distance also increases but approaches a limit.

Convergent Circle Analysis for Second-Order Model
For the second-order case, simulations demonstrate that for the second-order kinematics model given in Equations ( 5) and (7), the model can converge when the forwardfollowing weight k f ∈ (0, 1].Similar to the first-order model, the agents converge to a fixed-center circle, and the agents form a regular n-sided polygon evenly distributed on the converging circle.When the model's step length v c and the forward-following weight k f are determined, the final converging circle's radius is also determined.Similar to the derivation for the first-order model, based on the schematic diagram of the converging circles of neighboring agents at two consecutive time steps during the convergence of the second-order model, as shown in Figure 6.As shown in Figure 6, the intelligent agents converge to the circle O, forming a regular n-sided polygon.Points x i t−1 , x i t , x i t+1 , x i−1 t , and x i+1 t represent the positions of agent i and its neighbors at times t − 1, t, and t + 1, respectively, which are located on the convergent circle.Thus, v i t is collinear with x i t−1 x i t .∠ x i t Ox i t+1 represents the angle that an agent rotates on the convergent circle within a one-time step, denoted as θ.Since the distance the agents move in each time step is constant, equal to v c , it follows that x i t x i t+1 = v c .In addition, due to the agents' movement on the circle O, the distance x i t x i t+1 corresponds to a chord, and θ corresponds to the central angle subtended by this chord.Hence, the formula for the convergent circle's radius matches the one in Equation ( 8).Denoting the central angle corresponding to the chord x i t x i+1 t as α, we have Then, denote the angle between the velocity direction vector v i t of agent i at time t and the forward-following vector δ i f of agent i at time t as β 0 .Using the second-order control protocol from Equations ( 5) and ( 7), we can derive Further, we have where β 0 is the angle between v i t and δ i f , and β 2 is the angle between u i t and δ i f .Clearly, Using the property of the sum of interior angles in a triangle, we obtain Furthermore, at time t + 1, the velocity direction of agent i has shifted to vector u i t ; therefore, θ = β.In addition, considering that v i t and δ i f are unit vectors, Equation ( 17) can be simplified as ).
Equation ( 19) is an implicit equation that can be solved numerically to obtain the angle θ by using a solver.Then, the calculation of the radius of the convergent circle and the convergence adjacent distance follows the same approach as in Equation (13).
By selecting a step length v c = 1, the relationship between k f the convergent circle's parameters can be obtained, as shown in Figure 7.In Figure 7a, the relationship between k f and the reciprocal of the convergent circle radius is depicted, while Figure 7b shows the relationship between k f and the reciprocal of the adjacent distance during convergence.As can be observed, similar to the first-order model, the convergent circle size decreases with an increasing forward-following weight k f .Furthermore, under the same k f , more agents result in a larger convergent circle.

Convergent Speed Analysis of the First-Order Model
Convergence speed is a significant metric in multi-agent formation tasks, particularly in multi-robot formation.Setting control parameters reasonably and dynamically is essential to ensure the model converges quickly to the target circle.To quantify the relationship between convergence speed and model parameters, this section comprehensively analyzes convergence speed under various parameters.
Simulations were conducted by varying parameters and the number of agents.The number of agents n and parameter k f were selected from the Cartesian product of the following two sets: n ∈ [3, 4, 5, 6, 7, 8, 9, 10] and k f ∈ [0.55, 0.60, 0.65, 0.70, 0.80, 0.90, 1].Each set of parameters was repeated 20 times.The convergence speed was calculated using the following formula: where denotes the convergence threshold, in this section, = 0.001.e 0 represents the initial average absolute error of the group, and Step denotes the convergence step.
From Figure 8a,b, it can be observed that with an increase in the number of agents, the overall trend of convergence speed is decreasing, while with the increase in k f , the overall trend of convergence speed is increasing.In Figure 8c,d, taking the natural logarithm of the convergence speed shows that as the number of agents increases, the convergence speed decreases exponentially.

Convergent Speed Analysis for the Second-Order Model
In the subsequent part, we will perform a convergence speed analysis for the secondorder model.Similar to the analysis of convergence speed for the first-order model, the model's convergence speed is calculated based on Equation (20), with set to 0.001.The number of agents n and the parameter k f are chosen from the Cartesian product of the following two sets: n ∈ [3, 4, 5, 6, 7, 8, 9, 10] and k f ∈ [0.55, 0.60, 0.65, 0.70, 0.80, 0.90, 1].Each set of parameters is repeated in 20 simulations.In Figure 9a,b, it can be observed that with an increase in the number of agents, the overall trend in the convergence speed is decreasing, while, with an increase in k f , the overall trend in the convergence speed is increasing.In Figure 9c,d, taking the natural logarithm of the convergence speed shows a consistent pattern with the first-order model.As the number of agents increases, the convergence speed decreases exponentially.For the three-agents case, the convergence speed does not strictly increase with k f but instead shows a slight decrease when k f > 0.9.

Circular Formation Based on Dynamics Model of Robotic Fish
Although the proposed first-order model and the second-order model can both achieve the circular formation of agents, their implementations are based on the kinematics models of agents and cannot be directly applied to the formation of robots, especially for robotic fish with nonlinear dynamics.For the implementation of robotic fish formation, the dynamics characteristics need to be considered.Therefore, based on the dual-joint robotic fish dynamics model established in [25], and combined with the second-order model control protocol, a formation control scheme in a circular topology is designed.
It should be noted that the proposed method belongs to an indirect approach to circular formation, which does not directly specify the size of the target circle but determines it through model parameters.The size of the convergent circle is related to the forwardfollowing weight k f , the number of agents, and the step v c .When the number of agents is fixed, k f determines the angle θ (corresponding to angular velocity) that agents rotate on the circle during each update, and v c determines the distance (corresponding to linear velocity) that agents move during each update.For the kinematics model, given the number of agents, the desired convergent circle radius, and any three out of k f , v c , and θ, the remaining parameter can be calculated using Equations ( 13) and (19).However, for the dynamics model, its linear velocity and angular velocity are constrained and cannot directly correspond to the parameters of the kinematics model.On the other hand, the dynamics model updates in units of time, with a time step of 0.01 s.The kinematics model updates in terms of a distance step, and there is no physical correspondence between "distance" and "time".Therefore, in order to make the circular formation based on the robotic fish dynamics model match the parameters of the kinematics model when the number of agents and k f are the same, it is necessary to design a circular formation controller for robotic fish based on the proposed second-order model.

Design of Circular Formation Controller
In the defined circular formation, the angular velocity of agents during convergence in the robotic fish dynamics model is denoted as ω rob = η ts θ, and the linear velocity is denoted as v rob = η ts v c .Here, η ts is defined as the time scaling factor, which signifies how many steps of the kinematics model correspond to one second of the dynamics model.Therefore, for the circular formation of robotic fish based on the dynamics model, controlling the swimming speed v rob and angular velocity ω rob of the robotic fish allows us to match the parameters of the kinematics model's convergent circle.
Due to the negligible roll and pitch movements of robotic fish when swimming in a plane, this paper focuses solely on motion control within a two-dimensional plane to simplify the complexity of the formation problem.In this context, the control of the robotic fish's planar motion is broken down into yaw control, velocity control, and angular velocity control.Figure 10 illustrates the control system diagram for planar motion.Velocity control aims to align the swimming speed of the robotic fish with the desired speed.The fish's swimming speed is influenced by the amplitude, frequency, and bias of caudal fin oscillations, all of which are nonlinearly coupled.To simplify the controller complexity, assume a constant caudal fin oscillation amplitude and adjust the frequency of the caudal fin oscillations for the velocity control.A proportional integral (PI) control method is employed as follows: where f t denotes the oscillation frequency of the caudal fin's central pattern generator (CPG), and k p f t and k i f t are controller parameters.v e represents the velocity error.
The yaw control ensures that the robotic fish can track a target direction within the horizontal plane.When the desired yaw angle remains constant, the fish's direction of motion should stabilize at the desired yaw angle.Yaw control is achieved using a proportional-derivative (PD) controller: where b t is the oscillation bias of caudal fin, and k pbt and k dbt are controller parameters.ψ represents the yaw angle, and ψ e denotes the yaw angular error, which corresponds to the azimuthal difference between δ i and the robot's swimming direction v i * t in Equation (7).For circular formation, since the desired yaw angle is continually changing, adjusting the yaw angle velocity to stabilize at ω rob involves dynamically modifying k dbt through integral control action to eliminate angular velocity residuals: where k dbt −1 represents k dbt at the previous time step, k c is the controller parameter, and ∆t denotes the control cycle.
The velocity and yaw controller enables the robotic fish to form a convergent circle and maintain stable movement along the circular path.By automatically adjusting the parameter k dbt through yaw control, the angular velocity error can be eliminated, ensuring that the convergent circle's parameters match the calculation by Equations ( 13) and (19).
Figure 11 shows snapshots of the paths of five robotic fish in a circular formation with parameter k f = 0.3, a kinematics model step size of 0.5, and the scaling factor η ts = 0.6.

Simulation Results of Changing Speed and Scaling Factor
This part will discuss whether the convergent circle of the dynamics model matches the theoretical value.The proposed circular formation control algorithm will be validated from four aspects: the number of robotic fish, the forward-following parameter k f , the step v c , and the scaling factor η ts .
Firstly, verify whether the size of the convergent circle matches the formula given in Equation ( 13) when the step v c and scaling factor η ts vary.The steps are set to [0.3, 0.35, 0.4] m/s, and the scaling factors η ts are in the range of [0.8, 1, 1.2].k f is set to 0.5, and the number of robotic fish is three.The simulation results are depicted in Figure 12.The dashed lines in the figure represent the nearby convergent distance calculated according to Equation (13).The results from Figure 12 indicate that the scaling factor η ts does not affect the size of the convergent circle, while the step v c is directly proportional to the size of the convergent circle.

Simulation Results of Changing k f
Then, we will study the relationship between the size of the convergent circle and the theoretical value when k f varies.The swimming speed of the robotic fish is set to 0.5 m/s, and k f is varied within the range of [0.1, 0.3, 0.5, 0.7, 0.9, 1], corresponding to the angular velocities ω ∈ [5.36 • , 19.27  13.
The dashed lines in the figure represent the nearby convergence distance calculated based on Equation (13).These lines demonstrate that as k f varies, the size of the convergent circle in the robotic fish formation still conforms to Equation (13).

Simulation Results of Changing Number of Agents
Finally, we investigate the relationship between the size of the convergent circle and the theoretical values when the number of robotic fish varies.For this analysis, we set k f to 0.5, the step v c to 0.5 m/s, and vary the number of robotic fish as [3,5,7,9].The convergence behavior is depicted in Figure 14.In brief, by varying the number of robotic fish, the value of k f , the step v c , and the scaling factor η ts , the results demonstrate that the proposed control approach for robotic fish formation is equivalent to the second-order kinematics model in terms of forming the convergent circle.

Circular Formation Experiments of Biomimetic Robotic Fish
Using the two-joint biomimetic robotic fish developed in [25], a circular formation experiment was conducted with multiple robotic fish.The robotic fish takes a black koi fish as the bionic object, which is about 25 cm long, 5 cm wide, 9 cm high and weighs about 340 g.In this multi-robot circular formation experiment, the parameters were set as follows: k f = 0.2, implying that the angular velocity of the robotic fish during convergence was approximately 12 • ; the scaling factor was set to 1; the kinematics model step size was 0.12 for three robotic fish, 0.14 for four robotic fish, and 0.13 for five robotic fish, and the swimming velocity was 0.12 m/s, 0.14 m/s, and 0.13 m/s, respectively.The experimental scenario is depicted in Figure 15, and Figure 16 illustrates the paths of the three robotic fish during the experiment, with the horizontal and vertical axes measured in meters.According to Equation ( 13), for the three-robotic-fish formation, the convergent adjacent distance is 1.02 m, and the convergent angle is 60 • .The experimental results of the three-robotic-fish circular formation are shown in Figure 17.After 30 s, the average error in the adjacent distance is 0.12 m, and the average error in the convergent angle is 10.19 • .For the four-robotic-fish case, the convergent adjacent distance is 1.19 m, and the convergent angle is 90 • .The experimental results of the four-robotic-fish circular formation are shown in Figure 18.After 30 s, the average error in the adjacent distance is 0.23 m, and the average error in the convergent angle is 15.13 • .For the five-robotic-fish formation, the convergent adjacent distance is 1.11 m, and the convergent angle is 108 • .The experimental results of the multi-robot circular formation are shown in Figure 19.After 30 s, the average error in the adjacent distance is 0.38 m, and the average error in the convergent angle is 22.69 • .The results verify the effectiveness of the circular formation approach.

Discussion
This study is inspired by the milling motion observed in fish schools and similar behaviors in animals.For the task of circular formation in multi-robot fish systems, we proposed a circular formation algorithm based on a circular communication topology and the following rules.Specifically, we introduced first-order and second-order circular topology control protocols based on kinematics models, enabling multiple agents to converge uniformly to a target circle.Furthermore, we designed controllers for the robotic fish based on their kinematics models, enabling the simulation of circular formation with multiple robotic fish.The fish swarm demonstrated convergence to the target circle across various parameter settings, consistently maintaining the same target circle size as the kinematics model.Furthermore, we conducted experiments involving formation with three, four, and five robotic fish, confirming the effectiveness of our algorithm on the robotic fish platform.Through mathematical analysis, simulations, and experiments, we not only showcased the proposed algorithm's efficacy in addressing circular formation tasks for multi-robot systems but also suggested that the phenomenon of milling, commonly observed in the animal kingdom, is likely a result of creatures adhering to the following rule.This may lead to the spontaneous or externally induced formation of a stable circular communication topology, ultimately resulting in similar milling patterns across different species.
In Section 4, it can be observed that the errors in the convergence adjacent distance and angles are relatively large.This is partly due to the inherent under-actuation of the robotic fish and inconsistencies in the manufacturing of individual fish, which affect the control efficiency of the collective formation.Improving the consistency of the robotic fish and the control performance could lead to a smaller error in circular formation [26][27][28].
In the future, we will pursue research from two aspects.Firstly, we will delve further into the underlying mechanisms of animal collective behaviors, such as the bait-ball formation observed in fish schools.These widespread biological phenomena undoubtedly have principles that can be abstractly modeled.Based on the detailed study of the fish lateral line system [29], modeling, analyzing, and abstracting are carried out to validate conclusions drawn from the study of animal behaviors.Secondly, we will explore potential applications for circular formation tasks.Animals lack efficient communication and typically rely on passive perception rather than active communication for coordination.However, robots often construct wireless communication networks, offering higher communication effi-ciency.Animal milling behaviors are based on passive perception, which does not directly match the communication patterns of robots.Moreover, fish schooling behavior is a defense mechanism against predators, utilizing rapid rotations to lower the predators' efficiency.However, robots do not need to replicate such behaviors in their task scenarios.Seeking suitable application scenarios for the algorithm proposed in this paper remains a significant challenge, which may provide more insight into pursuit-evasion problems of multi-agent systems [30].

) a 1 a 5 a 4 a 3 a 2 Figure 2 .
Figure 2. Communication topology of the first-order model.

) a 1 a 5 a 4 a 3 a 2 Figure 3 .
Figure 3. Communication topology of the second-order model.

Figure 4 .
Figure 4. Schematic diagram of the convergent circle of the first-order model.

Figure 5 .
Figure 5. Relationship between convergent circle and parameters in the first-order model.(a) Relationship between convergent circle radius, k f , and the number of agents.(b) Relationship between convergent adjacent distance, k f , and the number of agents.

Figure 6 .
Figure 6.Schematic diagram of convergent circle of second-order model.

Figure 7 .
Figure 7. Relationship between convergent circle and parameters in the second-order model.(a) Relationship between the radius of convergent circle, k f , and the number of agents.(b) Relationship between convergence adjacent distance, k f , and the number of agents.

Figure 8 .Figure 9 .
Figure 8. Relationship between the number of agents and k f with convergence speed for firstorder model.(a) Convergence speed of the adjacent distance.(b) Convergence speed of the polygon angle.(c) Natural logarithm of the convergence speed of the adjacent distance.(d) Natural logarithm of the convergence speed of the polygon angle.

Figure 10 .
Figure 10.Control diagram of the robotic fish circular formation.

Figure 11 .
Figure 11.Circular formation simulation of five robotic fish.

Figure 12 .
Figure 12.Error diagrams of changing speed and scaling factor.

Figure 13 .
Figure 13.Error diagram of changing k f .
Adjacent distance errors.

Figure 14 .
Figure 14.Error diagram of changing the number of agents.

Figure 15 .
Figure 15.Experimental scenario of three robotic fish executing circular formation.

Figure 16 .
Figure 16.Path diagram of three robotic fish executing circular formation.

Figure 17 .
Figure 17.Results of the circular formation experiment of three robotic fish.

Figure 18 .
Figure18.Results of the circular formation experiment of four robotic fish.
Average angular error.

Figure 19 .
Figure 19.Results of the circular formation experiment of five robotic fish.