1 Introduction↩︎

Multimodal microrobots capable of working collaboratively in large autonomous teams have the potential to enhance our technological capabilities in complex environments, such as agricultural fields, confined spaces inside earthquake-collapsed buildings, narrow cavities within coral reefs, and deep waters. During the past few years, the robotics community has made significant technological progress in aquatic[1]–[10], terrestrial[10]–[18], and aerial microrobotics[19]–[24]. However, fundamental challenges in energy storage, manufacturing, control, actuation, sensing, and embedded intelligence must be addressed before the vision of autonomous collaborative microrobotics can become a reality. One major remaining challenge for mm-to-cm–scale locomoting robots is the ability to generate large forces relative to their size and weight to achieve high maneuverability. Regarding cm-scale aquatic propulsion, we have recently developed two new devices that utilize soft tails to produce thrust through fluid–structure interaction (FSI), driven by actuators powered by shape-memory alloy (SMA) wires. Namely, the single-tail propulsor created to drive the \(45\)-mg VLEIBot swimmer presented in[1] and shown in Fig. 1 (a), and the dual-tail propulsor developed to achieve the two-dimensional (\(2\)D) controllability on the surface of water of the \(90\)-mg VLEIBot⁺ swimmer also presented in[1] and shown in Fig. 1 (b). Together, these innovations enabled the creation of the \(900\)-mg autonomous VLEIBot⁺⁺ swimmer shown in Fig. 1 (c)[2]. Although the functionality of these insect-scale thrusters has been empirically demonstrated, the characteristics of the forces produced by the two different designs have not yet been investigated.

Here, we present a combined theoretical and experimental framework to characterize the forces produced by the two propulsors described above. Specifically, we consider a model for the soft tail driving a propulsor using reactive theory[25], which provides valuable insights for design, and for parameter analysis and selection. Also, we compare predictions obtained through this model with measurements collected employing a mm-scale μN-resolution force sensor designed and built using the techniques presented in[26] and the references therein. With this sensing device, for the tested single-tail propulsor operating at a frequency of \(1\) Hz, we measured maximum and cycle-averaged force values with multi-test means of \(0.45\) mN and \(2.97\) μN, respectively. For the dual-tail propulsor operating at a frequency of \(1\) Hz, we measured maximum and cycle-averaged force values with multi-test means of \(0.61\) mN and \(22.6\) μN, respectively. Comparisons of the experimental data with the results obtained through the reactive-force model indicate that this mathematical description is reasonably accurate in predicting the thrust generated by a single-tail propulsor, while failing at predicting the force produced by the dual-tail device. These observations suggest that the derived reactive model of force generation is applicable to situations in which a slender single tail interacts with a surrounding fluid in isolation. However, it fails at capturing the generation of hydrodynamic forces produced by the interaction of wake structures associated with two or more propulsive tails. We hypothesize that this issue can be remedied by incorporating a convective term into the model, which in this study we neglected a priori because during the performance of force tests, the tail does not move forward relative to the inertial frame of reference.

The analyses and experimental results presented in this paper represent the first systematic study of the mechanism of hydrodynamic force generation through FSI phenomena at the insect scale. The rest of the paper is organized as follows. Section [SECTION02] presents the derivation of the proposed reactive model of force generation. Section [SECTION03] describes the design and fabrication of the μN-resolution sensor used during force characterization. Section [SECTION04] presents the collected experimental data. Last, Section [SECTION05] summarizes the presented research.

2 Reactive-Force Model↩︎

In this section, we present the theoretical framework used to model the forces produced by VLEIBot-type propulsors. Employing the large-amplitude elongated-body approach presented in[25], we model the lateral undulation of a propulsive thin tail as a traveling wave with the characteristics shown in Fig. 2 (a). As seen in Figs. 2 (a) and (b), the triplet of unit vectors \(\boldsymbol{\mathcal{N}}=\{\boldsymbol{n}_1,\boldsymbol{n}_2,\boldsymbol{n}_3\}\) denotes the inertial frame of reference, fixed to the planet Earth, used to describe the kinematics of the tail and analyzing the dynamics of the fluid in its close proximity. Similarly, the triplet of unit vectors \(\boldsymbol{\mathcal{B}}=\{ \boldsymbol{b}_1,\boldsymbol{b}_2,\boldsymbol{b}_3\}\) denotes the body-fixed frame of reference. In agreement with these definitions, the directions of \(\boldsymbol{b}_1\), \(\boldsymbol{b}_2\), and \(\boldsymbol{b}_3\) coincide with those of \(\boldsymbol{n}_1\), \(\boldsymbol{n}_2\), and \(\boldsymbol{n}_3\) at any given instant, and \(\boldsymbol{b}_1\) is collinear with the midline of the tail when in the extended condition—i.e., the planform of the tail is perpendicular to the \(\boldsymbol{n}_1\)–\(\boldsymbol{n}_2\) plane. Consistent with these idealizations, we assume that the tail does not twist during undulation and, therefore, that the forces acting along the \(\boldsymbol{b}_3\) direction are negligible. To track the forces acting on the tail—all assumed inside the \(\boldsymbol{n}_1\)–\(\boldsymbol{n}_2\) plane—we define the Lagrangian coordinate \(s \in \left[0,l\right]\), according to which the arc length \(s\) is the distance from the origin of \(\boldsymbol{\mathcal{B}}\)—along the orthogonal projection of the tail onto \(\boldsymbol{n}_1\)–\(\boldsymbol{n}_2\)—to a vertical infinitesimal element of volume, and \(l\) is the total length of the tail.

