Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences
You have accessResearch articles

Reversible signal transmission in an active mechanical metamaterial

Alexander P. Browning

Alexander P. Browning

Mathematical Sciences, Queensland University of Technology, Brisbane, Australia

ARC Centre of Excellence for Mathematical and Statistical Frontiers, Queensland University of Technology, Brisbane, Australia

Google Scholar

Find this author on PubMed

,
Francis G. Woodhouse

Francis G. Woodhouse

Mathematical Institute, University of Oxford, Oxford, UK

Google Scholar

Find this author on PubMed

and

    Abstract

    Mechanical metamaterials are designed to enable unique functionalities, but are typically limited by an initial energy state and require an independent energy input to function repeatedly. Our study introduces a theoretical active mechanical metamaterial that incorporates a biological reaction mechanism to overcome this key limitation of passive metamaterials. Our material allows for reversible mechanical signal transmission, where energy is reintroduced by the biologically motivated reaction mechanism. By analysing a coarse-grained continuous analogue of the discrete model, we find that signals can be propagated through the material by a travelling wave. Analysis of the continuum model provides the region of the parameter space that allows signal transmission, and reveals similarities with the well-known FitzHugh–Nagumo system. We also find explicit formulae that approximate the effect of the time scale of the reaction mechanism on the signal transmission speed, which is essential for controlling the material.

    1. Introduction

    Mechanical metamaterials are artificially constructed and have mechanical properties defined by their structure [1]. Simple metamaterials consist of a one-, two- or three-dimensional array of elements connected by links [13] that may be elastic [47], magnetic [8,9] or electrostatic [4]. Mechanical metamaterials are highly tuneable [1012] and by altering the structure of these elements, and the properties of the links, materials have been developed that selectively transmit signals [13,14], behave as logic gates [5,15] or buckle after the application of an external stimulus [2]. There are many recent studies that experimentally realize simple mechanical metamaterials [6,14,1619]. An advantage of these designs is that they are often well suited to using three-dimensional printing technology [3,5,16,18,20]; however, a common theme among existing metamaterials is that they generally require an external source of energy to be provided in order to power their functions [5,21]. Many existing technologies can be thought of as static or inactive in the sense that they are limited by a fixed initial energy state, and are only able to respond to a finite number of stimuli before the manual introduction of external energy.

    Recent mathematical and experimental work has examined the properties of a class of one-dimensional bistable metamaterials [4,5]. These systems comprise elements, each of which consists of a mass connected to an external wall by a set of elastic elements that produce a bistable elastic potential. Individual elements are arranged in a one-dimensional lattice and are interconnected by linear springs (figure 1a). These elements may be tuned so that the elastic potential energy function is asymmetric, resulting in both a high and low potential energy stable configuration for each element (figure 1b). The system can, therefore, be designed so that an external stimulus, for example the change of a single node from the high to low potential energy stable state, can trigger a change in element configuration through the entire lattice [5]. This change is the transmission of a mechanical signal powered by stored elastic potential energy.

    Figure 1.

    Figure 1. (a) Schematic of the one-dimensional metamaterial [5]. Each element (red), of mass m, has a natural spacing of Δ. Each mass is connected to an external wall by a bistable elastic element (green) and to neighbouring elements by linear elastic springs (blue). The inset shows both stable states for each element, ui(t) =  ± δ where i is the mass index. (b) An example of bistable potential function (red) and its derivative (blue) for ai =  − 0.5 and δ = 1 in the inactive model, or the active model at t = 0. (bd) The effect of the biological reaction mechanism, which resets the potential function. Note that (d) corresponds to a reflection of (b) about ui = 0. For each plot, the stable steady states (dashed grey) and unstable steady state (dotted grey) are shown. (e,f ) Signal propagation through this material showing (e) the displacement of each element, ui(t), and (f ) the biological reaction, ai(t). The signal was initiated by moving the first element from a displacement of −δ to δ at t = 0, and was retransmitted by moving the last element from a displacement of δ to −δ at t = 400. We note the original signal could also be initiated from the right boundary by moving the last element from a displacement of −δ to δ, or indeed anywhere in the domain. Parameters used are m = 1 g, k = 1 g m s−2, γ = 1 g s−1, Δ = 0.002 m, δ = 1 m, ϵ = 0.01/s, η = 2, v = 1 g/(m2 s2) and N = 101 masses. (Online version in colour.)

    A limitation of this mechanical regime is that the system must receive an external energy input before the transmission of an additional signal [3,5]. Stiffness grading has been shown to overcome this limitation by exploiting a symmetric potential function [3]; however, these techniques may not allow the propagation of waves in systems with non-zero damping. The material in its current state is reset by manually moving each element back into its high potential energy configuration [5]. Our study introduces a theoretical biologically inspired mechanism that automatically resets each element to a high potential energy state, allowing the transmission of further signals. Many recent studies introduce the idea of manipulating biological subsystems in materials [2224], or discuss behaviours that arise in active matter systems, where biological systems exert mechanical forces [25]. Some systems are biologically inspired [21,23,24] where properties of the metamaterial are designed to mimic a biological phenomenon, and some systems exploit the properties of biological subsystems to produce new behaviours in the material [22]. One possibility for our mechanism is to exploit actin filaments in eukaryotic cells [2629] to convert energy provided to the cells as nutrients through chemical hydrolysis [30] into mechanical energy which can reset the bistable elements to a high potential energy configuration. The application of actin filaments in nanotechnology is well studied [28] and their exploitation in metamaterials has been previously suggested [29,30].

    The biological reaction mechanism we introduce is designed to react to changes in the displacement of individual elements and respond by inducing elastic potential energy back into each element. Mathematically, the effect of this process is to reset the potential function so that each element eventually reverts to a high potential energy state, as shown in figure 1bd. We present a mathematical characterization of this reaction mechanism and explore how the time scale of the reaction mechanism affects the ability of the system to transmit signals. We find that signal transmission through a coarse-grained description of the material takes the form of a travelling wave. Using a travelling wave model, we find explicit formulae that bound the parameter space for which signal transmission can occur, and approximate the effect the time scale of the reaction mechanism has on the signal transmission speed. We lay the foundation for future work on this system where the metamaterial can be tuned to produce useful new behaviours. The results we provide quantify the trade-off between the signal transmission speed and the time scale of the biological response, which are essential for control and tuning of the material. For clarity, throughout this work, we refer to the system without the reaction mechanism as the inactive system, and the system that includes the biological reaction mechanism as the active system.

    In §2, we present a mathematical model that describes the discrete active mechanical system. Following this, we take a continuous limit of the discrete model [4] with which we qualitatively explore the effect of the reaction mechanism on the ability of the system to transmit mechanical signals. We find evidence of travelling waves in the continuous model, where the wavespeed corresponds to the signal transmission speed. In §3, we solve for the wavespeed and shape in the case where the reaction mechanism is excluded. This analysis is then extended to explore the effect that the reaction mechanism has on the wavespeed and shape by taking a singular perturbation expansion (§§3a) and applying an energy conservation argument (§§3b). Finally, in §4, we discuss and summarize our results, and outline future work involving our active metamaterial.

    2. Mathematical model

    The metamaterial presented by Raney et al. [5] consists of N bistable elements of mass m, interconnected by N − 1 linear springs of stiffness k; each of the bistable elements is attached to two external walls by a pair of elastic elements, with a separation of Δ. A schematic of this physical system is shown in figure 1a. Denoting the displacement of the ith mass relative to the mean of its two steady states as ui(t), the state of the discrete system is governed by

    0=md2uidt2k(ui+1uiΔ)+γduidt+Vui,i=1,0=md2uidt2k(ui+12ui+ui1)+γduidt+Vui,i=2,,N1and0=md2uidt2k(uiui1+Δ)+γduidt+Vui,i=N,}2.1
    where γ is a damping parameter and V (ui, ai) describes the potential energy of the bistable elements that connect each mass to the external wall [4,5].

    In this study, we choose V (ui, ai) to be a quartic [4] defined in terms of its derivative,

    Vui=v(ui2δ2)(uiai),2.2
    where v describes the stiffness of the bistable elements and relates to the size of the energy gap between high and low potential energy configurations, ui =  ± δ are the stable fixed points of V (ui, ai) and ui = ai∈[ − δ, δ] is the unstable fixed point which governs the symmetry of V (ui, ai) (figure 1b). For this choice of potential energy function, ui = sign(ai)δ corresponds to an element in the high potential energy configuration and ui =  − sign(ai)δ corresponds to an element in the low potential energy configuration (figure 1b).

    In our study, we allow the symmetry of the potential function to vary in reaction to changes in the displacement, ui(t), by allowing ai = ai(t) and enforcing

    daidt=ϵ(uiηai).2.3
    Physically, ai represents a biological subsystem that receives energy from external sources and induces it into the material at a rate proportional to ϵ. The parameter η > 1 determines the extrema of the reaction parameter so that ai(t)∈[ − δ/η, δ/η]. This implementation means that if a user transmits a signal through the material by changing the displacement of a node (figure 1e), and waits a period of time of O(ϵ1) after the signal reaches the end of the domain, all elements of the system will have reverted to their high potential energy states. This resetting process of ai(t), and by extension V (ui, ai), is shown for a single element in figure 1bd, and throughout the material in figure 1f . The discrete system is similar to other fast–slow bistable systems, such as the FitzHugh–Nagumo model [3133]. In this context, we consider that the displacement function, ui(t), undergoes an excitable excursion in phase space in response to an external stimulus, and the variable representing the biological response, ai(t), behaves as a linear recovery variable. In addition to the results in figure 1e, we reproduce results in the electronic supplementary material, in the case that the second signal is initiated too early, so propagation cannot occur.

    If the length of the material is large relative to the separation of each mass, Δ, we can describe the material with a coarse-grained continuous model [34]. To derive a continuous description of the system described by equation (2.1), we consider a material of fixed length, L = (N − 1)Δ, and take the limit N so that Δ0. Following this, we define field functions u(x, t) and a(x, t) that describe ui(t) and ai(t), respectively, for x = (i − 1)Δ∈(0, L). When taking a continuous limit of the discrete system we require that the macroscopic quantities in the discrete model remain O(1) for physical reasons [4]. To do this, we replace unit quantities m, v and γ with density quantities, ρ = m/Δ, v^=v/Δ and γ^=γ/Δ, respectively, and scale the connecting spring force, k^=Δk.

    Dividing equation (2.1) by Δ and taking the limit Δ0 results in the continuous model,

    0=ρ2ut2k^2ux2+γ^ut+v^(u2δ)(ua)=02.4
    and
    0=atϵ(uηa).2.5
    A no flux boundary condition ∂u/∂x = 0 is applied at x = 0 and x = L.

    A key aspect of this study is to investigate the signal transmission speed through the parameter space, particularly as ϵ increases. We therefore non-dimensionalize equations (2.4) and (2.5) by scaling t=Tt^, x=Xx^, u=Uu^ and a=Aa^, where hat notation represents dimensionless variables. Choosing U = A = δ, T2=ρ/γ^ and X2=k^T2/ρ gives

    0=2u^t^22u^x^2+u^t^+ν(u^21)(u^ηa^)2.6
    and
    0=a^t^κ(u^ηa^),2.7
    where κ = ≥0, ν=v^δ2ρ/γ^2>0, η > 1 and x^(0,L^). The behaviour of the system can now be studied through the three-dimensional parameter space (ν, η, κ), where ν > 0 is the relative strength of the potential function; η > 1 describes the steady-state locations of a^(x^,t^); and κ≥0 is the relative time scale of the reaction mechanism. In this non-dimensional regime, the stable states are located at u^=±1, where the high potential energy state is always given by sign(a^).

    Previous studies provide evidence to suggest that the transmitted energy is input independent [35], so we do not expect the initial condition to affect the transmission speed in the centre of the domain. A signal is initiated in the continuous model using a Heaviside initial condition for the displacement, where

    u^(x^,0)={1,0x^<Q,1,Qx^L^,2.8
    and the biological response is kept as it was before the signal was initiated by setting
    a^(x^,0)=1η.2.9

    In figure 2a, numerical solutions to the continuous model show that the transition of the displacement variable, u^(x^,t^), is carried by a wave which appears to approach a constant shape and speed. Figure 2b demonstrates the slow biological response, where a^(x^,t^) also undergoes a slow transition in response to changes in u^(x^,t^). Results in figure 3 illustrate the dependence of the transmission speed on the time scale of the response, κ. These results indicate a negative monotonic relationship between κ and the transmission speed. Full details of the numerical scheme used to solve the continuous model are provided in the electronic supplementary material.

    Figure 2.

    Figure 2. Numerical solutions to the non-dimensional continuous model (equations (2.6) and (2.7)) showing (a) the displacement field, u^(x^,t^), and (b) the biological reaction, a^(x^,t^). Parameters used are ν = 4, η = 2 and κ = 0.01. The signal is initiated using the initial condition described by equations (2.8) and (2.9) with Q = 1. (Online version in colour.)

    Figure 3.

    Figure 3. The solution to the continuous model, equation (2.6), shown as kymographs where colouring represents (ad) the displacement function, u^(x^,t^), and (eh) the biological response, a^(x^,t^), for increasing values of κ. Parameter values used are ν = 4 and η = 2. Each signal is initiated using the initial condition described by equations (2.8) and (2.9) with Q = 1. (Online version in colour.)

    3. Travelling wave model

    Figure 2 suggests that signals are propagated through the system by waves which appear to approach a constant shape and speed. Motivated by this, we now look for a travelling wave solution to the continuous model by extending the domain to represent a material of infinite length, so that x^R+ [3,4,32,33,3642]. We define the wavespeed, c, and without loss of generality enforce c > 0 by investigating travelling waves that are initiated at the left boundary and move in only the positive x^ direction. In reality, we expect symmetry in these solutions as a travelling wave may also travel in the negative x^ direction if the signal is initiated in the centre of the domain.

    Substituting travelling wave variables f(z)=u^(x^,t^) and h(z)=a^(x^,t^), where z=x^ct^R, into equations (2.6) and (2.7), and dividing by c2 − 1, gives the travelling wave model,

    0=d2fdz2+c1c2dfdzν1c2(f21)(fh),f(±)=13.1
    and
    0=dhdz+κc(fηh),h()=1η,3.2
    which, in physical coordinates, corresponds to a wave profile connecting u^=1 to u^=1 as t.

    In figure 3, numerical solutions of the continuous model show that increasing the response speed κ typically reduces the wavespeed from a maximum which occurs when κ = 0. We, therefore, denote a fast wave as a travelling wave solution at, or near, κ = 0, and a slow wave as a travelling wave solution in the limit, or near, c = 0. Removing the active component of the system from the model by setting κ = 0 (and defining c0 as the corresponding wavespeed) gives

    0=d2f0dz2+c01c02df0dzν1c02(f021)(f0+1η),f0(±)=13.3
    and
    0=dh0dz,h0()=1η,3.4
    where f0(z) will correspond to the solution to the fast wave. Under certain parameter transformations, equation (3.3) is analogous to the well-studied bistable equation that arises from analysis of the FitzHugh–Nagumo model [3133,38]. The solution of equations (3.3) and (3.4) is therefore
    f0(z)=tanh(μ0z)3.5
    and
    h0(z)=1η,3.6
    where
    μ0=ν2η2+ν2andc0=1(η2/2ν+1).3.7
    For physically realistic parameter combinations, equation (3.7) suggests that c0 < 1.

    (a) Perturbation solution for the fast wave

    An exact solution in the limiting case κ = 0 (equation (3.5)) allows the formulation of a perturbation solution for 0 < κ≪1 [38,39,43]. Figure 4 shows an approximation to the wave profiles f(z) and h(z) for κ = 0.01. A fast transition region is seen in f(z) around z = 0 (figure 4a), suggesting behaviour necessitating a singular perturbation analysis [39,43]. For zO(μ01), the solution appears to match a solution for κ = 0, since, at this scale, h(z) is approximately constant (figure 4b). In figure 4a, we show that h(z) is not constant but rather a slow reaction of O(κ1), which transitions h(z) =  − 1/η to h(z) = 1/η as z. This observation further suggests a singular perturbation analysis, since the behaviour of the response mechanism h(z) varies significantly for 0 < κ≪1 compared with κ = 0.

    Figure 4.

    Figure 4. Approximation of the travelling wave profile using a solution to the continuous model for (a) the displacement function, f(z), and (b) the biological reaction, h(z). Parameters used are ν = 4, η = 2 and κ = 0.01. In (a), the horizontal axis is scaled to show both outer regions where zO(κ1), which is the rate of change of h(z). In (b), the horizontal axis is scaled to show the inner region where zO(μ01), which is the rate of change of f(z). (Online version in colour.)

    We propose a three-part perturbation solution about κ = 0 and define independent variables zO(μ01) to correspond to an inner region, and Z = κz for zO(κ1) to correspond to two outer regions. These three regions are shown in figure 4a. This regime requires the fast process, which occurs in the inner region, to be much faster than the slow process, which occurs in the outer region. That is, we require κμ0.

    Our aim in looking for a perturbation solution is to determine the effect of small perturbations in κ on c. To do this, we pose a perturbation expansion for the wavespeed through the entire domain,

    c=c0+c1κ+O(κ2),3.8
    where c0 is given by equation (3.7).

    In the inner region, we pose a perturbation solution of the form

    f(z)=f0(z)+f1(z)κ+O(κ2)3.9
    and
    h(z)=h0(z)+h1(z)κ+O(κ2),3.10
    where f0(z) and h0(z) correspond to the shape of the fast wave at κ = 0, given by equation (3.5). Full details of the perturbation solution in the inner region are given in the electronic supplementary material. In summary, the solution is given by
    f(z)=tanh(μ0z)+f1(z)κ+O(κ2)3.11
    and
    h(z)=1η+1c0μ0η(log(cosh(μ0z))μ0z+log(2))κ+O(κ2),3.12
    where f1(z) and c1 are defined by the solution of a second-order boundary value problem, which can be solved numerically. Full details of this numerical scheme are given in the electronic supplementary material. These solutions are shown in figure 5, and the wavespeed correction, c1, is summarized for various parameter combinations in electronic supplementary material, table S3.
    Figure 5.

    Figure 5. Comparison of the wave shape from a numerical solution to the continuous model (coloured solid curves) with the perturbation solution (black dotted and dotted-dashed curves) for (a) the displacement function, f(z), to O(1) and (b) the reaction mechanism, h(z), to O(κ). The inset (c) shows the crossover between the inner and outer perturbation solutions to h(z) in further detail, showing an excellent match to the wave shape from the continuous model, and between the inner and outer solutions. (Online version in colour.)

    In the outer region, we denote solutions to the slow system using uppercase variables, F(Z) and H(Z). We expect the outer solution to apply for ZO(1)zO(κ1) and as z,Z. In the outer region, equation (3.1) becomes

    0=κ2d2FdZ2+cκ1c2dFdZν1c2(F21)(FH),F(±)=1,3.13
    which has the solution
    F(Z)=sign(Z).3.14
    This agrees with the sharp transition region seen at this scale in figure 4a. Substituting equation (3.14) into equation (3.2) we see that
    0=dHdZ1c(sign(Z)η+H),H()=1η.3.15
    Full details of the perturbation solution in the outer region are given in the electronic supplementary material. In summary, the solution to equation (3.15) is given by
    H(Z)={1η2ηexp(Zc0)2c1ηc02Zexp(Zc0)κ+O(κ2),Z<0,1η,Z>0.3.16

    Figure 5bc shows a comparison between solutions for the reaction mechanism in the inner and outer regions (given by equations (3.11) and (3.16), respectively) and an approximation of the wave shape formed from the numerical solution to the continuous model. We see an excellent match between both the continuous model and the perturbation solutions, as well as between the inner and outer regions of the perturbation solution. These results provide excellent information about the behaviour of h(z) for 0 < κ≪1, and are particularly important as the behaviour for 0 < κ≪1 varies significantly from the behaviour at κ = 0.

    (b) Energy transport

    To determine information about the slow wave, and to obtain more information about the fast wave, we follow Nadkarni et al. [4,41] to derive an integrability condition—that is, a necessary condition for the existence of a solution with given boundary conditions—to investigate the transported energy. Multiplying the travelling wave model, given by equation (3.1), by df/dz and integrating gives

    0=d2fdz2dfdzdz+c1c2(dfdz)2dzν1c2(f21)(fh)dfdzdz.3.17

    For a parameter regime where the transition wave exists, the velocity will vanish in the far field so that f1 and df/dz0 as z±. Therefore, some components of equation (3.17) vanish,

    d2fdz2dfdzdz=[12(dfdz)2]=0and(f21)(fh)dfdzdz=dfdz(f21)hdz.
    Under the assumption |c|≠1, equation (3.17) becomes an integrability condition,
    c(dfdz)2dz=νdfdz(f21)hdz.3.18
    It is useful to note that, for f =  − tanh(σz), for some σ, which occurs as κ0 with σ = μ0,
    (f21)dfdzhsech4(σz)h.3.19
    This suggests that only the component of h(z) near z = 0 is important in gaining any approximation from the integrability condition, provided f(z) has the form of a hyperbolic tangent function. Since a travelling wave connecting f(z) = 1 to f(z) = − 1 as z will always have a sigmoidal form (figure 6), we expect this observation to apply for all regions of the parameter space where a travelling wave exists.
    Figure 6.

    Figure 6. Approximate numerical solutions to the travelling wave model, formulated from numerical solutions to the continuous model, for (a) the displacement function, f(z), and (b) the biological reaction, h(z). Results are shown for various values of the response parameter, κ, with ν = 4 and η = 2. Additionally, (a) shows the exact solutions to the travelling wave model, equations (3.1) and (3.2), for κ = 0 (black dashed-dot) and κ = κ* = 1.6 (black dashed). (Online version in colour.)

    In the inactive model, where κ = 0, Nadkarni et al. [4] show that the integrability condition reduces to

    Ek=cΔV2γ,3.20
    where Ek represents the total kinetic energy per density transported by the transition wave, and ΔV represents an energy gap or difference in the potential energy between the high and the low potential energy states. This result can be used to find an upper bound for c for κ > 0: since |h(t)|≤1/η and ΔV is a monotonically decreasing function of h(t), the available kinetic energy in the system is always bounded above by that which occurs when h(t)≡1/η, which is the case for κ = 0. This suggests that c(κ)≤c0 < 1 and substantiates numerical evidence seen in figure 3, which suggests a decreasing monotonic relationship between c and κ.

    In the following subsections, we apply the integrability condition to obtain approximations to c(κ) while holding ν and η constant. We also find the region of the parameter space that allows signal transmission. The advantage of these approximations is that they avoid numerical solutions to the continuous model to approximate the wavespeed. Numerical solutions to this model are not computationally inexpensive and it is generally difficult to obtain and verify the results.

    (i) Energy transport in the fast wave

    For κ≪1, it is reasonable to assume that f ≈  − tanh(μz), where μ depends on κ and μμ0 as κ0. By assuming f(z) has a similar form to f0(z) for κ≪1, it is reasonable to use the perturbation solution (equation (3.12)) as an approximation for h(z). Since the integrability condition (equation (3.18)) depends only on the component of h(z) near z = 0, we use only equation (3.12). Allowing μ0 = μ(0) and c0 = c(0), where μ = μ(κ) and c = c(κ) depend on κ, we assume that

    h(z)1η+κcμη(log(cosh(μz))μz+log(2)),κ1.3.21

    Equation (3.21) corresponds to the smoothed piecewise-defined function where growth is equal to 2κ/() for z < 0, and zero for z > 0. This function corresponds exactly to an O(1) approximation to h(z), which uses f(z) =  − sign(z). Substituting equation (3.21) into the integrability condition (equation (3.18)) provides the relationship

    4c2μ23=2ν(6cμ5κ)9η.3.22
    This result does not allow c and μ to be determined independently, so we consider a far-field expansion of the travelling wave [36], that is, an expansion about ϕ ≈ 0, where ϕ = eμz,
    f(z)=tanh(μz)=e2μz1e2μz+1=ϕ21ϕ2+1=1+2ϕ2+O(ϕ4).3.23
    Substituting ϕ = eμz and equation (3.23) into the travelling wave equation for h(z) (equation (3.2)) gives
    f(z)1+2e2μz,|z1|1andh(z)1η+2κη(κ+2cμ)e2μz,|z1|1.}3.24
    Substituting equation (3.24) into equation (3.1) provides a far-field relationship,
    0=νηcμν+2μ2(1c2).3.25

    To obtain a useful analytical expression for c and μ, we pose a perturbation solution to equations (3.22) and (3.25) around κ = 0, such that

    c=c0+c~1κ+O(κ2)3.26
    and
    μ=μ0+μ~1κ+O(κ2),3.27
    where c0 and μ0 are given by equation (3.7). We note that c~1 approximates c1, the gradient of c(κ) at κ = 0, where c1 is defined exactly by the solution of the perturbation problem. Substitution into equations (3.22) and (3.25) gives
    c~1=5ν(4c02μ0+c04μ0)24μ02(2c0ημ0ν).3.28
    We compare this estimate of c1, which is calculated numerically by solving the boundary value problem that comes from the singular perturbation expansion, with that estimated from the continuous model in electronic supplementary material, table S3. In figure 7, we show that, for κ≪1,
    c(κ)c0+c~1κ3.29
    matches numerical results for c(κ).
    Figure 7.

    Figure 7. Comparison between the wavespeed as a function of κ, from the continuous model (circles), with: the analytical expression at κ = 0, given by equation (3.7) (diamonds); the analytical expression for κ* at c = 0, given by equation (3.37) (squares); the analytical approximation given by equation (3.29), which applies for κ≪1 (dashed-dotted curve); the analytical approximation given by equation (3.41), which applies for c≪ 1 (dashed curve); the combined approximation given by equation (3.42), which applies everywhere (dotted curve); and the approximation given by equation (3.44), which applies everywhere (solid curve). (Online version in colour.)

    (ii) Energy transport in the slow wave

    We denote f*(z) and h*(z) as the shape of the slow wave, which occurs as c0 and κκ, where κ* is yet to be determined. In addition, we expect the curve c(κ) to be perpendicular to the κ axis at κ* to maintain continuity in the symmetry of the problem where solutions with a negative wavespeed are equally valid. To determine a governing equation for the slow wave, we take c0, so that equations (3.1) and (3.2) become

    0=d2fdz2ν(11η)f(f21),f(±)=13.30
    and
    0=fηh,3.31
    which have the analytical solution
    f(z)=tanh(μz)andh(z)=1ηtanh(μz),}3.32
    where
    μ=ν2(11η).3.33

    Direct substitution of f*(z) and h*(z) into the integrability condition (equation (3.18)) causes all terms to vanish, so higher order behaviour of O(c) as c0 is important. To allow for this, we note that equation (3.2) can be solved, for all κ, to give

    h(z)=κcηzf(s)exp[κc(zs)]ds.3.34
    Equation (3.34) shows that h(z) depends on the product of f(s) and a function that decays rapidly as z moves away from s. Expanding f(s) in a Taylor series about z gives
    h(z)κcηz[n=0(sz)nn!dnfdzn]exp(κc(zs))ds=1ηn=0cnκndnfdzn.3.35

    Assuming that f∼tanh(μ(κ)z) provided c≪1, where μ(κ*) = μ*, and truncating the infinite series given by equation (3.35) after n = 3, the integrability condition (equation (3.18)) gives the relationship

    4κ33=16ν(7κ28c2μ2)105η.3.36
    For η and ν fixed, c0 only as κκ. Substituting c = 0 into equation (3.36) gives
    κ=45νη.3.37
    Under the assumption that c(κ) is monotonically decreasing, the result in equation (3.37) provides an analytical expression for the region of the parameter space where we expect signal transmission. In figure 7, we show that this expression matches the numerical results. Interestingly, equation (3.37) depends only on the ratio of the other parameters, ν/η. The scales of the horizontal axes for figure pairs figure 7a,c and b,d have been chosen to highlight this.

    To gain information about the shape of c(κ) near κ = 0, we pose a perturbation solution around κ = κ* to the system governed by equation (3.36) and the far-field relationship, equation (3.25), where

    c=c+c12κκ+O(κκ)3.38
    and
    μ=μ+μ12κκ+O(κκ),3.39
    where μ*, given by equation (3.33), and c* = 0 apply at κ = κ*. Substitution of these expansions into equations (3.25) and (3.36) gives
    c12=75(η1),3.40
    so that, for κ ≈ κ*,
    c(κ)7(κκ)5(η1).3.41
    Figure 7 shows that equation (3.41) matches numerical results for c(κ) for a surprisingly wide range of κ < κ*, especially for larger η. This result is particularly important as we have not constructed a full perturbation solution to the travelling wave model about c = 0.

    (iii) Combined approximation

    We can combine the approximations to c(κ) from the slow and the fast wave to obtain a curve that behaves like equation (3.29) for κ0 and like equation (3.41) as κκ. To do this, we propose a form like

    c(κ)α1+α2κ+α3κκ,3.42
    where we choose
    α1=(c0+2c~1κ),α2=c0κ+2c~1andα3=2(c0+c~1κ)κ,3.43
    where c0 is given by equation (3.7); c~1 is given by equation (3.28); and κ* is given by equation (3.37). In figure 7, we show that this combined approximation provides a reasonable approximation to c(κ) for κ∈[0, κ*].

    (iv) Whole domain ansatz

    Results in equations (3.5) and (3.32) show that f(z) is described exactly by a hyperbolic tangent at both κ = 0 and κ = κ*, and figure 6a suggests that f(z) remains sigmoidal. Therefore, it may be reasonable to approximate f(z) ≈ tanh(μ(κ)z) for κ∈[0, κ*], where μ(κ) depends on κ. Additionally, assuming that μ(κ) decays linearly from μ = μ0 at κ = 0 to μ = μ* at κ = κ*, we obtain an approximation to f(z) for all κ,

    f(z)tanh(μ(κ)z)andμ(κ)=μμ0κκ+μ0.3.44
    Substituting the approximation for μ(κ) given by equation (3.44) into the far-field matching condition (equation (3.25)) gives an approximation to c(κ) which can apply for κ∈[0, κ*]. We show that this approximation is reasonably accurate throughout κ∈[0, κ*] in figure 7; however, it is clear from numerical results in figure 6, when κ = 0.4, that the solution does not have the symmetry of a hyperbolic tangent function.

    4. Discussion and conclusion

    Currently, the inactive metamaterial described mathematically by Nadkarni et al. [4] and experimentally realized by Raney et al. [5] is able to transmit mechanical signals by the release of stored potential energy. A limitation of this design is that mechanical energy must be manually introduced into the system before additional signals can be transmitted. Our study presents a novel biologically inspired metamaterial that incorporates a theoretical biological mechanism that harvests energy to reset the system to a high potential energy state, allowing the transmission of additional signals. Energy may be induced into the active metamaterial through a biological process, such as actin filaments in eukaryotic cells [2628]. That said, our analysis does not necessarily require this mechanism to have a biological origin: the reaction mechanism may also represent a mechanical system where energy is added through other electrochemical [44], photovoltaic, thermodynamic [20] or pneumatic [21] subsystems.

    By finding evidence of travelling wave solutions, we are able to analyse limiting behaviour describing the signal transmission speed and wave shape. We provide a detailed analysis to qualitatively and quantitatively understand the effect of our reaction mechanism on signal transmission abilities of the material. Our main results consist of a set of analytical approximations that quantify the signal transmission speed as a function of the parameters which describe the physical properties of the material. Results in figure 7 show that the approximation we develop to apply through the whole domain, given by equation (3.42), provides an excellent match to the numerical results, particularly for large η. In addition, our approximation for the wavespeed near the slow wave, given by equation (3.41), is able to provide excellent information about the shape of c(κ) as c0, which we find is difficult to obtain numerically. This approximation is also able to provide a region of the parameter space for which signal transmission can occur, given by κ < 4ν/(5η) (equation (3.37)). This understanding of the effect of our mechanism on the signal transmission speed is useful as it allows our active metamaterial to be tuned to produce desirable new behaviours. For example, our results allow quantification of the trade-off between signal transmission speed and the response time, which is essential for controlling the material. These insights are also essential for building a material containing a biological mechanism that induces energy into the system. Decreasing ν and η in the same proportion increases the transmission speed at the cost of increased sensitivity to noise-induced misfiring, but may be essential if the energy budget is small.

    A key aspect of our study is to follow Nadkarni et al. [4] by representing the bistable potential energy function as a quartic (equation (2.2)). This approach leads us to obtain numerous analytical approximations that characterize the effect of the biological mechanism on the transmission speed which, although qualitatively reliable, may not always be quantitatively appropriate for particular systems [5]. In fact, the analytical expression for the transmission speed is a result of the similarity between our model and the well-studied bistable equation [38]. These choices mean that our system has mechanical and algebraic properties that are similar to other bistable systems, such as the FitzHugh–Nagumo model [31,33]. That said, we do not assume that the time scale of the response is significantly slower than the time scale of the excitement, as is often the case in analysis of such models. Indeed, our aim is to develop an intelligent biomechanical material that has tuneable properties. In some sense, it is desirable that the response is as fast as possible to allow for a short period of time between signal reception and retransmission. Future work may examine the role of heterogeneities in the properties of the material [3,45]. Such features could allow the material to selectively transmit signals by creating energy barriers that interact with signals of certain properties [13].

    The travelling wave analysis we conduct assumes a material of infinite length over a large period of time. However, applications of our material will have a finite length and may have properties not suitable for a continuum model. For example, in a material where the spacing between elements is not significantly different from the length of the material, a discrete travelling wave analysis may be more appropriate. The discrete problem is known to be substantially more difficult than the continuous problem [38], so the limiting transmission speed our analysis provides may still be useful. Furthermore, the inclusion of our biologically inspired mechanism can be incorporated into passive metamaterials of higher dimensions to enable new behaviours and the travelling wave analysis can be extended to investigate two-dimensional signal propagation. In the electronic supplementary material, we produce results which show the transmission of concurrent signals (electronic supplementary material, figures S1 and S2) and interacting signals initiated from both ends of the material (electronic supplementary material, figure S3). Further analysis is needed to examine the behaviour of these types of interacting waves [46] and the material's ability to transmit oscillatory or concurrent signals.

    To conclude, we have presented a novel, biologically inspired, active metamaterial that can reversibly transmit mechanical signals. This work provides an analytical expression that describes the mechanical properties of the material required for signal transmission. We also provide numerous approximations that quantify the effect of the mechanical properties, and the time scale of the biological response, on the transmission speed. This work demonstrates how a new class of biologically inspired metamaterials are able to produce useful new functionalities. The type of analysis we present is invaluable for tuning and controlling the active metamaterial.

    Data accessibility

    Additional data are provided in the electronic supplementary material. Key algorithms used to generate results are available on Github at github.com/ap-browning/rspa-2019.

    Authors' contributions

    All authors conceived and designed the study; A.P.B. performed the analysis and numerical simulations, and drafted the article; all authors provided comments and gave final approval for publication.

    Competing interests

    We declare we have no competing interests.

    Funding

    This work is supported by the Australian Research Council (DP170100474) and a Royal Society exchanges grant (no. IE160805).

    Acknowledgements

    We thank Kevin Burrage and Ian Turner for their helpful discussions. We also thank the three anonymous referees for their comments.

    Footnotes

    Electronic supplementary material is available online at https://dx.doi.org/10.6084/m9.figshare.c.4563902.

    Published by the Royal Society. All rights reserved.

    References