Spontaneous Synchronization of Beating Cilia : An Experimental Proof Using Vision-Based Control

This article investigates the formation of spontaneous coordination in a row of flexible 2D flaps (artificial cilia) in a chamber filled with a high viscous liquid (Re = 0.12). Each flap is driven individually to oscillate by a rotary motor with the root of the flap attached to its spindle axle. A computer-vision control loop tracks the flap tips online and toggles the axle rotation direction when the tips reach a pre-defined maximum excursion. This is a vision-controlled implementation of the so-called “geometric clutch” hypothesis. When running the control loop with the flaps in an inviscid reference situation (air), they remain in their individual phases for a long term. Then, the flaps are studied in the chamber filled with a highly viscous liquid, and the same control loop is started. The flexible flaps now undergo bending due to hydrodynamic coupling and come, after a maximum of 15 beats, into a synchronous metachronal coordination. The study proves in a macroscopic lab experiment that viscous coupling is sufficient to achieve spontaneous synchronization, even for a symmetric cilia shape and beat pattern.


Introduction
A wide range of biological systems use synchronization in their movement patterns [1,2], ranging from small-scale unicellular organisms to larger scale sperms and microswimmers [3,4].At low Reynolds numbers, this coordination is crucial for the propulsion of microswimmers or the generation of transport on ciliated walls.Motile cilia are found in many different tissues, from the brain [5] to the lung and the oviduct, and in many organisms, from Chlamydomonas [6] and Volvox [7,8] algae to Paramecium.It is the beat coordination of the invidiual cilia in the array that plays an essential role in the locomotion of sperm, the cleaning of breathing air, and the movement of oocytes in the fallopian tube [9][10][11][12].Metachronal synchronization results when cilia beat with a certain constant phase shift between each other.Similar coordination was observed not only for natural cilia but also for artificial cilia and flagella in microswimmer propulsion [12][13][14][15].Many studies have proposed ways to use the cilia propulsion mechanism in nanomechanics and microfluidic devices.Inspired by the progress made in the fabrication of small-scale flexible structures, the hydrodynamic interaction of artificial cilia has become again of interest to optimize such engineered systems.At low Reynolds numbers, the inertial terms in the Navier-Strokes equation can be neglected [3], simplifying the momentum equation to the linear Stokes equation with the pressure p, the velocity v and the dynamic viscosity µ of the fluid.As the right-hand side in Equation ( 1) is time-independent, it suggests that for symmetric beating patterns it may not be possible to spontaneously come into a synchronized beating pattern [16,17].However, the flexibility of the cilia adds another degree of freedom to the equation [18,19].Typically, the beating pattern of natural cilia is non-symmetric in time; in the power stroke, the cilium is moved forward while stretched out straight, and in the return stroke it is retraced back near the wall.This was adopted by computer-simulation studies for 2D arrays of cilia [20,21], together with the "geometric clutch" hypothesis that was first introduced by Lindemann [22,23].The results suggested that hydrodynamic interactions are sufficient to achieve spontaneous synchronization.Synchronization of model flagella was studied for macroscopic rotating helices [16] and colloidal spheres driven by optical traps [24,25].However, up to now simple 2D symmetric flaps have not been tested on their possible self-synchronization, to the best of our knowledge.Often in such 2D systems, synchronization was just imposed by the control of the actuators to achieve maximum efficiency of the transport [26].Therefore, the process of spontaneous synchronization has not been studied in such experiments under well-defined boundary conditions so far.The aim of the present work is to introduce such an experiment and to test the conditions that lead to spontaneous synchronization.Therefore, we use silicon rubber flaps as artificial, macroscopic cilia that are driven by motor-controlled rotating axles in which the flaps are clamped in at one end and are forced to undergo a symmetric beating motion at a predefined frequency.The integration of the fluidic system inside a feedback control loop is done by vision-controlled implementation of the "geometric clutch" hypothesis.Bending deflections of the flaps' tips are tracked while applying geometric thresholds for the toggling of the rotating direction of each individual motor separately.Hereby, we will experimentally prove that viscous coupling can force a row of individually beating flaps from random initial conditions into a synchronized beating pattern, forming a metachronal wave.The paper is structured as follows: Section 2 describes the methods, and results are shown in Section 3. Finally, discussion and conclusions are given in Section 4.