The position of a point on the centerline of the tail, relative to the origin of \(\boldsymbol{\mathcal{N}}\), at the Lagrangian coordinate \(s\) and time \(t\) is given by \(\boldsymbol{r}(s,t) = \left[x(s,t)~y(s,t)~0\right]^T\). For modeling purposes, we assume that the tail moves from right to left at a constant speed \(U\)—information embedded in \(x(s,t)\)—and that the associated traveling wave travels from left to right, as shown in Fig. 2 (a). Under these assumptions, \(\boldsymbol{\mathcal{B}}\) is also an inertial frame. The tail material (AirTech A\(4000\)R\(14417\)) is assumed inextensible and, as a consequence, \(\left\|\partial{\boldsymbol{r}}/\partial{s}\right\|_2 = 1\) and \(\left[ds/dt \right]_{\text{tail}} = 0\), for all \((s,t)\). Thus, the velocity of a point on the centerline of the tail can be found as \[\begin{align} \boldsymbol{v}_{\text{b}}(s,t) = \frac{D \boldsymbol{r}}{Dt} = \frac{\partial \boldsymbol{r}}{\partial t} + \frac{\partial \boldsymbol{r}}{\partial s} \left[ \frac{ds}{dt} \right]_{\text{tail}} = \frac{\partial \boldsymbol{r}}{\partial t}. \label{EQ01} \end{align}\tag{1}\] Also, we define a local Lagrangian frame, \(\boldsymbol{\mathcal{L}}\), with its primary axis aligned with the unit tangent vector of the tail centerline, \(\boldsymbol{u}_{\text{t}}(s,t)\), and its secondary axis aligned with the unit normal vector of the tail centerline, \(\boldsymbol{u}_{\text{n}}(s,t)\), as shown in Fig. 2 (a). Using these definitions, (1 ) can be rewritten as \[\begin{align} \boldsymbol{v}_{\text{b}}(s,t) = v_{\text{b},\text{n}}(s,t) \boldsymbol{u}_{\text{n}}(s,t) + v_{\text{b},\text{t}}(s,t) \boldsymbol{u}_{\text{t}}(s,t). \label{EQ02} \end{align}\tag{2}\] Next, following Lighthill’s approach, we write the velocity of a parcel of fluid at \((s,t)\) in \(\boldsymbol{\mathcal{L}}\) coordinates as \[\begin{align} \boldsymbol{v}(s,t) = v_{\text{n}}(s,t) \boldsymbol{u}_{\text{n}}(s,t) + v_{\text{t}}(s,t) \boldsymbol{u}_{\text{t}}(s,t), \label{EQ03} \end{align}\tag{3}\] where, \(v_{\text{n}} = \boldsymbol{v} \cdot \boldsymbol{u}_{\text{n}}\) and \(v_{\text{t}} = \boldsymbol{v} \cdot \boldsymbol{u}_{\text{t}}\).

In this approach, we assume that the fluid is inviscid and that cannot penetrate the surface of the tail. From the no-penetration condition, it immediately follows that \(v_{\text{n}} = v_{\text{b},\text{n}}\); therefore, both \(v_{\text{n}}\) and \({\partial v_{\text{n}}}/{\partial t}\) can be estimated from the kinematics of the tail’s centerline. From the inviscid-fluid condition, it follows that, in general, \(v_{\text{t}} \neq v_{\text{b},\text{t}}\), which implies that the no-slip condition does not apply as a parcel of fluid can slide along the surface of the tail; therefore, \(v_{\text{t}} = \left[ds/dt\right]_{\text{fluid}}\) cannot be directly estimated from the kinematics of the tail’s centerline.

