Razor clam to RoboClam: burrowing drag reduction mechanisms and their robotic adaptation

The adage 'clear as mud' is often used to describe the difficulty of visually investigating burrowing animals in situ. To explore the soil mechanics at play during E. directus burrowing, we constructed a 2D ant farm, or Hele-Shaw cell [19], to see through the substrate surrounding a digging animal and measure deformation of the granular medium with particle image velocimetry (PIV) [9]. This work showed that the uplift and contraction motions of E. directus' valves agitate the surrounding soil (figure 2(A)), creating a state of localized fluidization. This event occurs at a faster timescale than that required for the soil to naturally fail and landslide toward the animal (figure 2(B)).

The discontinuity at the failure surface, shown by θ_f in figure 2(B), enables fluidization to occur—as the valves contract beyond the point when incipient failure is induced in the soil, the particles in the failure zone are free to move with the pore fluid while the particles outside the failure zone remain stationary. The relevant Reynolds number of the pore fluid flow, $\mathrm{Re} = \frac{\rho _f v_v d_p}{\mu _f}$, calculated from E. directus' valve velocity v_v (determined from the valve contraction angular velocity and valve width), particle diameter d_p, and the pore fluid density ρ_f and viscosity μ_f, varies between 0.02 and 56, depending on particle size (0.002–2 mm [6, 8, 10]), animal size (10–20 cm, from experimental observation), and valve contraction velocity (v_v ≈ 0.011–0.028 m s⁻¹ [8]). As this range of Reynolds numbers falls primarily within the regime of Stokes drag (the transition to form drag occurs at a Reynolds number of approximately 100) [20], the characteristic time for a particle to reach the pore fluid velocity can be estimated through conservation of momentum:

$$\begin{equation} m_p \frac{{\rm d} v_p}{{\rm d} t} = 6 \pi \mu _f d_p (v_v - v_p) \rightarrow t_{{\rm char}} = \frac{{d_p}^2\rho _p}{36\mu _f}, \end{equation} \tag{ 1 }$$

where m_p is the mass of a soil particle, $\frac{{\rm d} v_p}{{\rm d} t}$ is the acceleration of the particle, and ρ_p is the density of the particle. For the 1 mm soda lime glass beads used in our experiments, t_char = 0.075 s. This timescale is considerably less than the ≈ 0.2 s valve contraction time measured by Trueman [8], indicating that soil particles surrounding E. directus can be considered inertialess and are advected with the pore fluid during valve contraction.

The discontinuity created by the failure surface is critical to achieving fluidization, as without it substrate particles would follow the fluid flow field, which is incompressible and governed by ${\nabla }{\cdot }\overline{v}=0$. No divergence in the flow field creates no divergence between particles, and thus no unpacking. However, in the presence of a finite failure zone, E. directus' contraction motion reduces the volume between the valves, which draws pore water into the region surrounding the animal. This pore water mixes with the failed soil to create a locally unpacked, fluidized zone.

The previous subsection relates to our measurements of E. directus in a two-dimensional experimental setup; the following discussion is focused on describing the mechanics that govern how the failure, and fluidized, zone forms around the animal in three dimensions.

As E. directus initiates contraction of its valves, it reduces the level of stress acting between the valves and the surrounding substrate, causing failure in the soil. Figure 3(A) shows a Mohr's circle representation [21] of the effective stress states at equilibrium, before contraction (circle a), and during the initiation of contraction, which brings the soil to one of two states of active failure, caused either by an imbalance between radial and vertical stresses (circle b), or an imbalance between radial and hoop stresses (circle c). Effective stress is the actual stress acting between soil particles, neglecting hydrostatic pressure of the pore fluid, and is denoted in this paper with a prime. The term 'active' corresponds to the reduction (rather than increase) of one of the principal stresses to induce failure [6].

Zoom In Zoom Out Reset image size

As the soil begins to fail, it will tend to naturally landslide downward at a failure angle θ_f. At this point, the shear stresses in the soil are equal to its shear strength. This condition is shown in figure 3(A), with the applied stress circles b and c tangent to the failure envelope, which lies at the same angle as the friction angle of the soil φ, a property commonly measured during a geotechnical survey. The failure angle is the transformation angle between the principal stress state and the stress state at failure. This angle can also be determined by connecting the tangency point on the failure envelope, the horizontal effective stress at failure $\sigma ^\prime _{hf}$, and the principal stress axis (figure 3(A)), and is given by

$$\begin{equation} \theta _f= \frac{\pi }{4}+\frac{\varphi }{2}. \end{equation} \tag{ 2 }$$

Equation (2) was used to plot the failure angle in figure 2(B), with the friction angle of the substrate measured as 25°.