Fluidic Chamber with Artificial Cilia
The experimental setup is designed to mimic a row of five artificial cilia (n = 5) with individual beating actuations at their roots.Each of the flexible cilia with a length L = 36 mm and width W = 20 mm is made of a silicon rubber sheet of thickness T = 2 mm (material polydimethylsiloxane (PDMS), Young's modulus E = 5.96 MPa).The Young's modulus was qualified with an impulse-response test as described in Favier et al. [27].The obtained value is well within the range documented for standard silicon rubber (Young's modulus 1-50 MPa) [28].
Each of the flaps (indicated with index i from 1 to 5) is clamped with one short end to the spindle axle of a rotary DC servo motor.The individual DC motors with the flaps are fixed on a support frame in a horizontal row such that the flaps' free ends point vertically down with an interspacing between each of the spindles of ∆s = 20 mm.The row is inserted from top in the center of a closed chamber with a squared (60 mm × 60 mm) cross section, see Figure 1.The chamber is made out of transparent Perspex to allow optical access to the system.Both lateral sides of the chamber can be closed and liquid can be filled such that the flaps are fully submerged up to the top wall of the chamber.As working liquid, pure glycerol is used with a density of ρ = 1260 kg•m −3 and a viscosity of µ = 0.95 Pa•s at room temperature.
For the present studies, the flap beating cycle is initiated with a constant beating frequency in which A is the preset beat amplitude A = 8.3 mm, which is set to ensure that each two neighboring flaps do not just touch each other at a maximum angle of α max = ±13 the rotary velocity of the DC motor .α = ∂α/∂t = 0.3 rad/s.Small adjustments of the individual amplitudes A i are necessary in the range of ∆A i < 0.1 mm to take into account small uncertainties in the DC motors controller response between programmed angular speed and the actual measured ones.Therefore, a calibration procedure is done for the individual amplitudes A i to end up with the same frequency f beat for each individual flap i.The characteristic parameter of the flow is the Reynolds number, which is defined with the preset amplitude A and the tip speed v reading For reference measurements with the control loop in an inviscid environment, the chamber is emptied to ensure that no viscous coupling affects the flaps' motion.
For reference measurements with the control loop in an inviscid environment, the chamber is emptied to ensure that no viscous coupling affects the flaps' motion.

Computer Vision Control
The control loop presented herein requires the online tracking of the tips of the cilia and the feedback of the position to the main control processor and the motors (see Section 2.3).Therefore, a high-speed camera (ProcImage 500-Eagle high-speed camera, 1280 × 1024 px 2 , Photon Lines Ltd, Bloxham, UK) is monitoring the tips of the five flaps recognized by a fluorescent marker glued onto each tip.Two mirrors are used to facilitate the illumination and imaging of the system, as shown in Figure 2. A LED light source (IL-106G Green LED Illuminator, HARDSoft, Krakow, Poland) is adjusted to continuous illumination mode.The fluorescent material glued onto each tip scatters the emitted light back to the camera.An optical filter (optical edge filter, cut-off wavelength 550 nm, Novasoft, Aarhus, Denmark) is used to reduce stray light from the back of the frame, which is painted in black.The markers are arranged in a zigzag pattern to help identify the tips separately as single objects in the image processing phase, even when the flaps are getting close to each other.The scene is recorded with a rate of 505 fps at 700 × 225 px 2 image format, and a typical image of the 5 visible marker dots is shown in Figure 3.Each image is processed online to extract the tips coordinates using image binarization and centroid detection.The resulting coordinates of the tip markers are then fed to the main control processor at a sampling frequency of fS = 15 Hz.The current implementation allows a maximum observation period of approximately 1 min, corresponding to a total number of 18-20 cycles.

Computer Vision Control
The control loop presented herein requires the online tracking of the tips of the cilia and the feedback of the position to the main control processor and the motors (see Section 2.3).Therefore, a high-speed camera (ProcImage 500-Eagle high-speed camera, 1280 × 1024 px 2 , Photon Lines Ltd, Bloxham, UK) is monitoring the tips of the five flaps recognized by a fluorescent marker glued onto each tip.Two mirrors are used to facilitate the illumination and imaging of the system, as shown in Figure 2. A LED light source (IL-106G Green LED Illuminator, HARDSoft, Krakow, Poland) is adjusted to continuous illumination mode.The fluorescent material glued onto each tip scatters the emitted light back to the camera.An optical filter (optical edge filter, cut-off wavelength 550 nm, Novasoft, Aarhus, Denmark) is used to reduce stray light from the back of the frame, which is painted in black.The markers are arranged in a zigzag pattern to help identify the tips separately as single objects in the image processing phase, even when the flaps are getting close to each other.The scene is recorded with a rate of 505 fps at 700 × 225 px 2 image format, and a typical image of the 5 visible marker dots is shown in Figure 3.Each image is processed online to extract the tips coordinates using image binarization and centroid detection.The resulting coordinates of the tip markers are then fed to the main control processor at a sampling frequency of f S = 15 Hz.The current implementation allows a maximum observation period of approximately 1 min, corresponding to a total number of 18-20 cycles.
is recorded with a rate of 505 fps at 700 × 225 px 2 image format, and a typical image of the 5 visible marker dots is shown in Figure 3.Each image is processed online to extract the tips coordinates using image binarization and centroid detection.The resulting coordinates of the tip markers are then fed to the main control processor at a sampling frequency of fS = 15 Hz.The current implementation allows a maximum observation period of approximately 1 min, corresponding to a total number of 18-20 cycles.

Feedback Loop System
A mechatronic system is designed to actuate, monitor, and control the row of flaps.The system shown in Figure 4 is composed of three subsystems: (i) the main plant, which consists of the chamber, the flaps, and the five motors; (ii) the supervisor control system for monitoring the flaps and controlling the servo motors with the desired speed and direction; and (iii) the optical recording system, which consists of the high-speed camera and the illumination.The energy flows into and out of the control system as follows: The five flaps are controlled to start beating from random initial positions (t = 0).The motors' speeds and beating amplitudes are set by the controllers to the same calibrated values.Recall that the steady-state speed of the DC motor is proportional to the motor torque and therefore to the torque imposed at the root of the flaps.Then, the motion is started and continuously monitored by tracking the tip markers.Each motor is advised to reverse direction once the attached flap's tip reaches its predefined maximum position.

Feedback Loop System
A mechatronic system is designed to actuate, monitor, and control the row of flaps.The system shown in Figure 4 is composed of three subsystems: (i) the main plant, which consists of the chamber, the flaps, and the five motors; (ii) the supervisor control system for monitoring the flaps and controlling the servo motors with the desired speed and direction; and (iii) the optical recording system, which consists of the high-speed camera and the illumination.The energy flows into and out of the control system as follows: The five flaps are controlled to start beating from random initial positions α i (t = 0).The motors' speeds and beating amplitudes are set by the controllers to the same calibrated values.Recall that the steady-state speed of the DC motor is proportional to the motor torque and therefore to the torque imposed at the root of the flaps.Then, the motion is started and continuously monitored by tracking the tip markers.Each motor is advised to reverse direction once the attached flap's tip reaches its predefined maximum position.
This implementation in the control loop is equivalent to the geometric clutch hypothesis introduced by Lindemann [22,23], which suggests that the cilium motor reverses the beating direction once a certain terminal position of the cilia tip is reached.Herein, the bending deformation is induced by the viscous coupling of the surrounding liquid and neighboring flaps and causes the tip of the flap to deviate from the simple zig-zag trajectory prescribed by the rotating axle.This allows for phase variations in the system due to force-induced bending deformation of the flexible cilia. of the control system as follows: The five flaps are controlled to start beating from random initial positions (t = 0).The motors' speeds and beating amplitudes are set by the controllers to the same calibrated values.Recall that the steady-state speed of the DC motor is proportional to the motor torque and therefore to the torque imposed at the root of the flaps.Then, the motion is started and continuously monitored by tracking the tip markers.Each motor is advised to reverse direction once the attached flap's tip reaches its predefined maximum position.This implementation in the control loop is equivalent to the geometric clutch hypothesis introduced by Lindemann [22,23], which suggests that the cilium motor reverses the beating direction once a certain terminal position of the cilia tip is reached.Herein, the bending deformation is induced by the viscous coupling of the surrounding liquid and neighboring flaps and causes the tip of the flap to deviate from the simple zig-zag trajectory prescribed by the rotating axle.This allows for phase variations in the system due to force-induced bending deformation of the flexible cilia.

Reference Tests in Inviscid Environment
This step is required to guarantee that the measured time-variant phase differences among neighboring beating flaps in the presence of the liquid can be regarded as a reaction to the hydrodynamic forces of the surrounding fluid and the viscous coupling among the flaps and are notimposed by timing issues raised by uncertainties in the image processing and motor control.Therefore, the control loop is first run in an inviscid environment with the liquid chamber fully emptied.Theoretically, in the reference case the phase relationship between the flaps in the row should remain exactly the same for infinite number of beats after the start of the control loop.All uncertainties such as those in the imaging processing to detect the tips centroids and in the calibrated motor speeds can accumulate to an observable variability of the phase, which can then be used to quantify an upper limit within the observation period.The performance in air is then compared to the control loop run with the flaps in the liquid.The coordinates of the tips are by-passed to a recorder during the motion tracking while the control loop is running.

Results
Figure 5 shows the recorded position trajectories for both situations: the control loop in reference situation and the control loop including strong viscous coupling.
Indeed, the results for the reference situation demonstrate that the phase relationship between the flaps in the row remains approximately constant over the number of recorded cycles.A small variation of flap i = 3 is seen after 15 cycles, which is less than 10% phase drift.All flaps perform a zig-zag-type oscillatory motion while they remain in straight shape.In contrast, for the system with high viscous coupling, the phase differences are time-variant.The flexible flaps start bending as a reaction to the stronger fluid forces, which leads to a deviation of tip location relative to the position prescribed by the angle of the motor with the flap in straight shape.Maximum values of about 2 mm were observed, which corresponded to roughly 25% of the beat amplitude.The system comes after about 15 cycles into a coordinated synchronous motion when the phase lag among the tips remains approximately constant in space (between pairs) and in time.The straight black lines in Figure 5 bottom show this observation clearly.The lines are passing through the instants where the individual beats reach their peak position.Only flap i = 1 appears to be lagging behind the other tips.There are some possible explanations for this observation, which are discussed in the next section.
the control loop run with the flaps in the liquid.The coordinates of the tips are by-passed to a recorder during the motion tracking while the control loop is running.

Results
Figure 5 shows the recorded position trajectories for both situations: the control loop in reference situation and the control loop including strong viscous coupling.).Note that we selected for visibility reasons in the plot for air the initial phase shift between pairs of cilia to be at the extreme value of 180°, which results in out-of-phase beating, which is completely outside the range for a metachronal synchronization.Therefore, the black lines in the upper plot do not indicate any synchronization, as there is no interaction between the flaps.).Note that we selected for visibility reasons in the plot for air the initial phase shift between pairs of cilia to be at the extreme value of 180 • , which results in out-of-phase beating, which is completely outside the range for a metachronal synchronization.Therefore, the black lines in the upper plot do not indicate any synchronization, as there is no interaction between the flaps.
For further characterization of the metachronal coordination between neighboring beating flaps n i and n i−1 , we define the cycle ratio CR from the phase differences of reversal times T peak in the corresponding beat cycles: Figure 6 shows the temporal evolution of the phase differences between each couple of neighboring flap tips for flap i = 2-5.The non-dimensional phase differences CR n3-5 of the four synchronized flaps are found to converge to 15-25% of the cycle period between each neighboring couple.This is close to the phase lag observed in natural cilia to obtain optimum transport conditions [3,21].For comparison, the black line shows the behavior of one pair of flaps CR n4 for the reference case in air.
Figure 6 shows the temporal evolution of the phase differences between each couple of neighboring flap tips for flap i = 2-5.The non-dimensional phase differences CRn3-5 of the four synchronized flaps are found to converge to 15-25% of the cycle period between each neighboring couple.This is close to the phase lag observed in natural cilia to obtain optimum transport conditions [3,21].For comparison, the black line shows the behavior of one pair of flaps CRn4 for the reference case in air.