Reactive theory as presented in[25] is based on the notion of added mass, which is an idealization of the volume per unit length of incompressible fluid displaced when the modeled tail undulates while submerged in the surrounding medium. Here, we describe the added mass per unit length along the tail as \[\begin{align} m(s) = 2 \rho \beta h(s), \label{EQ04} \end{align}\tag{4}\] where \(\rho\) is the density of the fluid; \(\beta\) is the thickness of the tail material; and, \(h(s) = 0.694 \sqrt{l-s}\) defines the shape profile of the parabolic tail considered in this paper and shown in Fig. 2 (b). In this approach, resistive forces associated with fluid flow produced by viscosity are assumed negligible. Thus, as explained in the Appendix, we only consider the forces acting on the tail that arise from the reaction to accelerating the corresponding added mass. Namely, using Newton’s laws and including Lighthill’s assumptions, we obtain that the total instantaneous reactive force exerted on the undulating tail is given by \[\begin{align} \boldsymbol{F}_{\text{r}}(t) = -\int_0^l m(s) \left( \frac{\partial v_{\text{n}}}{\partial t} + v_{\text{t}} \frac{\partial v_{\text{n}}}{\partial s} \right) \boldsymbol{u}_{\text{n}}(s,t) ds. \label{EQ05} \end{align}\tag{5}\] Then, by extracting the component along the direction of locomotion, \(\boldsymbol{b}_1\), we obtain the magnitude of the signed instantaneous thrust force, \[\begin{align} F_{\text{th}} = -\int_0^l m(s) \left( \frac{\partial v_{\text{n}}}{\partial t} + v_{\text{t}} \frac{\partial v_{\text{n}}}{\partial s} \right) \boldsymbol{u}_{\text{n}} \cdot \boldsymbol{b}_1 ds, \label{EQ06} \end{align}\tag{6}\] and its corresponding time average, \[\begin{align} \bar{F}_{\text{th}} = -\int_0^l m(s) \left \langle \left( \frac{\partial v_{\text{n}}}{\partial t} + v_{\text{t}} \frac{\partial v_{\text{n}}}{\partial s} \right) u_{\text{n},1} \right \rangle ds, \label{EQ07} \end{align}\tag{7}\] where \(u_{\text{n},1} = \boldsymbol{u}_{\text{n}} \cdot \boldsymbol{b}_1\) and \(\langle \cdot \rangle\) denotes the time average of the corresponding function over an undulation cycle. The integrals specified by (6 ) and (7 ) can be separated into an unsteady term and a convective term. The unsteady term represents the reaction force exerted on the tail due to direct normal acceleration of the fluid parcel in contact with it; this term can be significant even when the undulating tail does not travel relative to the inertial frame—i.e., \(U=0\)—and the instantaneous curvature of the tail is small. In contrast, the convective term can be considered negligible for \(U=0\) and a small curvature; additionally, it is difficult to estimate \(v_{\text{t}}\) from video footage without the use of advanced techniques such as particle image velocimetry (PIV). Therefore, for the experiments and numerical estimations presented in Section [SECTION04], we consider that \(v_{\text{t}}=0\) because thrust is measured with the undulating tail attached to a static sensor.

As discussed in[25], over an undulation cycle—through a traveling wave like the one depicted in Fig. 2 (a)—the actuator delivers a finite amount of mechanical energy at the connection with the tail. Therefore, assuming that the elastic energy stored in the tail material is negligible, the cycle-averaged power carried by a traveling wave along the tail can be estimated as \[\begin{align} P_{\text{w}}(s) = \frac{1}{2}m(s)\left\langle v_{\text{n}}^2(s,t) \right\rangle v_{\text{w}}, \label{EQ08} \end{align}\tag{8}\] where \(v_{\text{w}}\) is the traveling velocity of the wave, which is determined by the actuation frequency. Using basic kinematics, we conclude that to maximize the stroke amplitude during operation, \(\left\langle v_{\text{n}}^2(s,t) \right\rangle\) must increase and, due to energy dissipation, \(P_{\text{w}}(s)\) is expected to decrease along \(s\). This observation indicates that a high-performance propulsion system of the type considered here must necessarily combine high-work-density (HWD) actuation and a soft-tail design with a decreasing \(h(s)\) along \(s\). As already mentioned, the model specified by (6 ) and (7 ) neglects resistive forces such as skin-friction and boundary-layer dissipation, which become more dominant as the Reynolds number decreases. Therefore, to describe force generation at low Reynolds numbers (\(\text{Re} < 1\)), resistive-force models are necessary. The propulsors studied in the presented research operate at Reynolds numbers in the range from about \(300\) to \(900\), and it is not immediately obvious that the approach presented here is suitable for this flow regime. To assess this and other issues, we developed the μN-resolution force sensor discussed in the next section.

3 A New  -Resolution Force Sensor↩︎

3.1 Dynamic Modeling for Design↩︎