To describe soil failure in three dimensions, we propose a simplified model of E. directus as a cylinder with contracting radius that is embedded in saturated soil (figure 3(B)). To neglect end effects, the length of the cylinder is considered to be much larger than its radius. When E. directus initiates valve contraction, it induces incipient failure without moving the substrate particles; as this relaxation in pressure can be considered quasi-static and elastic [6], stresses due to inertial effects can be ignored and the total radial and hoop stress distribution in the substrate can be described with the following thick-walled pressure vessel equations [22], which have been modified to geotechnical conventions (with compressive stresses positive) and to reflect an infinite soil bed in lateral directions:

$$\begin{equation} \sigma _r = \frac{R_0^2 (p_i - p_0)}{r^2} + p_0 \end{equation} \tag{ 3 }$$

$$\begin{equation} \sigma _\theta = -\frac{R_0^2 (p_i - p_0)}{r^2} + p_0, \end{equation} \tag{ 4 }$$

where σ_r is total radial stress, σ_θ is total hoop stress, R₀ is E. directus' expanded radius (before contraction), p_i is the pressure acting on the animal's valves, and p₀ is the natural lateral equilibrium pressure in the soil. It is important to note that these equations still hold if there is a body force acting in the z-direction, such as in soil. In this case, the pressure vessel equations describe the state of stress within annular differential elements stacked in the z-direction. The total vertical stress is given as

$$\begin{equation} \sigma _z = \rho _t g h, \end{equation} \tag{ 5 }$$

where h is the clam's depth beneath the surface of the soil, ρ_t is the total density of the substrate (including solids and fluids), and g is the gravitational constant. It should be noted that there are no shear stresses within the soil in principal orientation, as τ_rz = τ_θz = 0 because E. directus is modeled with a high aspect ratio (L ≫ R₀) and τ_rθ = 0 because of symmetrical radial contraction.

The undisturbed horizontal effective stress in the substrate is determined by subtracting hydrostatic pore pressure u from the natural lateral equilibrium pressure:

$$\begin{equation} \sigma ^{\prime }_{h0} = p_0-u. \end{equation} \tag{ 6 }$$

The undisturbed horizontal and vertical effective stresses can be correlated through the coefficient of lateral earth pressure

$$\begin{equation} K_0 = \frac{\sigma ^{\prime }_{h0} }{\sigma ^{\prime }_{v0}}, \end{equation} \tag{ 7 }$$

which is a soil property that can be measured through geotechnical surveys [6, 23]. By also knowing the void fraction of the soil and the particle and fluid density, ρ_p and ρ_f respectively, p₀ can be determined as

$$\begin{equation} p_0 = K_0 \sigma ^{\prime }_{v0} + u = K_0 g h (1 - \epsilon ) (\rho _p - \rho _f) + \rho _f g h. \end{equation} \tag{ 8 }$$

Failure of the substrate will occur when p_i is lowered to a point where the imbalance of two principal effective stresses produces a resolved shear stress that equals the shear strength of the soil. This resolved failure shear stress can be created by an imbalance between radial and vertical stresses (figure 3(A), circle b) or radial and hoop stresses (figure 3(A), circle c). From the geometry of either circle and the failure envelope defined by φ, the relationship between stresses at failure for either mechanism is

$$\begin{equation} \frac{\sigma ^{\prime }_{rf}}{\sigma ^{\prime }_{vf}} = \frac{\sigma ^{\prime }_{rf}}{\sigma ^{\prime }_{\theta f}} = \frac{1 - \sin \varphi }{1 + \sin \varphi } = K_a, \end{equation} \tag{ 9 }$$

where the subscript f denotes the stresses at failure and K_a is referred to as the coefficient of active failure. It is important to note that this failure analysis is also valid for cohesive soils. The difference between cohesive and granular soils when plotted on a Mohr's circle is that the failure envelope does not pass through (0,0), as cohesive stresses give soil shear strength even when no compressive stresses are applied. At sufficient depths the failure envelope can be approximated as running through (0,0) for any soil type, as compressive stresses due to gravity will dominate cohesive stresses.

Soil failure due to an imbalance between radial and vertical stresses will occur when the applied radial effective stress equals the radial stress at failure. The radial location of the failure surface in this condition, $R_{f_{rv}}$, can be found by combining (3) for radial stress with (6), (7), and (9), and realizing that the vertical effective stress at failure and equilibrium is unchanged, namely

$$\begin{eqnarray*} \sigma ^{\prime }_r \big |_{r = R_{f_{rv}}} &=& \sigma ^{\prime }_{rf} \nonumber \\ \frac{R_0^2 (p_i - p_0)}{R_{f_{rv}}^2} + p_0 - u &=& \frac{K_a}{K_0} (p_0 - u) \nonumber \end{eqnarray*}$$