Discussion and Conclusions
Previous numerical models of the cilia beating patterns have simulated the spontaneous emergence of metachronal waves due to hydrodynamic interactions [18,21].Included in these theoretical models is often the geometrical clutch hypothesis, which was introduced by Lindemann [22,23] to explain how the cilia motor is reversing the direction.It is assumed that once a cilium The gray shaded region is the range of cycle ratios reported in literature for metachronal waves of ciliated walls to achieve optimum fluid transport [21].

Discussion and Conclusions
Previous numerical models of the cilia beating patterns have simulated the spontaneous emergence of metachronal waves due to hydrodynamic interactions [18,21].Included in these theoretical models is often the geometrical clutch hypothesis, which was introduced by Lindemann [22,23] to explain how the cilia motor is reversing the direction.It is assumed that once a cilium experiences certain bending, it retracts back following the beat stroke.Up to now, there has been no experimental proof of this spontaneous synchronization for rows of artificial 2D cilia under controlled conditions.In our study, we have verified the previous hypotheses in macroscopic experiments by documenting the spontaneous synchronization in a row of flaps due to their hydrodynamic interaction.This has been made possible by (i) introducing system flexibility through the use of silicon rubber flaps as artificial cilia, which start bending in the highly viscous environment at Re = 0.12 and (ii) introducing online computer vision as the feedback methodology within the control loop that mimics the geometric clutch hypothesis.Since the silicon flaps bend as a reaction to the hydrodynamic coupling, it was important to optically track the flap tips for control of the beat reversal commands.This additional degree of freedom in the cilia tip motion relative to the imposed motor force allows coordination to take place.This load response has been measured and discussed for natural cilia and flagella in [29,30].Our system is therefore equivalent to the key property of natural cilia required for synchronization, as they are able to change speed (or tip speed) in response to hydrodynamic load.Hydrodynamic coupling is herein the source driving the system to synchronization, as nothing else has changed between the two control experiments than the fluid.This agrees with the conclusions given in [31].
The results show that a system of 5 individually beating flaps can form spontaneous synchronization after several cycles, while the same system in an inviscid environment remains fixed in its initial phase relationships.The observed convergence to phase differences of 15-25% of cycle duration agrees with the range found in natural cilia coordination for optimum flow transport [3,21].This hints at a possible relaxation of the control system to a minimum of dissipation in the liquid [21].Additional flow measurements in our experiment are planned to gain details of the flow field during synchronization.As seen from the results, the synchronization is not perfect, as the differences in phase do not converge to the exact same values.Flap number i = 1 could not reach the phase relationship to the next as close as the others during the recording period of total number of 18 cycles.This could be related to the non-symmetric boundary effects for flap number i = 1 and i = 5.For both, the beating in direction away from the inner neighboring flaps is less influenced by viscous coupling than it is for the inner ones.As a result, the phase synchronization therefore may drift towards lock-on to either of the end flaps.A definite answer to the boundary effects can only be found by testing a chamber with a circular row of flaps, which eliminates these effects.However, variations in the phase shifts were also observed in the numerical simulations, even for imposed periodic boundary conditions [21].Real-time control is, in general, sensitive to time constrains and system internal delays, which may lead to aliasing effects [32].It takes 0.067 s from the detection of the terminal position to the action of the flap reversing the beat, given by the sampling frequency.At maximum tip speed, this delay equates to a possible variability of the amplitude of ∆A/A i = 0.06.This uncertainty introduces a time variant parameter in the system synchronization.However, this random effect can be fairly neglected compared to the phase shift introduced by the bending of the flaps due to hydrodynamic interactions, as proven by comparison to the reference measurements in air, see also Figure 6.In addition, during the calibration stage, the five motors were found to have slightly different speeds.Nevertheless, synchronization can even emerge in the presence of different intrinsic beat frequencies and oscillator noise as long as hydrodynamic coupling is strong enough.Fluctuations have been observed for natural systems, too, see ref. [6].
A number of additional results is shown here to address several questions of importance regarding statistics and reproducibility of the results.Figure 7 shows another experiment in the viscous liquid at random initial conditions.The results prove that the system again gets into metachronal synchronization after about 15 cycles.Interestingly, while in Figure 5 the metachronal wave runs from left to right (or from flap i = 5 to i = 1), it is the opposite for the other experiment shown in Figure 7.As the beating profile of the cilia is in principal time-reversible and the cilia shape is symmetric, and all cilia are uniformly spaced, we assume that small deviations from symmetry in the motor response or in the shape of the flaps or the initial conditions play a role in selection of direction.Sometimes, we even could observe spontaneous symmetry breaking at times larger than 1 min.However, we could not judge these observations as physically relevant, as they only appeared at much later times, and possible accumulations and aliasing affects as discussed above may then play a role.The number of cycles to get into synchronization depends on the initial conditions, and we observed in several experiments a variation of this time-scale.Note, that we quantify this time of synchronization as the instant when all cilia pairs get into the phase relationship within the range marked as a gray-shaded box in Figure 6.When the cilia start with a constant phase shift representing already a metachronal state, this time-scale is zero.On the other hand, all other experiments showed a maximum of equal or less than 15 cycles until synchronization was reached in the present system at the given properties of the The results prove that the system again gets into metachronal synchronization after about 15 cycles.Interestingly, while in Figure 5 the metachronal wave runs from left to right (or from flap i = 5 to i = 1), it is the opposite for the other experiment shown in Figure 7.As the beating profile of the cilia is in principal time-reversible and the cilia shape is symmetric, and all cilia are uniformly spaced, we assume that small deviations from symmetry in the motor response or in the shape of the flaps or the initial conditions play a role in selection of direction.Sometimes, we even could observe spontaneous symmetry breaking at times larger than 1 min.However, we could not judge these observations as physically relevant, as they only appeared at much later times, and possible accumulations and aliasing affects as discussed above may then play a role.The number of cycles to get into synchronization Fluids 2018, 3, 30 9 of 12 depends on the initial conditions, and we observed in several experiments a variation of this time-scale.Note, that we quantify this time of synchronization as the instant when all cilia pairs get into the phase relationship within the range marked as a gray-shaded box in Figure 6.When the cilia start with a constant phase shift representing already a metachronal state, this time-scale is zero.On the other hand, all other experiments showed a maximum of equal or less than 15 cycles until synchronization was reached in the present system at the given properties of the liquid.We expect also a dependency on the Reynolds-number, which is left open for future studies with liquids of varying viscosity.
For further illustration of the variability in the motion pattern, we calculated the variance in the zig-zag type motion by superposition of the individual "teeth" in the zig-zag pattern.For reference, we fixed the tip of each tooth at the position of cycle time equal to 0. From the data superposed in this way, we calculated the probability function to be within a certain radial distance to the mean of the triangular shape using a box-counting method.The resulting distribution is shown as color plot in Figure 8 for flap i = 2. Dark red color means a probability of 1 that this position is reached in all successive cycles in a similar time relative to the tip of the teeth (defined by the search radius).As demonstrated, the air experiment shows a rather sharp-edged triangle with small variations between individual beats.In comparison, the synchronization process in the viscous environment causes larger deviations, especially at the flap reversals, where neighboring flaps get closer to each other.Further evidence of this variability is demonstrated by comparing the normalized auto-correlation profiles of the flap motion as given in Figure 9.The normalized autocorrelation of a perfect zig-zag motion pattern would result in extrema of 1 and −1 at time-lags of even and odd multiples of half of the cycle ratio.For the air experiments, the positive peaks are all above a correlation value of 0.95, shown by the "+" type symbols in the plot.In comparison, in the viscous case successive cycles underlie the documented small variability during the synchronization process; therefore, the correlation peaks decrease with larger time-lags until correlation is fully lost after 8-9 cycles.Similar profiles can be shown also for other flaps in the row.Peaks of normalized autocorrelation function of the zig-zag motion profile for a selected flap in air (red symbol "+") and glycerol (blue symbol "x") for a total period of 1 min.For illustration, the continuous auto-correlation profile is added for the case of glycerol as a solid line in blue.
In the future, further improvements will be possible on the control loop timing to increase the sampling frequency and thus allow longer periods of study.This can be achieved by embedding the PC monitoring and control logic into the high-speed camera processor itself.It will allow or the Figure 9. Peaks of normalized autocorrelation function of the zig-zag motion profile for a selected flap in air (red symbol "+") and glycerol (blue symbol "x") for a total period of 1 min.For illustration, the continuous auto-correlation profile is added for the case of glycerol as a solid line in blue.
In the future, further improvements will be possible on the control loop timing to increase the sampling frequency and thus allow longer periods of study.This can be achieved by embedding the PC monitoring and control logic into the high-speed camera processor itself.It will allow or the extension of the studies to a wider range of Reynolds-numbers and parametric variations of the geometry of the cilia.One major difference between biological cilia and the flaps in our study is their 2D shape with rectangular cross-section.Nevertheless, we expect similar synchronization for cylinders as long as we are in the Stokes-flow regime where the shape of the body loses importance in viscous forces.Furthermore, the chamber offers optical access to capture the flow details using Particle Image Velocimetry simultaneously while the control loop is running.It might help to further investigate the above-claimed hypothesis of the relaxation of the control system towards a minimum of dissipation in the liquid.