From simple experiments, we estimate that the instantaneous hydrodynamic forces generated by the microswimmers shown in Fig. 1 must exhibit oscillatory components with sizes as small as \(6\) μN over the range of \(0\) to \(1.0\) mN. Therefore, the sensor we presented in[26]—with a resolution of \(7.48\) μN—is not well suited for the task at hand. This limitation prompted us to improve the resolution of existing μN-resolution force sensors for microrobotic applications. The resulting sensor presented here has a resolution of \(0.274\) μN and an operational range of \(-6.85\) to \(6.85\) mN. To develop this new device, we one more time followed the analysis and design principles for dual-cantilever structures (DCSs) discussed in[26], [27], and references therein. A lumped model for a structure of this type is depicted in Fig. 3 (b) (a), in which the target plate is necessary to measure the deformation of the distal end of the DCS using a capacitive sensor and, by convention, the upward and downward directions are considered to be positive and negative, respectively. According to the idealization shown in Fig. 3 (b) (a), a static mapping between the force value, \(F\)—defined as the signed magnitude of the corresponding vector force—and the displacement of the structure’s target plate, \(\delta\), follows Hooke’s law, \(\delta = F/k_{\text{tot}}\), where the elastic constant can be calculated as \(k_{\text{tot}} = 24EI/l_{\text{b}}^3\). In this description, \(E\) is the Young modulus of the structure’s material (Invar-\(36\)); \(l_{\text{b}}\) is the length of the DCS, as defined in Fig. 3 (b) (a); and, \(I\) is the moment of inertia of each composing beam’s cross-sectional area given by \(I = wd^3/12\), where \(w\) is the width of the beam, as defined in Fig. 3 (b) (a), and \(d\) is the thickness of the Invar-\(36\) material. Using basic beam theory, the DCS can be modeled as an equivalent mass-spring-damper system, whose response to a periodic excitation of the form \(F(t) = F_0 e^{j\omega_0t}\), with \(j = \sqrt{-1}\), has a magnitude given by \[\begin{align} \left|\delta(\omega_0) \right| = \frac{F_0}{k_{\text{tot}}\left[\eta^2 + (1-\omega_0^2/\omega_{\text{n}}^2)^2\right]^{\frac{1}{2}}}, \label{EQ09} \end{align}\tag{9}\] in which \(F_0\) is a known constant amplitude; \(\omega_0\) is a known constant frequency; \(\eta = 0.007\) is a hysteretic loss coefficient used to account for energy losses, provided by the manufacturer of the Invar-\(36\) material; and, \(\omega_{\text{n}}=\sqrt{k_{\text{tot}}/m_{\text{tot}}}\) is the undamped natural frequency of the equivalent dynamics of the structure, with the equivalent mass \(m_{\text{tot}}\). During development, through an iterative design process, we kept the sensor dimensions reported in[26] unchanged—i.e., \(w = 3\) mm, \(d=0.153\) mm—and varied the length, \(l_{\text{b}}\), which we finally selected to be \(16\) mm. The normalized magnitude response, \(\left| \delta(\omega_0)\right|/\delta(0)\), in decibels (dB) for the resulting design across the frequency set \(\left\{ \omega_0/2\pi \in \left[10^0, 10^3\right]\right\}\) Hz is shown in Fig. 3 (b) (b). As seen, the natural frequency of the system is \(288\) Hz, which is sufficiently high to obtain a wide useful measurement band for the desired application and flexibility to achieve a high resolution, precision, and accuracy. In this case, we predict a measurement error of less than \(1\)% for excitation frequencies below \(28.7\) Hz.

3.2 Fabrication of the Sensor↩︎

The fabrication process of the sensor, which consists of three steps, is graphically explained in Fig. 3 (a). In Step \(1\), the Invar-\(36\) pieces of the DCS are premachined using a \(3\)-W UV laser (Photonics Industries DCH-\(355\)-\(3\)), and folded along perforated lines to create \(3\)D structures; in Step \(2\), following the clip-brazing method introduced in[26], small clips of brazing wire—with a \(80/16/4\) wt.% Ag-Cu-Zn composition—are placed at the interfaces between pieces and folding lines, and then melted with a butane torch, to create a stable rigid monolithic structure; lastly, in Step \(3\), the \(1\)-cm-diameter target plate, also made of Invar-\(36\), depicted in Fig. 3 (b) (a) is mirror polished and connected to the DCS as depicted in Fig. 3 (b) (a), using conductive silver epoxy (MG Chemicals \(8331\)D). After fabrication, to prevent oxidation, we mechanically polished and coated the structure with oil (WD-\(40\)). As seen in Fig. 3 (b), before performing a force characterization experiment, the tested propulsor is glued to the lower distal end of the DCS, and the target plate is precisely aligned beneath the cylindrical probe of the displacement capacitive sensor (Physik Instrumente D-\(510.021\)) employed to measure the instantaneous value of \(\delta\).