yielding the dimensionless radius for radial–vertical stress imbalance-induced failure:

$$\begin{equation} \frac{R_{f_{rv}}}{R_0} = \left[ \frac{p_i - p_0}{\big(\frac{K_a}{K_0} - 1 \big) (p_0 - u)} \right]^\frac{1}{2}. \end{equation} \tag{ 10 }$$

If soil failure is caused by an imbalance between radial and hoop stresses, the radial location of the failure surface, $R_{f_{r\theta }}$, can be found by combining (3) and (4) with (6) and (9):

$$\begin{eqnarray*} \sigma ^{\prime }_r \big |_{r = R_{f_{r \theta }}} &=& \sigma ^{\prime }_{rf} \nonumber \\ \frac{R_0^2 (p_i - p_0)}{R_{f_{r \theta }}^2} + p_0 - u &=& K_a \left[ -\frac{R_0^2 (p_i - p_0)}{R_{f_{r \theta }}^2} + p_0 - u \right] \nonumber \end{eqnarray*}$$

yielding the dimensionless radius for radial–hoop stress imbalance-induced failure:

$$\begin{equation} \frac{R_{f_{r \theta }}}{R_0} = \left[ \frac{(K_a + 1) (p_i - p_0)}{(K_a - 1) (p_0 - u)} \right]^\frac{1}{2}. \end{equation} \tag{ 11 }$$

The dominant failure mechanism in the soil surrounding a contracting cylindrical body is determined by the type of failure (radial–vertical or radial–hoop) that results in the largest failure surface radius. The ratio of failure radii for both mechanisms can be calculated by combining (10) and (11) into

$$\begin{equation} \frac{R_{f_{r v}}}{R_{f_{r \theta }}} = \left[ \frac{K_a - 1}{(K_a + 1) \big(\frac{K_a}{K_0} - 1 \big)} \right]^\frac{1}{2}. \end{equation} \tag{ 12 }$$

Figure 4(A) shows values of (12) for the full range of possible soil types, with 0.19 ⩽ K_a ⩽ 0.52 and 0.31 ⩽ K₀ ⩽ 1 [6, 23]. For the region of radial–hoop stress-dominated failure (above the dashed line in figure 4(A)), $\frac{R_{f_{r v}}}{R_{f_{r \theta }}} \approx 1$, indicating that the failure radius from both mechanisms is approximately the same size. For the region of radial–vertical stress-dominated failure (below the dashed line in figure 4(A)), $\frac{R_{f_{r v}}}{R_{f_{r \theta }}} > 1$ and rapidly rises to an asymptote for lower values of K₀ and higher values of K_a. The asymptote and the imaginary region in figure 4(A) correspond to natural stress imbalances that exceed the shear strength of the soil and are thus impossible. Since $R_{f_{r v}} \gtrapprox R_{f_{r \theta }}$, radial–vertical stress imbalance-induced failure is considered the dominant mechanism and is used to predict the failure zone size from the following scaling arguments.

Zoom In Zoom Out Reset image size

If during contraction p_i is assumed to be approximately zero, corresponding to complete stress release between E. directus' valves and the surrounding soil, and u ≈ 0.5p₀ because of the relative densities and volumetric fractions of the soil particles and interstitial water, (10) reduces to

$$\begin{equation} \frac{R_f}{R_0} \approx \left(\frac{2}{1-\frac{K_a}{K_0}}\right)^{\frac{1}{2}}. \end{equation} \tag{ 13 }$$

This facilitates a prediction of R_f using only two soil properties, K_a and K₀, both of which are commonly measured during a geotechnical survey [24].

Applying the range of possible K_a and K₀ values to (13) yields $1 < \frac{R_f}{R_0} < 5$ in most conditions (figure 4(B)). These results demonstrate that soil failure around a contracting cylindrical body is a relatively local effect, and for reductions of p_i to near zero, depth-independent. Equation (13) also does not depend on any soil cohesion terms, indicating that localized substrate failure and fluidization should be possible in both granular and cohesive soils.

To explore the transfer of localized fluidization burrowing to engineering applications, we designed and built RoboClam, a robot that replicates the digging kinematics of E. directus (figure 5(A)). The robot consists of a control platform with two pneumatic pistons that actuate an end effector in the same up-in-down-out motions as E. directus' valves (figures 5(B)–(F)). The end effector is the same size scale as E. directus (9.97 cm long and 1.52 cm wide) and is capable of the same contraction displacement of an adult animal (∼6.4 mm). Unlike E. directus, which uses its foot—a soft, dexterous organ—to move its valves up and down, RoboClam uses a simpler mechanical arrangement to actuate the end effector up and down with a pneumatic piston, which is capable of a much longer stroke than the animal's foot. Pneumatics were chosen as the main power source so RoboClam can be safely tested in real, undisturbed ocean substrates, which enables accurate comparisons to E. directus performance and mitigates container wall effects.