Fluids 2018, 3 ,
x FOR PEER REVIEW 3 of 11 ones.Therefore, a calibration procedure is done for the individual amplitudes Ai to end up with the same frequency fbeat for each individual flap i.The characteristic parameter of the flow is the Reynolds number, which is defined with the preset amplitude A and the tip speed v reading = ≈ 0.12

Figure 1 .
Figure 1.Schematic view of the closed working chamber: (a) front view of the chamber with the five flexible flaps at resting conditions and (b) side view of one flap clamped on the short side to the axle of one of the five motors.

Figure 1 .
Figure 1.Schematic view of the closed working chamber: (a) front view of the chamber with the five flexible flaps at resting conditions and (b) side view of one flap clamped on the short side to the axle of one of the five motors.

Figure 2 .
Figure 2. Schematic view of the optical set-up to record the tip motion of the artificial cilia.Figure 2. Schematic view of the optical set-up to record the tip motion of the artificial cilia.

Figure 2 .
Figure 2. Schematic view of the optical set-up to record the tip motion of the artificial cilia.Figure 2. Schematic view of the optical set-up to record the tip motion of the artificial cilia.

Fluids 2018, 3 , 11 Figure 3 .
Figure 3. Original image of the flap tips after binarization and color inversion for a single instant in the process of self-organization.For reference, we overlaid circles marking the centers of the centroids and the flap tip contours (dashed lines).The oscillatory motion is along the horizontal image axis, perpendicular to the flap span.The left-most flap is i = 1, and the right most is i = 5 in rising order.Note that the roots of the model cilia are uniformly spaced, as shown in Figure 1.The image shown here is what is seen by the camera as the tips of the cilia while they are in motion.Thus, the tip contours are not uniformly spaced in the image.