3.3 Calibration of the Sensor↩︎

Before operation, we calibrated the sensor using the same procedure described in [26]. Under this methodology, as depicted in Fig. 3 (b), a series of custom-made laboratory weights are sequentially hung from the DCS with a tested propulsor attached to it in order to find a relationship between weight and the voltage outputted by the capacitive displacement sensor. To mitigate the chances of human error, during each calibration test, we hang the corresponding laboratory weight almost instantaneously using a linear stage. This action smoothly excites the second-order dynamics of the DCS with a Heaviside-like force, thus generating a step response in the sensor reading. Two examples of such responses are shown in Fig. 3 (c); here, the signals in blue and red correspond to calibration weights with masses of \(10.9\) mg and \(49.5\) mg, respectively. The resulting calibration mapping used to interpret the data obtained through the experiments presented in Section [SECTION04] is shown in Fig. 3 (d). In this plot, each red data point indicates the mean and experimental standard deviation (ESD) of the steady-state values corresponding to the ten step responses obtained through ten back-to-back experiments performed for each calibration weight in the set \(\left\{0.107,0.2011,0.299,0.393,0.486\right\}\) mN. Also, in Fig. 3 (d), the blue line represents the least-squares fitting of the data. This linear function maps the probe voltage generated by the capacitive displacement sensor, \(V_{\text{p}}\), to the sensed force, \(F_{\text{s}}\), according to \(V_{\text{p}} = 1.46\cdot F_{\text{s}}\), with a coefficient of determination of \(0.999\). We performed the calibration procedure several times with both the single-tail and dual-tail propulsors connected to the sensor as depicted in Fig. 3 (b), and we obtained almost identical results. For this reason, we used the same mapping in Fig. 3 (b) to interpret all the force data discussed in Section [SECTION04].

4 Experiments for Force Characterization↩︎

4.1 Experimental Setup and Data Processing↩︎

The experimental setup used for characterizing the forces produced by both the single-tail and dual-tail tested propulsors is shown in Fig. 4 (a). Here, a Mathworks Simulink Real-Time host–target system, equipped with an analog-digital/digital-analog (AD/DA) board (National Instruments PCI-\(6229\)), is used to generate, process, and record signals at a rate of \(5\) kHz. During operation, the target computer generates the pulse-width modulation (PWM) signal required for actuation, with prescribed frequency (\(f\)) and duty-cycle (DC) parameters. This signal is power-amplified by a MOSFET-based circuit (YYNMOS-\(4\)), which opens and closes the conductive pathways from a power supply to the actuation mechanism of the tested propulsor. To obtain the desired amplitudes of actuation, the voltage of the power supply is adjusted according to the procedure discussed in[7]. As seen in the inset of Fig. 4 (a), during a force characterization experiment, the tested propulsor remains submerged from just above the leading edge of its soft tail in a pool filled with water. During the experiments, we tested both studied propulsors at all the combined \(f\)–DC pairs in the sets \(f\in\{1:1:4\}\) Hz and \(\text{DC}\in\{1:10\}\) %. To ensure similar operational tail strokes across different experiments and devices, we tuned the on-height voltages of the exciting PWM signals such that the current flowing through the SMA wires was \(250\) mA in all the cases. The photographs in Figs. 4 (b) and (c) show the single-tail and dual-tail tested propulsors connected to the force sensor, before operation.

During operation, the voltages outputted by the capacitive displacement sensor are filtered through a signal conditioner (Physik Instrumente E-\(852\) PISeca) before being sampled by the AD/DA board of the host–target system and recorded. For the data presented in this section, after collection, the signals were digitally processed to remove high-frequency electrical noise as well as low-frequency sensor drifting possibly caused by changes in the environmental temperature and humidity. More specifically, the raw data was filtered through a low-pass finite-impulse response (FIR) zero-phase digital filter of order \(5 \times 10^3\) with a cutoff frequency of \(0.01\) Hz to extract the low-frequency drift. This extracted drift was then subtracted from the raw signal to obtain a drift-free signal. This drift-free signal was then filtered through another low-pass FIR zero-phase digital filter of order \(2\times10^3\) with a cutoff frequency of \(50\) Hz to cancel high-frequency noise, which yields the final filtered data used for analysis. We define the cycle-peak force, \(F_{\text{p}}\), as the maximum force value produced by a propulsor during an actuation cycle, and the cycle-averaged force, \(F_{\text{a}}\), as the average of the force values measured over the same actuation cycle. For each independent force test, we used five steady-state actuation cycles to compute the mean of the corresponding five \(F_{\text{p}}\) values, which we denote as \(\bar{F}_{\text{p}}\); similarly, for each independent force test, we used five steady-state actuation cycles to compute the mean of the corresponding five \(F_{\text{a}}\) values, which we denote as \(\bar{F}_{\text{a}}\). For the purposes of data analysis, for each \(f\)–DC pair, we performed five independent force tests and we computed the means and ESDs of the \(\bar{F}_{\text{p}}\) and \(\bar{F}_{\text{a}}\) values.