Zoom In Zoom Out Reset image size

The end effector is sealed within a rubber boot to prevent soil particles from jamming the valve expansion/contraction mechanism (figure 5(G)). Displacement of the mechanism is accomplished with a sliding wedge that moves the two valves of the end effector in and out. The wedge and valves are exactly constrained (figure 5(H)), with contact lengths/widths greater than two in order to prohibit jamming during any part of the stroke [25]. The wedge was designed to intersect the center of pressure of the valves regardless of its position, assuming the center of pressure is located approximately at the center of the valves. The geometry and exact constraint of the end effector mechanism allow its efficiency to be characterized by measuring the coefficient of friction between the wedge and valves, which was found to be 39% [26].

As RoboClam burrows, the control system tracks the total amount of energy input to the system by integrating the forces acting on the pneumatic actuators over their displacement. Because we were able to characterize the frictional losses in the end effector and the actuators, as well as account for potential energy changes, the energy expended for soil deformation while burrowing could be tracked [26]. This is an important facet of the machine's design, as overall energy consumed is device-dependent; the purpose of creating RoboClam was to test the new method of using a machine to burrow via localized fluidization. After the method is understood and characterized, other burrowing machines can be designed for optimized efficiency.

RoboClam is controlled using a GA. The control space for the robot is composed of eight independent parameters: up/in/down/out pneumatic pressure and the time or displacement associated with each motion. A GA was used because it can handle many independent variables in an optimization problem. The random mutation and recombination of traits used by a GA may allow it to find a global minimum, even in situations in which other optimization methods would not [27]. During the ocean tests reported in this work, MATLAB's⁴ built-in GA was used, with a population of 10–20 individuals running for 10–20 generations. In laboratory tests reported, RoboClam ran from customized GA software written in Python⁵.

A GA is analogous to evolution in nature [27]. In the case of RoboClam, the GA begins a test by randomly generating a population of individuals, where each individual is a set of instructions (movement times/displacements and pressures) to run the machine. The GA runs the robot with each individual and records its digging performance. Individuals are then compared by a metric to be optimized—a 'cost' function that relates to digging efficiency for RoboClam. Individuals in the population with the lowest cost are allowed to interbreed (mix parameters) with each other, high-cost individuals are killed off, and new individuals are added to the population to form the next generation. The process repeats for many generations, ideally continually decreasing the cost to a global minimum, which appears as a cost asymptote.

In quantifying optimized burrowing efficiency, two factors are relevant: the overall energy expenditure per depth of the robot β, where $\frac{E}{\delta } = \beta$; and the power law relationship α between energy expended and depth, where $\alpha = \frac{\ln E}{\ln \delta }$ with the energy-depth data centered about (0, 0). Fluidized soil exhibits α = 1, whereas static soil exhibits α = 2. The cost function used for the first set of tests reported in this work (figure 6) was

$$\begin{equation} {\rm cost}_{{\rm ocean}} = \beta \alpha . \end{equation} \tag{ 14 }$$

These tests were conducted in real E. directus habitat off the coast of Gloucester, MA. One hundred and twenty five tests were run during the trial, representing 13 generations of the GA (with the final generation incomplete). The controlled digging parameters were up and down time, in and out displacement, and the pressures associated with each movement. The robot was moved to a new location for each test, greater than ten characteristic lengths away from the previous test, as to maintain virgin soil conditions.

Zoom In Zoom Out Reset image size

We found that low values of β resulted in low-energy burrowing for depths of 20–30 cm, but with relatively high α (≈2, indicating non-fluidized substrate). At greater depths these burrowing techniques would not provide any energetic advantage over existing techniques used to penetrate static soil [7]. As a result, the cost function used for the second set of data reported in this work (figure 7) was

$$\begin{equation} {\rm cost}_{{\rm lab}} = \alpha , \end{equation} \tag{ 15 }$$

as α ≈ 1 demonstrates localized fluidization—the aim of our experiments. These tests were conducted in 1 mm soda lime glass beads saturated with water, the same substrate used in our visualization experiments of E. directus. The beads were contained in a 0.125 m³ (33 gallon) steel drum with the robot mounted on top. The diameter of the drum was greater than 10X the width of the end effector. Three hundred and seventy one individual tests, constituting ten separate trials, were performed with this setup and are reported in this work.

Zoom In Zoom Out Reset image size