Figure 3 .
Figure 3. Original image of the flap tips after binarization and color inversion for a single instant in the process of self-organization.For reference, we overlaid circles marking the centers of the centroids and the flap tip contours (dashed lines).The oscillatory motion is along the horizontal image axis, perpendicular to the flap span.The left-most flap is i = 1, and the right most is i = 5 in rising order.Note that the roots of the model cilia are uniformly spaced, as shown in Figure 1.The image shown here is what is seen by the camera as the tips of the cilia while they are in motion.Thus, the tip contours are not uniformly spaced in the image.

Figure 4 .
Figure 4. Schematic view of the integrated mechatronics (a) and the feedback loop chart (b).

Figure 4 .
Figure 4. Schematic view of the integrated mechatronics (a) and the feedback loop chart (b).

Figure 5 .
Figure 5. Position trajectories of the five flap tips during the feedback control loop; (a): flaps in reference case (air); (b): flaps in highly viscous liquid.The black lines indicate the phase relationships between maxima of neighboring flaps over time.The position corresponds to the centroid coordinate of the flap tip marker in pixel units along the horizontal image axis in Figure 3.The colors indicate the different flaps (red: flap i = 1, blue: flap i = 5).Note that we selected for visibility reasons in the plot for air the initial phase shift between pairs of cilia to be at the extreme value of 180°, which results in out-of-phase beating, which is completely outside the range for a metachronal synchronization.Therefore, the black lines in the upper plot do not indicate any synchronization, as there is no interaction between the flaps.