4.2 Results and Discussion↩︎

Figs. 5 (a) and (b) show the measured instantaneous forces respectively generated by the single-tail and dual-tail tested propulsors, using a \(1\)-Hz PWM actuation signal with a DC of \(5\)%. In these cases, it can be observed that the dual-tail propulsor generates noticeably higher peaks of positive force, while the single-tail propulsor generates deeper valleys of negative force. These two phenomena combined indicate that the average force in the single-tail case must be lower than in the dual-tail case. Figs. 5 (c) and (d) show the instantaneous forces measured for the single-tail and dual-tail tested propulsors, using a \(2\)-Hz PWM actuation signal with a DC of \(10\)%. In these cases, it can be observed that the dual-tail propulsor generates drastically higher peaks of positive force, and both propulsors generate visually comparable valleys of negative force. These observations indicate that the average force in the single-tail case must be considerably lower than in the dual-tail case. Figs. 5 (e) and (f) show sequences of video frames, taken at \(0.1\)-s intervals, of full actuation cycles corresponding to the first second of experimental data presented in Figs. 5 (a) and (b), respectively.

The data obtained through the force characterization experiments are summarized in Fig. 6. Here, each data point in Fig. 6 (a) indicates the mean and ESD of the \(\bar{F}_{\text{p}}\) values corresponding to five back-to-back force tests performed using the single-tail propulsor and a given \(f\)–DC pair; and, each data point in Fig. 6 (b) indicates the mean and ESD of the \(\bar{F}_{\text{a}}\) values corresponding to the same five back-to-back force tests associated with the data shown in Fig. 6 (a). Similarly, each data point in Fig. 6 (c) indicates the mean and ESD of the \(\bar{F}_{\text{p}}\) values corresponding to five back-to-back force tests performed using the dual-tail propulsor and a given \(f\)–DC pair; and, each data point in Fig. 6 (d) indicates the mean and ESD of the \(\bar{F}_{\text{a}}\) values corresponding to the same five back-to-back force tests associated with the data shown in Fig. 6 (c). These data clearly show that, in all single-tail cases, the means of the \(\bar{F}_{\text{p}}\) values obtained with actuation frequencies of \(3\) and \(4\) Hz are significantly lower than those measured at frequencies of \(1\) and \(2\) Hz; while, in most but not all single-tail cases, the means of the \(\bar{F}_{\text{a}}\) values at \(1\) Hz are noticeably higher than those corresponding to the other three actuation frequencies. By comparison, in all the dual-tail cases, both the means of the \(\bar{F}_{\text{p}}\) and \(\bar{F}_{\text{a}}\) values obtained with actuation frequencies of \(3\) and \(4\) Hz are significantly lower than those measured at frequencies of \(1\) and \(2\) Hz.

Furthermore, it can be observed that while the means of the \(\bar{F}_{\text{p}}\) values are generally comparable between the single-tail and dual-tail cases, the means of the \(\bar{F}_{\text{a}}\) values are significantly larger in the dual-tail case. As seen, among all single-tail cases, the maximum means of the \(\bar{F}_{\text{p}}\) and \(\bar{F}_{\text{a}}\) values both occur at a frequency of \(1\) Hz and DC of \(5\)%, which were measured to be \(0.45\) mN and \(2.97\) μN, respectively. By comparison, among the dual-tail cases, the maximum mean of the \(\bar{F}_{\text{p}}\) values—\(0.61\) mN—occurs at an actuation frequency of \(1\) Hz and DC of \(5\)%, while the maximum mean of the \(\bar{F}_{\text{a}}\) values—\(22.6\) μN—occurs at an actuation frequency of \(1\) Hz and DC of \(7\)%. Additionally, it is interesting to note that, at an actuation frequency of \(2\) Hz, the dual-tail propulsor performs drastically better than its single-tail counterpart; specifically, among all the \(2\)-Hz cases, the maximum mean of the \(\bar{F}_{\text{p}}\) values for the dual-tail propulsor—\(0.56\) mN, occurring at a DC of \(9\)%—is \(4.31\) times larger than that corresponding to the single-tail at a DC of \(9\)%.

