Blindfolded or disoriented people have the tendency to walk in circles rather than on a straight line even if they wanted to. Here, we use a minimalistic walking model to examine this phenomenon. The bipedal spring-loaded inverted pendulum exhibits asymptotically stable gaits with centre of mass (CoM) dynamics and ground reaction forces similar to human walking in the sagittal plane. We extend this model into three dimensions, and show that stable walking patterns persist if the leg is aligned with respect to the body (here: CoM velocity) instead of a world reference frame. Further, we demonstrate that asymmetric leg configurations, which are common in humans, will typically lead to walking in circles. The diameter of these circles depends strongly on parameter configuration, but is in line with empirical data from human walkers. Simulation results suggest that walking radius and especially direction of rotation are highly dependent on leg configuration and walking velocity, which explains inconsistent veering behaviour in repeated trials in human data. Finally, we discuss the relation between findings in the model and implications for human walking.
In human walking, the body is alternately supported by the left and right leg, with a double support phase during the load transfer from one leg to the other. The bipedal spring-loaded inverted pendulum (BSLIP) model  is a mimimalisic gait template  which mimics key features of human walking: besides the alternating body support, each leg shows an M-shaped ground reaction force similar to what is observed in humans. Owing to its simplicity, it fails to reproduce other characteristics of human walking, such as typical contact times or force amplitudes, especially for very slow and fast walking . However, we decided to use this model for our present modelling study, as it provides a good compromise between similarity to human data and simplicity. An important idea that comes with simplicity is the hope that results derived from a simple template still persist in more elaborate models, and thus might provide insights into fundamental biomechanical relations.
The BSLIP model originally focuses only on the sagittal plane, whereas human walking happens in three dimensions. This is especially apparent when we change walking direction. But we do not only turn voluntarily. When we are blindfolded or have no orientation cues, we walk rather in circles than in straight lines even if we wanted to [4–6]. Interestingly, whether or not the turning direction is a function of differences in leg lengths has been debated for decades [4,5,7–9].
In this study, we extend BSLIP to three dimensions. We are interested if stability persists, and how differences in the legs lead to walking in circles. We transfer the three reference patterns (A,B,C) from the original manuscript  into a three-dimensional gait by adding an azimuthal angle to the leg. These three points correspond to three different kinds of periodic solutions of BSLIP. Solutions of type A show a symmetric M-shaped force profile, whereas solutions of type B have an asymmetric M, with the first force hump having larger amplitude. Type C solutions show a comparatively shallow force profile, typically accompanied with low force magnitude, long double support phases and high cadence.
2.1. Model description
The model (figure 1) consists of a point mass on top of two massless springs. The motion alternates between single and double support phases. Double support starts when during single support the swing leg touches ground. Single support occurs when during double support the trailing leg's force becomes zero. Then, the trailing leg becomes the swing leg. Note that in contrast to similar BSLIP definitions [1,10,11], flight phases cannot occur in this model definition. We say that a stride begins with the onset of single support of leg 1, that is, the take-off event of leg 2.
Each leg exerts force only along its direction from tip to the centre of mass (CoM). The magnitude is , with and li(t) denoting the ith leg's stiffness, rest length and momentary length, respectively. Here, the rest length l0 is 1 m. We restrict the leg force to positive values (pushing), that is, if the magnitude of the leading leg's force becomes less than 0 during stance, it is set to 0. A negative force of the trailing leg triggers its transition to swing. This prevents flight phases with two legs in swing phase, and enforces an alternating touchdown sequence of the legs. During swing, the leg does not exert force. The equations of motion for the CoM are 2.1with g denoting gravity (9.81 ms−1) and m denoting body mass. The values for k, l0 are selected according to Geyer et al. .
In this study, we apply two different leg placement policies. Primarily, and if not stated otherwise, we use a body-oriented (BO) leg orientation scheme: during swing, the leg is aligned with respect to the current sagittal plane which is spanned by the current CoM velocity vector and gravity. The leg's orientation with respect to the current sagittal plane is given by the polar angle α and the azimuthal angle β (figure 1). The values for α and β are kept constant for every stride, but can differ between legs. By construction, this model is invariant to rotations in the horizontal plane, because it does not refer to any external direction except gravity. We refer to this model as the BO-model.
For comparison, we use a world-oriented (WO) leg placement policy. Here, the only difference is that the sagittal plane which is used for leg alignment is spanned by time-invariant coordinate axes of a world frame, and thus is itself time-invariant. This model is not invariant to rotations of the horizontal plane, as a rotation of the horizontal plane would affect the sagittal plane and thus the swing leg orientation. We refer to this model as the WO-model.
We compare the stability of periodic straight walking solutions which share the same motion and differ only in the leg placement policy.
2.2. Relative periodic coordinates
Right steps typically accelerate the CoM to the left, whereas left steps do the opposite. To return to an initial lateral position and velocity, we thus need to perform two steps. Here, a stride consists of a right step followed by a left step. More precisely, in our model, we are seeking two-step relative periodic solutions, that is, solutions that have the same state after two steps up to a symmetry operation (here: translation and rotation in the horizontal plane). In our case, a solution is relative periodic if the reduced state xc, defined as , is periodic. l1 denotes the stance leg's current length, |v| denotes the magnitude of the velocity and ζ denotes the azimuthal angle between instantaneous sagittal plane and the stance leg (figure 1).
Given xc, the state in Euclidean coordinates is determinable up to the model's orientation and position in the horizontal plane. This implies that rotations of the model in the horizontal plane do not alter the relative periodic coordinates. The one missing dimension of the five-dimensional relative periodic coordinates, compared with the six-dimensional coordinates of the full state, corresponds exactly to the invariance of the model with respect to rotations in the horizontal plane. As a consequence of this state definition, strides whose final state (measured with respect to the foot contact point) are identical to the initial state up to a rotation in the horizontal plane have the same periodic coordinates at the beginning and the end. This allows us to search for relative periodic solutions as fixpoints of return maps in relative periodic coordinates.
2.3. Finding relative periodic solutions
We select the model parameters to be m = 80 kg, E = 816 J and following the original paper of Geyer et al. . We then transfer the three highlighted periodic gaits from Geyer et al. to the three-dimensional walking model: configuration A: k1,2 = 14 kN m−1, α1,2 = 69°, B: k1,2 = 14 kN m−1, α1,2 = 72.5° and C: k1,2 = 20 kN m−1, α1,2 = 76.5°. Note that we slightly altered the polar angles α for B and C to obtain more robust solutions.
To find a relative periodic solution for given parameters in the BO-model, we seek an initial state xc which results in the identical state after two steps of the model. We formulate a difference function D by D(xc) = xc − S(xc), with S(xc) denoting the model's state after two steps starting from xc. Zeros of this difference function D now correspond to relative periodic solutions. We use the Newton–Raphson algorithm to find zeros of D. Because we look for periodic solutions of a given energy, we allow variations in initial conditions that leave only the system energy unchanged. We do so by removing |v| from xc for the Newton–Raphson algorithm, and computing |v| prior to every simulation run from the remaining states of xc, the model parameters and the required energy.
For periodic solutions of each symmetric configurations in the BO-model, we compute a corresponding solution in the WO-model. This is achieved by rotating the BO-model's initial conditions in the horizontal plane such that after one stride, the lateral (z) coordinate of the BO-model's CoM remains unchanged. We use the same initial conditions, but we have to adapt the angle β, because, at touchdown, the sagittal plane now points into another direction. For the WO-model, the corresponding azimuthal angle β for the swing leg orientation is computed as the sum of the BO-model's azimuthal angle and the angle between CoM velocity and xy-plane at touchdown of the BO-model's stride. This policy yields the same ground contact points (foot locations) in both models. We successfully verified the periodicity of these solutions in the WO-model and the equivalence of the trajectories to the corresponding trajectories of the corresponding motion of the BO-model.
2.4. Stability analysis
In this study, we refer to stability as asymptotic stability. We apply a Poincaré map  analysis to compute stability. The Poincaré map was computed using small disturbances of the initial conditions, to obtain the derivative of the final state with respect to the initial conditions (the Jacobian). If all eigenvalues of the Jacobian have magnitude less than one, then the model is locally asymptotically stable.
The BO-model's and WO-model's full states each have six dimensions, thus their Jacobians of the full state have six eigenvalues. The location of the model's CoM is measured with respect to the current (single) stance leg's contact point. We expect the following properties of the Jacobian: (i) there is one direction that corresponds to a perturbation along the trajectory. Such perturbations do not affect the final state at all, which corresponds to an eigenvalue of zero. (ii) As the model cannot change its energy, we expect that energy changes owing to perturbations stay exactly the same after one stride. This corresponds to an eigenvalue of 1. (iii) The BO-model is invariant with respect to rotations in the horizontal plane; after a stride, the rotation will be exactly preserved. This ‘neutrality to walking direction’ corresponds to another eigenvalue of 1.
We use the six-dimensional return map only to compare the BO-model with the WO-model. Further, we deviate from textbook stringency and say a solution is asymptotically stable if all but two eigenvalues have magnitude less than one, and two eigenvalues equal 1.
In the BO-model, we can alternatively use the five-dimensional relative periodic coordinates xc to compute the return map, which take neutrality to rotations in the horizontal plane into account. Further, we take into account that the model preserves the total energy, by allowing only initial conditions for return map computation that correspond to the same energy. This is accomplished by computing the coordinate |v| of the relative periodic coordinates from the other coordinates, system parameters and desired system energy. We now vary only the four remaining dimensions of the relative periodic coordinates xc to compute the return map, which yields a four-dimensional map that does not account for variations of system energy or rotations in the horizontal plane. In other words, the two neutral directions that we expect in the full state return map—one related to changes of energy, the other related to orientation in the horizontal plane—are removed from this map. Note that the BO-model itself is still neutral to energy perturbations and rotations in the horizontal plane. Using this reduced return map, we call a periodic solution asymptotically stable if and only if all eigenvalues of this four-dimensional map have magnitude strictly less than one.
Eigenvalues tell us how quickly small perturbations decay. They contain no information about the magnitude of perturbations that can be mitigated, that is, the size of the basin of attraction. In the two-dimensional model, Rummel et al.  showed that the stability, indicated by the magnitude of the largest eigenvalue of the return map, is not correlated to the robustness as indicated by the size of the basin of attraction. Computing the robustness is beyond the scope of this study.
2.5. Influence of leg parameters
The BO-model's walking motion, projected to the horizontal plane, does not exactly follow a circle (figures 2 and 4). Still, we can define and compute a radius of the circular walking motion. For this, we take the angle γ by which the CoM trajectory has turned after one stride, and the distance covered d (figure 2). The walking radius r then is obtained by elementary geometry: Note that this radius depends on the choice of the Poincaré section which separates strides. If the Poincaré section corresponds to more outer or more inner points, the radius will increase or decrease accordingly.
To study the BO-model's parameter dependence, we introduce the difference parameters Δk, Δα, Δβ and Δl0, which denote the differences between both legs in the corresponding leg parameters. The actual model parameters are then computed by k1,2 = k0 ± 1/2Δk, etc. Here, k1,2 denotes the stiffness for leg 1 or 2, respectively, and k0 denotes the reference parameter from solution A, B or C. The new parameter β0 is antisymmetric, that is, β2,0 =−β1,0. We further varied the total system energy E for selected solutions, to investigate the effect of energy changes on the periodic motions.
We compute the influence of leg parameters on the walking radius twofold: first, we compute the corresponding relative periodic solutions for A, B and C with symmetric leg parameters (Δk = Δα = Δβ = Δl0 = 0) as described above. For symmetry reasons, the corresponding motions show a straight walking motion. Starting from these symmetric solutions, we compute the derivative of the walking curvature (which is the reciprocal of the walking radius) with respect to difference parameters. We choose the walking curvature, because the walking radius is infinite at the reference solution. For every variation in the difference parameters, we compute the corresponding relative periodic solution of the BO-model and calculate its walking curvature.
Additionally, we provide a map of the walking radius and the largest eigenvalue for a subset of different leg parameters in the BO-model, namely with varying Δk and Δα. Starting with the straight walking patterns A, B and C, we set β1,2 = ±0.05 rad or β1,2 = ±0.1 rad, vary Δα in steps of 0.1°, and vary Δk in steps of 50 Nm−1. We then search a stable periodic solution, starting from the nearest periodic solution found until then. When no stable relative periodic solution is found, we increase Δα and repeat the variation of Δk. This mapping yields a qualitative and quantitative dependency of the walking radius from variations in leg parameters. However, this mapping technique does not probably discover every stable periodic solution that exists.
We found stable relative periodic solutions in the BO-model when configurations A, B and C were modified with a leg orientation pointing into the lateral direction, as described above (β1,2 ∈±0.05, ±0.10 rad). This yielded walking patterns with step widths between 2.3 and 8.5 cm, which is less than the preferred step widths of human walking (≈13 cm ). For every configuration, we modified Δk, Δα, ΔL0 and Δβ until no solution could be found or instability was detected. For each modulated parameter, we found stable solutions for a certain range of parameter variation, indicating structural stability of the model.
To compare the leg placement policies of the BO-model and the WO-model, we transferred those stable straight walking solutions of the BO-model to the WO-model (see Methods section). We verified that the trajectories remained unchanged when we transferred them from the BO-model to the WO-model. Subsequently, we computed the six-dimensional full state return map for these periodic solutions in the WO-model, and found that none of these solutions remained stable (at least one eigenvalue was strictly larger than 1) when leg placement policy was exchanged. In fact, we did not observe any stable-walking configuration in the WO-model.
The dependence of walking curvature from difference parameters is given in table 1 for selected configurations in the BO-model. Asymmetry in any single leg parameter leads to walking in circles. This implies that asymmetry in one leg parameter can be compensated for by asymmetry in other parameters (figure 3).
The dependence of the walking radius and stability on the difference parameters Δα and Δk is displayed in figure 3. The locally predicted relation between Δα and Δk for straight walking patterns from table 1 is in good agreement with simulation results for larger parameter changes throughout all stable solutions we found. We further see that stability is not necessarily reduced by asymmetric parameters but can strongly increase. A similar observation in the two-dimensional model was reported by Merker et al. .
For the highlighted configurations in figure 3, we computed the dependence of the walking curvature from energy in the BO-model. We found that for A2, A3, B1 and C1, an energy increase was accompanied by an increase in walking radius, indicating a more straight walking motion. Configuration A1, which showed straight walking for E = 816 J, did not remain straight but started veering to the right for higher energies, and veering to the left for lower energies.
Our results show that BSLIP preserves stability properties when it is extended to three dimensions, and a BO leg placement policy is applied (BO-model). Further, we find that asymmetries in the limbs typically lead to walking in circles, as one would expect. Certain combinations of leg parameters can result in straight walking, despite asymmetric leg configurations and resulting asymmetric gait pattern. For most configurations, the BO-model predicts a walking radius between about 10 and a few 100 m (figure 3). A comparison with empirical data [6,8] shows that both range and magnitude of predicted walking radius are comparable to data from human walkers, which tend to walk slightly more straight.
We find that asymmetric gait patterns of the BO-model are not inferior to symmetric ones in terms of stability. This is in line with findings from Merker et al. , who showed for the sagittal plane model that ‘considerable differences between contralateral legs can be tolerated and may even provide advantages to the robustness of the system dynamics’.
In our model, we primarily choose to align the leg with respect to an instantaneous sagittal plane (BO-model), spanned by gravity and the CoM velocity vector, instead of a world frame oriented sagittal plane (WO-model). In the WO-model, we could not expect that our model is neutral to its orientation in the world frame, which however is a prerequisite for walking in circles. Indeed, our simulations confirmed that the examined straight walking solutions are stable only if the leg is aligned with respect to velocity (i.e. in the BO-model) instead of a world frame (i.e. in the WO-model). From very similar running templates , we know that a fixed leg alignment with respect to a world frame does not result in stable running in three dimensions, but that leg orientation with respect to the current velocity does . In a previous modelling study examining postural stability of the trunk in walking , we could also show that positioning the leg with respect to the trunk instead of the vertical axis might not only preserve but increase the stability of the model. These findings give us a consistent picture stating that for stable locomotion, the body reference frame could be a superior choice for leg positioning compared with a world-based frame. This finding has potential applications in control schemes for bipedal walking robots.
When humans have access only to body-related feedback, that is, in the absence of vision and possibly acoustic guidance, they cannot maintain a straight walking path but tend to walk in circles . Boyadjian et al.  provide convincing evidence that this can be attributed to asymmetries in limb functionality instead of systematic bias in perceived body orientation. However, there are some studies that disagree [5,6,8]. In human walking, asymmetry is the norm rather than the exception , and our model shows that persistent functional asymmetry (i.e. in the BO-model, asymmetry in parameters) goes along with veering. Despite the simplicity of our model, all parameters have a physiological interpretation. While this is apparent for leg orientation angles α and β as well as for leg rest length l0, also leg stiffness k can be interpreted as strength of the leg.
There are further factors that influence veering, which cannot be captured directly by the model. For example, Kallie et al.  report that asymmetries in the motor noise [18–21] can also lead to walking in circles. In addition, distraction , sex , drugs [24,25] and deviations in stepping frequency  have been reported to have significant influence on veering. Our model cannot account for these factors, as well as for some physiological asymmetries such as different leg masses or different feet lengths. Still, some of those factors might be represented indirectly in the model, for example a longer foot could correspond to an increased angle of attack α, or the effect of a drug (e.g. pain medication for chronic pain patients) could reduce differences in exposed functional leg strength, represented by Δk, during walking.
Human walking is not exactly periodic but there is stride to stride variability, with long-range correlations over multiple strides . This can be modelled as a time-varying bias in leg parameters, resulting in a potential change of veering direction over time as found in experimental data [6,28]. Further, slight velocity changes in repeated trials might result in changes of the veering direction, as shown in our BO-model. This is especially evident for solutions with little or no veering (e.g. A1 in figures 3 and 4). Only for persons with persistent gait asymmetries such as some stroke survivors , we would expect an unique veering direction, as reported by Turton et al. .
What can we learn from this modelling study for human walking? Our model allows us to investigate the effects of various influences to the walking pattern in terms of stability and veering. We can derive predictions, such as (i) smaller step width results in less veering, (ii) the veering direction can change with changes in velocity, especially if little veering was observed before. Besides analysing and understanding normal human walking, gait models such as the proposed three-dimensional BSLIP might be used to derive rehabilitation strategies when symmetric gait cannot be achieved by providing suitable asymmetric reference patterns. We believe that template models such as the proposed three-dimensional BSLIP provide a valuable starting point for such models.
This work was supported by the BALANCE project. The BALANCE (Balance Augmentation in Locomotion, through Anticipative, Natural and Cooperative control of Exoskeletons) project is partially funded under grant no. 601003 of the Seventh Framework Programme (FP7) of the European Commission (Information and Communication Technologies, ICT-2011.2.1).
- Received June 4, 2014.
- Accepted July 2, 2014.
- © 2014 The Author(s) Published by the Royal Society. All rights reserved.