Figure 5 .
Figure 5. Position trajectories of the five flap tips during the feedback control loop; (a): flaps in reference case (air); (b): flaps in highly viscous liquid.The black lines indicate the phase relationships between maxima of neighboring flaps over time.The position corresponds to the centroid coordinate of the flap tip marker in pixel units along the horizontal image axis in Figure 3.The colors indicate the different flaps (red: flap i = 1, blue: flap i = 5).Note that we selected for visibility reasons in the plot for air the initial phase shift between pairs of cilia to be at the extreme value of 180 • , which results in out-of-phase beating, which is completely outside the range for a metachronal synchronization.Therefore, the black lines in the upper plot do not indicate any synchronization, as there is no interaction between the flaps.

Figure 6 .
Figure 6.Evolution of the flap-to-flap phase relationship in the first 18 cycles, shown in terms of the cycle ratio CRn3-5 between pairs of flaps 2-5.For comparison, the black line shows the behavior of one pair of flaps CRn4 for the reference case in air.The curve fitting is done with a smoothing spline.The gray shaded region is the range of cycle ratios reported in literature for metachronal waves of ciliated walls to achieve optimum fluid transport[21].

Figure 6 .
Figure 6.Evolution of the flap-to-flap phase relationship in the first 18 cycles, shown in terms of the cycle ratio CR n3-5 between pairs of flaps 2-5.For comparison, the black line shows the behavior of one pair of flaps CR n4 for the reference case in air.The curve fitting is done with a smoothing spline.The gray shaded region is the range of cycle ratios reported in literature for metachronal waves of ciliated walls to achieve optimum fluid transport[21].