Similarly, among all the \(2\)-Hz cases, the maximum mean of the \(\bar{F}_{\text{a}}\) values for the dual-tail propulsor—\(14.08\) μN, occurring at a DC of \(10\)%—is positive and \(37.01\) times larger than that corresponding to the single-tail at a DC of \(10\)%, which is, additionally, negative. The data presented in Fig. 6 provide compelling evidence demonstrating that the tested dual-tail propulsor produces considerably higher average thrusts than its single-tail counterpart, which explains to a significant extent the improved locomotion performance achieved by the dual-tail VLEIBot⁺ and VLEIBot⁺⁺ swimmers compared to that measured for the single-tail VLEIBot swimmer. Note, however, that the presented force tests do not account for a periodic lateral motion experienced by the single-tail VLEIBot during swimming, a phenomenon likely produced by the periodic reactive torque acting on the robot’s head as the tail flaps in the surrounding fluid. Furthermore, since the front rigid parts of the tested propulsors are held fixed relative to the inertial frame of reference, the forces measured according to the method graphically described in Fig. 4 are not fully representative of the total thrust produced by a swimmer moving forward at a constant speed \(U \neq 0\).

Last, we evaluate the predictive power and limitations of the model specified by (6 ) with \(v_{\text{t}}=0\), using a custom algorithm implemented in MATLAB that extracts the curvature of the tested tail from video footage. Specifically, this algorithm reads each video frame, converts it into a binary image, removes unwanted features, and then extracts the instantaneous centerline of the tested tail. The extracted instantaneous curvature of the tail is then used to estimate the kinematic variables required by the model specified by (6 ), assuming \(v_{\text{t}}=0\). Figs. 7 (a) and (b) compare the instantaneous thrust forces obtained through direct measurement and model-based estimation corresponding to a \(1\)-Hz PWM excitation with a DC of \(5\)%, for the single-tail and dual-tail propulsors, respectively. It can be clearly observed in Fig. 7 (a) that, for the single-tail case, the model-based estimated and measured forces are similar in shape and comparable in magnitude. However, the model predicts a slight positive force after the initial spike; in contrast, the measured signal exhibits a negative force valley during the same interval of the actuation cycle. In the dual-tail case, it can be observed that the amplitude of the model-based estimate is considerably lower than that measured experimentally. We hypothesize that this discrepancy was partly caused by the introduction of noise during the image-processing-based extraction of the tail’s curvature from video footage. We believe, however, that the most important reason for the underestimation of the model-based thrust is that the operator specified by (6 ), with \(v_{\text{t}}=0\), neglects the convective contribution associated with the \(v_{\text{t}} {\partial u_{\text{n}}}/{\partial s}\) term. To see the logic behind this argument, note that, in the dual-tail case, \(v_{\text{t}}\) is not expected to be negligible even though \(U=0\) because of interactions between the wake structures generated by the two undulating tails of the propulsors. It is well known that when wake structures interact, constructive or destructive patterns form, which might enhance or attenuate the process of thrust production[28]. Besides these limitations, the results presented in this paper are significant as they suggest that direct measurement of thrust through static experiments can be leveraged to design single-tail propulsors. Nonetheless, the experimental data and associated analyses highlight the need for further investigation into the modeling of the involved hydrodynamic forces to account for multi-tail wake–structure interaction and other convective phenomena.

5 Conclusion↩︎

We presented a study of the characteristics of the thrust generated by two newly developed insect-scale propulsors that leverage the FSI between unactuated soft tails and their surrounding fluid to function. The first actuator studied is single-tailed, and the second is double-tailed. In this approach, a driving actuator generates a wave that travels along a tail to produce hydrodynamic forces useful for locomotion. To describe the mechanism of thrust generation, we performed an analysis based on reactive theory that provides insights for the design and development of high-performance insect-scale propulsors. To experimentally measure the instantaneous thrust generated by the two studied devices during static tests, we designed and fabricated a high-sensitivity μN-resolution force sensor. Employing this sensor, we discovered that the tested single-tail propulsor generates relatively small cycle-averaged force values, with the best case being \(2.97\) μN, which occurs at low actuation frequencies (\(\leq 2\) Hz); by comparison, we measured cycle-averaged force values of up to \(22.6\) μN for the dual-tail propulsor.

Interestingly, we also found that, at respectively \(0.45\) mN and \(0.61\) mN, the measured peak force values for the single-tail and dual-tail cases are very similar. Considering these data, we hypothesize that the dual-tail propulsor generates a significantly higher average thrust due to the interaction between the wake structures created by its two composing soft tails. This hypothesis is consistent with the finding that the reactive-force model presented in Section [SECTION02] estimates force values similar to those measured for the single-tail case; in contrast, a noticeable mismatch between the modeled and experimentally obtained propulsion forces occurs in the dual-tail case. This discrepancy is most likely due to the omission of the convective term in the formulation of (6 ), which seems to be necessary to account for interaction between wake structures. In the near future, we will refine the reactive-force model presented in Section [SECTION02] to account for wake-structure phenomena. Aligned with this objective, we will improve our understanding of the dynamics relating actuation output with thrust for locomotion, which is essential for developing advanced control schemes for insect-scale aquatic robots.

Appendix: Reactive-Force Model↩︎

In this appendix, we review the fundamental ideas introduced by Lighthill in[25] to model the total reactive force acting on a slender undulating tail. Here, we start by finding the material acceleration of a parcel of fluid at \(\left(s,t\right)\) as \[\begin{align} \boldsymbol{a}(s,t) = \frac{D\boldsymbol{v}}{Dt} = a_{\text{n}} \boldsymbol{u}_{\text{n}} + a_{\text{t}} \boldsymbol{u}_{\text{t}}, \label{EQ10} \end{align}\tag{10}\] with \[\begin{align} \begin{aligned} a_{\text{n}}(s,t) &= \frac{Dv_{\text{n}}}{Dt} + v_{\text{t}} \frac{\partial\theta}{\partial t} + v_{\text{t}}^2 \frac{\partial \theta}{\partial s}, \\ a_{\text{t}}(s,t) &= \frac{Dv_{\text{t}}}{Dt} - v_{\text{n}} \frac{\partial\theta}{\partial t} - v_{\text{n}} v_{\text{t}} \frac{\partial \theta}{\partial s}, \end{aligned} \label{EQ11} \end{align}\tag{11}\] where \(\theta(s,t)\) is the local tangent angle of the tail centerline from \(\boldsymbol{b}_1\). Next, we invoke Newton’s second and third laws to obtain the hydrodynamic force per unit length exerted by the fluid on the tail, \[\begin{align} \begin{aligned}\boldsymbol{f}(s,t) &= -m(s) \boldsymbol{a}(s,t) \\ &= -m(s) \left[\left( \frac{Dv_{\text{n}}}{Dt} + v_{\text{t}} \frac{\partial\theta}{\partial t} + v_{\text{t}}^2 \frac{\partial \theta}{\partial s} \right) \boldsymbol{u}_{\text{n}} \right.\\ & \left.+\left( \frac{Dv_{\text{t}}}{Dt} - v_{\text{n}} \frac{\partial\theta}{\partial t} - v_{\text{n}} v_{\text{t}} \frac{\partial \theta}{\partial s} \right) \boldsymbol{u}_{\text{t}} \right], \end{aligned} \label{EQ12} \end{align}\tag{12}\] where \(m(s)\) is the added mass per unit length, which in this case is specified by (4 ). Additionally, we consider that (i) the motion of the tail’s centerline is confined to the \(\boldsymbol{b}_1\)–\(\boldsymbol{b}_2\); (ii) \(\boldsymbol{\mathcal{B}}\) does not rotate relative to \(\boldsymbol{\mathcal{N}}\); (iii) the local curvature defined by \(\theta(s,t)\), and its variations over \(s\) and \(t\) are small; and, following Lighthill’s assumptions, (iv) reactive forces arise from the normal acceleration of the displaced fluid mass and act perpendicularly on the tail surface. Directly from these assumptions it follows that the reactive force per unit length is given by \[\begin{align} \boldsymbol{f}_{\text{r}}(s,t) = -m(s) \boldsymbol{a}_{\text{n}}(s,t) = -m(s) \frac{Dv_{\text{n}}}{Dt} \boldsymbol{u}_{\text{n}}(s,t), \label{EQ13} \end{align}\tag{13}\] in which the material derivative of the normal velocity component, \(v_{\text{n}}\), satisfies \[\begin{align} \frac{Dv_{\text{n}}}{Dt} = \frac{\partial v_{\text{n}}}{\partial t} + \frac{ds}{dt} \frac{\partial v_{\text{n}}}{\partial s} =\frac{\partial v_{\text{n}}}{\partial t} + v_{\text{t}} \frac{\partial v_{\text{n}}}{\partial s}. \label{EQ14} \end{align}\tag{14}\] Therefore, \[\begin{align} \boldsymbol{f}_{\text{r}}(s,t) = -m(s) \left( \frac{\partial v_{\text{n}}}{\partial t} + v_{\text{t}} \frac{\partial v_{\text{n}}}{\partial s} \right) \boldsymbol{u}_{\text{n}}(s,t), \label{EQ15} \end{align}\tag{15}\] which allows us to compute the total instantaneous reactive force acting on the tail as specified by (5 ).

References↩︎