Figure 7 .
Figure 7. Repetition of the control experiment in glycerol; compare this to Figure 5b.Note, again, emergence of metachronal coordination after about 15 cycles.Herein, the observed metachronal wave runs in opposite direction to the one documented in Figure 5, bottom.For explanation of the lines see the figure caption in Figure 5.

Figure 7 .
Figure 7. Repetition of the control experiment in glycerol; compare this to Figure 5b.Note, again, emergence of metachronal coordination after about 15 cycles.Herein, the observed metachronal wave runs in opposite direction to the one documented in Figure 5, bottom.For explanation of the lines see the figure caption in Figure 5.

Fluids 2018, 3 ,Figure 8 .
Figure 8. PDF of the individual "teeth" to overlap with the mean in the zig-zag motion profile for a selected flap in air (a) and glycerol (b) for a total period of 1 min.Blue values show a larger deviation from the mean, while red is aligned with the mean of all teeth.Note that all teeth are fixed with their maximum position at beginning of cycle time 0.

Figure 8 .
Figure 8. PDF of the individual "teeth" to overlap with the mean in the zig-zag motion profile for a selected flap in air (a) and glycerol (b) for a total period of 1 min.Blue values show a larger deviation from the mean, while red is aligned with the mean of all teeth.Note that all teeth are fixed with their maximum position at beginning of cycle time 0.

Figure 8 .
Figure 8. PDF of the individual "teeth" to overlap with the mean in the zig-zag motion profile for a selected flap in air (a) and glycerol (b) for a total period of 1 min.Blue values show a larger deviation from the mean, while red is aligned with the mean of all teeth.Note that all teeth are fixed with their maximum position at beginning of cycle time 0.

Figure 9 .
Figure9.Peaks of normalized autocorrelation function of the zig-zag motion profile for a selected flap in air (red symbol "+") and glycerol (blue symbol "x") for a total period of 1 min.For illustration, the continuous auto-correlation profile is added for the case of glycerol as a solid line in blue.