ISSN: 2320-2459

All submissions of the EM system will be redirected to **Online Manuscript Submission System**. Authors are requested to submit articles directly to **Online Manuscript Submission System** of respective journal.

**Christopher O Adeogun ^{*}**

Department of Physics, Mountain Top University, Ogun, Nigeria

- *Corresponding Author:
- Christopher O Adeogun

Department of Physics,

Mountain Top University,

Ogun,

Nigeria;

**Email:**oluwatobichris88@gmail.com

**Received:** 30-Apr-2023, Manuscript No. JPAP-23-97477; **Editor assigned:** 02-May-2023, Pre QC No. JPAP-23-97477 (PQ); **Reviewed: **16-May-2023, QC No. JPAP-23-97477; **Revised: **30-Jun-2023, Manuscript No. JPAP-23-97477 (R); **Published: **07-Jul-2023, DOI: 10.4172/2320-2459.11.4.004

**Citation:** Adeogun CO, Analysis of the Bipedal, Biphasic Forward Propulsion of the Center of Gravity of the Human Body using Theoretical Physics . J Pure Appl Phys. 2023;11:004.

**Copyright:** © 2023 Adeogun CO. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Visit for more related articles at** Research & Reviews: Journal of Pure and Applied Physics

We introduce an innovative method for the investigation of human gait, which is based on the visualisation of the vertical component of the movement of the center of mass during walking or running, in the space of the coordinates position, velocity, and acceleration of the center of mass. Collected data has been numerically approximated by the best fitting curve for a non-linear model. The resulting equation for the best fitting plane or curve in this space is a differential equation of second order. The model that we suggest is a Duffing equation with coefficients that depend on the height of a walker or runner and on the angular frequency of the oscillation. Statistics about the distribution of the Duffing stiffness depending on the speed is presented.

Duffing stiffness; Statistics; Differential equation of second order; Velocity; Centre of mass

Research on the mechanics of human gait can be of interest to different disciplines, for example sport science, medicine, and robotics. In this paper we discuss a model for the movement of the vertical coordinate of a person’s Centre of Mass (COM) during walking and running.

Human locomotion is an inherently complicated process requiring the complex integration of neural and musculoskeletal control in response to both internal and external forces. In an attempt to strip away complexity and gain an understanding of the fundamental principles underpinning human locomotion, simple mechanical models have been developed. The mechanical simplification of locomotion allows the identification of just a few key parameters that can be manipulated to examine cause and affect relationships and identify which features most influence the system^{ [1,2]}.

Blickhan suggested a linear spring-mass model for hopping in 1989^{ [3]}. Other papers followed, for example^{ [4-7]}. The motion of the center of mass is described by the equation mz_{tt}+Kz=-mg, where m is the body mass, z is the vertical deflection of the center of mass with the origin on the treadmill surface and the direction chosen upwards. The constant K is the stiffness, and g is gravitational acceleration. By z_{tt }we denote the second derivative of z, *i.e*., the vertical acceleration of the center of mass.

There have been different approaches on how to calculate leg stiffness. Blickhan’s approach uses the formula K=mω^{2}_{0}, where ω_{0} is the stride’s angular frequency of the oscillation, which, during gait, reflects the stride’s angular frequency.

The other approach for calculation of the leg stiffness is to find the ratio of F_{max}, the maximum value of the vertical ground reaction force, and ΔL, the absolute value of the leg compression,* i.e.*, K=F_{max}/ΔL. This definition of leg stiffness is used in several papers ^{[4-9]}. An overview of these two approaches is presented in^{ [6]}. There is a third approach to leg stiffness calculation based on the measurements of loss of mechanical energy by walking/running, W.

Leg stiffness K is found from the formula W= ½K(Δr)^{2}, where Δr is the shortening of a spring.1

In examining mechanical and metabolic determinants of the human walking gait, Kuo ^{[10]}, and ^{[11]} employed an anthropomorphic three dimensional, passive dynamic model, in which human legs were represented as rigid inverted pendulums with small point masses modelling each foot and a larger mass modelling the concentration of the COM at the pelvis. These studies drew on earlier models of a rigid swing leg during walking [12] and continued the view that walking and running were two distinct gaits that could not be described using similar mechanical models. This view, however, was discounted by Geyer, Seyfarth, Blickhan ^{[13]} who demonstrated that a compliant legged, spring-mass bipedal model consisting of two linear, equal and massless springs and a single COM point mass, as an extension of Blickhan’s one-dimensional model, reliably predicted ground reaction forces and COM behaviour in both human walking and running. Several subsequent studies further validated the efficacy of a bipedal spring mass model of walking^{ [4,12-15]}. While much of the twenty first century research in the field has adopted the bipedal spring-mass model and focused on adapting or adding selected elements to improve prediction accuracy for both walking and running gait mechanics, Blickhan’s spring mass model remains largely valid and has been applied, with modifications to suit certain parameters, in recent studies ^{[16,17]}.

Our mathematical model is based on the analysis of the collected data of the three dimensional movement of COM. In this paper we concentrate on the projection of the movement of COM on the vertical axis. We plan to collect more data and apply our new method for the investigation of projections on the horizontal axes the axis along the walking/running movement and the axis that is perpendicular to it.

We suggest a new approach for finding leg stiffness for the simple harmonic oscillation model, and then develop a more precise model assuming that the stiffness for a fixed speed is not a constant but depends on the displacement of the COM. The innovative idea of our method is to visualise the data for the vertical component of a motion of the COM, z(t), as a curve in the three dimensional space, (z, z_{t}, z_{tt}). Here z_{t }and z_{tt }are the first and the second derivatives of the function z(t), correspondingly. These curves are interesting by themselves as they are individual for each person and velocity. We provide examples of the variety of these curves.

In this paper we discuss the properties that are common for these curves. We apply the Fast Fourier Transform (FFT) to the data and filter out the smaller frequencies.

**Data recording**

We have collected data for six participants aged from 18 to 55 years, three males and three females. The study followed ethical protocols as per ethics requirements (HE19-239). We have measured the vertical coordinates z(t) of the COM for each participant walking or running on the treadmill. The markers we used were the Left and Right PSIS and ASIS, then we computed the average of all four. The data was collected for different integer velocities, at 100 frames per second, over 10 seconds, for each velocity, using an 8 camera, Qualisys Motion capture system with the COM reconstructed using a pelvic marker set within Visual3D.

Using MATLAB, we have visualised the data as curves in the three dimensional space (z, z_{t}, z_{tt}). We have filtered out high frequency oscillations, such as noise and individual features, using the fast Fourier transform function in MATLAB. For FFT threshold equal to 0.3, the curves are close to an ellipse. The direction of the motion along the curve is found, for example, in the following way. Find the point with the greatest z-coordinate. This is the highest point of the centre of mass during the gait cycle. The velocity at this point is equal to zero. Hereafter the movement goes down, *i.e*., the velocity z_{t }becomes negative. In both pictures of **Figure 1**, from those perspectives, the movement is clockwise.

**Data interpretation**

The pictures in the space with the coordinates position, velocity, and acceleration, are rich in information. For example, **Figures 2 and 3** show the data for walking (4 km/hr) and running (9 km/hr) of a participant, correspondingly. Here we look into the projection to the plane z, z_{tt},* i.e*., the horizontal axis shows the position of the COM (in meters) and the vertical axis shows the acceleration (in m/sec^{2}). The part CD of both curves corresponds to the phase of the gait when a foot touches the surface of the treadmill. In the case of walking it is a flat line; vertical acceleration is close to zero. In the case of running, acceleration is diminishing because of braking during initial foot contact. The phase DE corresponds to the propulsion during toe off, where acceleration (and consequently force) increases. On the segment EA acceleration diminishes, turns to zero when the the COM reaches its average position, and is minimal at the point A. The minimum acceleration for the walking curve is -2 m/sec^{2}, the minimum acceleration for the running curve is about -10 m/sec^{2}, *i.e*., close to the gravitation constant g. The acceleration is less during walking than during running because the body is always in contact with the ground whereas during running there is a flight phase. The arc AB on **Figure 3** corresponds to the flight phase of running. This part of the curve is more complicated than just a flat constant z_{tt}=-g, because it is smoothly connected with the rest of the curve. Some information that we get from these curves is common for all participants and walking/running speeds, but some features are individual for example not each participant has the “flight” component AB at running speeds, due to individualised transitions between walking and running gait patterns.

Our aim in the future is to collect data for more participants and compare individual properties of the curves, see Section 9. In this paper we concentrate on their common properties, and suggest three models, based on differential equations. We build our three models based on data, purely numerically and mathematically. The first model is well-known, it is a harmonic oscillation. Our novelty here is the method we develop for computing leg stiffness.

**Figure 2:** Data for walking, 4 km/hr. The horizontal axis shows the position and the vertical axis shows the acceleration of the COM. The part CD corresponds to the phase of the gait when a foot touches the surface, DE corresponds to the propulsion during toe off, during EA the COM moves upwards and the acceleration diminishes.

**Figure 3:** Data for running, 9 km/hr. The horizontal axis shows the position and the vertical axis shows the acceleration of the COM. The part CD corresponds to the phase of the gait when a foot touches the surface, DE corresponds to the propulsion during toe-off, during EA the COM moves upwards and the acceleration diminishes. The arc AB corresponds to the flight phase.

**Modelling gait as harmonic oscillator**

**Best fitting plane and the interpretation of its coefficients:** We begin with visualisation of the data. Consider a curve in the three dimensional space (z, z_{t}, z_{tt}) that represents smoothed gait data, as for instance the red line on **Figure 1** (FFT=0.3). This curve is not quite flat, but it is close to a flat ellipse with the centre at z=z_{0}, z_{t}=z_{tt}=0. The equation of a plane in three dimensional space (z, z_{t}, z_{tt)} that passes through the point (z_{0}, 0, 0) is

az_{tt}+bz_{t}+c(z−z_{0})=0, (1)

where a, b, and c are constants, and z-z_{0} is the vertical displacement of the average of the center of mass. The zero for z(t) is chosen to be on the level of the treadmill, and the z-axis direction is upwards.

We find the coefficients a, b, and c in equation (1) numerically, by finding the best fitting plane for the smoothed data curve. The equation for the plane is interpreted as a second order differential equation that describes a damped driven harmonic oscillator.

The coefficients in equation (1), in the context of the spring model, have the meaning of a=m, which is the mass of the oscillating object, b=ν, which is the viscosity of the spring, and c=K, which is the spring stiffness,

mz_{tt}+νz_{t}+K(z−z_{0})=0. (2)

If we consider a stable gait, we need to assume, in this linear model, that the viscosity is zero. Equation (2) turns to

mz_{tt}+K(z−z_{0})=0, (3)

Which is Hook’s Law that describes a movement of a spring with stiffness K^{ [18]}.

F(t)=−K(z(t)−z_{0}). (4)

**Where the gravitation constant is hidden**

We note that the gravitation constant is included in equation (2) implicitly. Namely, we can rewrite the equation as

mz_{tt}+νz_{t}+K(z−z_{st})=−mg. (5)

The coordinate z_{st }is the average of the vertical coordinate of the centre of mass of a standing body but in a walking/running posture. It is not exactly the same as the coordinate of COM in a standing position. The relation between z_{0} and z_{st} is calculated from Equations (2) and (5) and is z_{0}=z_{st}−mg/K.

**Modelling by a non-linear homogeneous differential equation**

**Best fitting curve, interpreted as a non-linear second order differential equation: **Non-linear gait dynamics has been discussed, for example, in^{ [19]}. The author remarks that stride to stride fluctuations, which are often considered to be noise, actually convey important information. With the aim to describe these fluctuations, we refine the method used for the harmonic oscillation model. We approximate the movement of the COM during walking or running by a Duffing equation,* i.e*., a homogeneous non-linear second order differential equation. We write

mztt+ν(t)z_{t}+K(t)(z−z_{0})=0, (6)

Where, unlike the harmonic oscillation model, we consider stiffness and viscosity to not be constant, but functions depending on time, K(t) and ν(t). This makes sense, as stiffness and viscosity depend on a phase of a stride. For example, the slope on **Figures 2 and 3** is changing, and the slope in the plane z, z_{tt} reflects stiffness:
$\frac{\mathrm{K(t)}}{m}$
at a time t. We divide both sides by m and consider the approximation of the functions v(t)=
$\frac{\mathrm{\nu (t)}}{m}$
and K(t)=
$\frac{\mathrm{K(t)}}{m}$
by cubic polynomials,

z_{tt}+(v_{1}z_{t}+v_{2}z^{2}_{t}+v_{3}z^{3}_{t})+(k_{1}(z−z_{0})+k_{2}(z−z_{0})^{2}+k_{3}(z−z_{0})^{3})=0. (7)

As the value z_{0} is not known beforehand, we look for the best fitting curves of the form of least squares:

z_{tt}+v_{1}z_{t}+v_{2}z^{2}_{t}+v_{3}z^{3}_{t}+k_{1}z+k_{2}z^{2}+k_{3}z^{3}+C=0. (8)

The best fitting curve in the chosen coordinates is a differential equation of second order. We visualize the solution of this differential equation as a curve in the same coordinate system as the data curve. See, for example, **Figures 4, 7, 8**, where the red curve is the data curve and the blue curve is the solution of the differential equation given by the best fitting curve. Initial conditions are the values z(0) and z_{t}(0). Equation 8 is equivalent to equation (7) after expanding brackets and renaming coefficients. From our data, the viscosity coefficients v_{1}, v_{2}, v_{3 }turned out to be zero if the FFT threshold is large enough, for instance if it is 0.3. This is not surprising because we have assumed no excitation force (homogeneous equation), and the participants have maintained stable walking/running. If the viscosity is non-zero, the curve is a stable spiral approaching z_{0}, as in **Figure 8**. Scaled mean squared error between the Duffing equation output and the observed data was computed. As the values z, z_{t} and z_{tt }have different units, each value was scaled by dividing by the difference between the minimum and maximum values on the corresponding axis.

**Finding and analysing the fixed points**

We are also interested in the fixed (equilibrium) points ^{[20]} for the curves described by equation (8). We first rewrite the second order differential equation (7) as a system of first order differential equations:

z_{t}=z_{1}, (9)

(z_{1})_{t}=−v_{1}z_{1}−v_{2}z^{2}_{1}−v_{3}z^{3}_{1}−k_{1}z−k_{2}z^{2}−k_{3}z^{3}−C. (10)

We remind that for stable gait v_{1}=v_{2}=v_{3}=0. The fixed points are found by solving the equilibrium equations (see for example (3))

z_{1}=0, (11)

−k_{1}z−k_{2}z^{2}−k_{3}z^{3}−C=0, (12)

*i.e*., the fixed points are the roots of the cubic equation k_{1}z+k_{2}z^{2}+k_{3}z^{3}+C=0. The solutions to equation (8) are plotted in the same coordinate system as the data curves, see for example **Figure 4.** The solution curves are stable if we set one fixed point equal to zero, *i.e*., C=0. Then the two other fixed points occur at z_{0} (the center of the closed curve, average coordinate of the center of mass during walking/running), and at h, which is approximately 2z_{0} for stable walking/running. The differential equation equation (8) becomes

z_{tt}+kz(z−z_{0})(h−z)=0. (13)

The connection between the coefficients of equation (8) and equation 13 is

k_{1}=−kz_{0}h, k_{2}=k(h+z_{0}), k_{3}=−k, k>0, and v_{1}=v_{2}=v_{3}=0. (14)

Numerical computations show that, for stable gait, the fixed point z=z_{0} is a center, while the fixed points z=0 and z=h are saddles.

**Interpretation of the parameters in the model**

Equation (13) shows that we can model the movement of the center of mass by a Duffing equation, up to a constant k, knowing only h and z_{0}. The values h and z_{0} are close to the height of a person and to the average coordinate of the COM in motion, correspondingly. The gravitation constant g is involved implicitly in the differential equation, in a similar way as in the harmonic oscillation model.

The coefficient k does not have the meaning of the square of the angular frequency ω²=
$\frac{\mathrm{stiffness}}{\mathrm{Mass}}$
, as in a linear case, but k behaves in a similar way: It increases with an increase of the walking/running speed. We call this constant the Duffing stiffness. The meaning of the coefficient k is found from the following consideration.

We rewrite equation (13) as

z_{tt}+k(−(z−z_{0})^{3}+(h−2z_{0})(z−z_{0})^{2}+z_{0}(h−z_{0})(z−z_{0}))=0. (15)

For z close to z_{0} the linear approximation at z_{0} is

z_{tt}+kz_{0}(h−z_{0})(z−z_{0})=0. (16)

Comparing this equation with Hook’s Law (equation (4)) we get an expression that relates the angular frequency ω with the coefficient k and the values z_{0 }and h:

ω^{2}≈kz_{0}(h−z_{0}). Hence, the Duffing stiffness

K≈ | ω² z _{0}(h-z_{0}) |
(17) |

Case h = 2z_{0}

Assuming that the center z_{0} is approximately in the middle between two saddle points, 0 and h, we get h=2z_{0}, and Equation (15) becomes

z_{tt}+kz^{2}_{0}(z−z_{0})−k(z−z_{0})^{3}=0, (18)

where k>0, and

ω^{2}

K≈ | ω²z _{0}² |
(19) |

Equation (18) is the Duffing equation for a softening oscillator ^{[3]},* i.e.*, the stiffness diminishes with the displacement.

**Example based on the collected data**

**Figure 4** shows the observed smoothed data curve (red) and the curve corresponding to the solution of Equation (13) (blue) for one of the participants, together with the fixed points. The differential equation that describes the movement of the centre of mass of this participant running at the speed 8 km/hr is z_{tt}+279(z−0.9)(1.8−z)=0. We have rounded |k_{3}|=k to an integer, and z_{0} and h to the first decimal place (**Figure 4**).

We can check that equation (4) and the approximation equation (17) are realistic by computing the number of strides in a second and comparing it with the real value.

From k=279, z_{0}=½ h=0.9 we can compute ω ≈
$\surd \frac{K}{\mathrm{z\u2080\xb2}}$
≈ 18.56, *i.e.,* the number of strides in a second is
$\frac{\omega}{\mathrm{2\pi}}$
≈ 2.95. This participant made 24 strides in 10 secs, *i.e*., the real number of strides in a second is 2.4.

**Statistical analysis**

**UNE data: ****Figure 5** (left) shows that the Duffing stiffness depends on the speed of walking/running for each participant. By visually inspecting the graph, we see some outlier points. To further examine these, we calculated the standard deviation of Duffing stiffness for each participant, for the walking data. Two participants had standard deviations >400, while the rest of the participants had standard deviations <100. Our observation of these two participants during the data recording showed that they were uncomfortable with some of the speeds, and they reported having no experience with treadmills.

Therefore, we excluded these two participants from all further analyses. **Figure 5** (right) shows the data for the remaining four participants.

The behaviour of the Duffing stiffness is different for walking (3-8 km/hr) and running (9-14 km/hr). Therefore, we separately fit data for walking and running speeds. We fitted a linear mixed effect model for the walking data. We included speed as the fixed effect and Duffing stiffness as the dependent variable. We allowed both the slopes and the intercept to vary across participants. The model showed a slope estimate of ̂β = 41.4, t=9.1, p=0.004. The intercept was -14.6, p=0.5.

For the running data, we performed an equivalent analysis. Here, the slope associated with speed was not significant, though trending in a positive direction, ̂β =13.9, t=3.0, p=0.05. The intercept was 245.6, t=4.0, p=0.02.

**Comparing UNE dataset against public datasets**

Next, we compare our data against public datasets for 42 walking participants ^{[21]}, and 29 running participants. As these sources were set up slightly differently than the UNE data, in order to match the data, we estimated the COM during walking or running, using the midpoint of the ASIS markers and subtracting the lowest value for the heels marker. **Figure 6 **shows the relationship between speed and Duffing stiffness.

First we removed one participant with only one data point. Then, to identify outliers, we calculated the standard deviation in Duffing stiffness for each participant and excluded outliers for the walking and the running datasets, separately. For the walking data ^{[22]}, the standard deviation ranged from 39.2 to 3183.5. We chose to remove all participants with a standard deviation >100. This left us with 12 participants for the walking data. We proceeded to fit the data with a Linear Mixed Effect model, akin to the UNE data. Again, we get a significant effect of speed on Duffing stiffness, ̂β=36.7, t=11.8, p<0.0001. The intercept is not significant, 24.0, p=0.09. For the running data^{ [21]}, we removed three out of 29 participants with standard deviations >100. The Linear Mixed Effect model showed a significant effect of speed, ̂β=13.5, t=14.0, p<0.0001. The intercept is 188.1, t=19.7, p<0.0001.

Finally, we compared the two datasets against each other, to determine whether there are any overall differences in Duffing stiffness, or whether there is difference in the relationship between speed and Duffing stiffness. We created two models, one for running and one for walking, including both of the trimmed datasets. In the Linear Mixed Effect model, Duffing stiffness acted as the dependent variable, and the fixed effects were speed, dataset UNE versus ^{[21,22]}, with UNE acting as a baseline. We allowed the slope and the intercept to vary across participants. For walking, the effect of speed was, again, significant, ̂β=36.7, t=12.4, p<0.0001. However, the effect of dataset was not significant, ̂β=43.1, t=1.5, p=0.2, nor was the interaction between dataset and speed, ̂β=5.3, t=0.9, p=0.4. Similarly, in the running dataset, the effect of speed was ̂β=13.5, t=13.2, p<0.0001. Neither the effect of dataset, nor the interaction of dataset and speed were significant, ̂β=27.5, t=0.9, p=0.4 and ̂β=3.4, t=1.1, p=0.3, respectively.

In summary, we found, overall, a significant positive relationship between speed and Duffing stiffness. For the running UNE data, the slope was not significant; however, when we combined the two datasets, the running slope was significant, and we found no interaction. Thus, the lack of significance in the UNE running data may be a result of low statistical power. We found no main effect of dataset, nor an interaction between dataset and speed. Thus, we find no evidence of a difference across datasets.

**Identifying and examining outliers**

From a practical perspective, an interesting aspect is participants whose data deviates from the fitted model. Here, we defined outliers based on standard deviation. This is the simplest method, which can easily be applied by a sport scientist. In our analyses, we drew a somewhat arbitrary threshold, where we treated all participants with a standard deviation >100 as outliers. In practice, this can act as a sign for the practitioner to follow up on the walker or runner if their standard deviation is high.

The reasons for high standard deviations could be different. For example, a typical problem for uncomfortable speeds less than 3 km/hr is an additional loop, as on the red line in **Figure 7.** The Duffing equation does not take into account the loop, as the blue approximation curve demonstrates. The other example of an outlier is illustrated in **Figure 8**, where the curve brakes down into two parts, corresponding to the left or to the right leg. Red lines represent the data, blue lines are the solutions of the differential equation (**Figure 8**).

A third reason for outliers that we noticed is an instability of walk, when each step varies in the amplitude and in the average height.

**Duffing equation with viscosity and excitation force**

**Motivation: **Approximation by a homogeneous equation with zero viscosity does not take into account the asymmetry caused by damping and excitation forces.

For example, **Figure 9** shows a smoothed data curve for a running participant (14 km/hr, FFT threshold=0.3) in the phase plane (z, z_{t}). This curve is close to an ellipse, and the absolute value of the slope of the main axis of the ellipse, CA, is equal to v=ν/m

if we assume that the viscosity is a constant. (If the viscosity ν=0, the slope is vanishing, and the corresponding axis is horizontal.) The symmetry with respect to z_{t}=0 is disturbed by the slope.

Now we compare this slightly asymmetrical typical data curve with the symmetrical energy level curves, see **Figure 10**. Solutions of equation (13) with different initial values of z(z_{t}=0) give a set of energy level curves. The red curve is the data curve, the movement occurs in the direction ABCD. Energy is gained twice in each cycle (gait) in a phase of a “step”, BC, and in a phase of a “fall” (in general not free fall), DA. The maximum of the energy occurs at points A and C, the minimum occurs when at points B and D. There is a natural desire to find an approximation that considers the viscosity and the restoring force.

**Looking for the best fitting curve for a non-homogeneous differential equation**

We were looking for the best fitting curve in different forms, as for example z_{tt}+vz_{t}+C_{1}z+C_{2}z^{2}+C_{3}z^{3}=f cos(Ωt−ϕ), where the term vz_{t} is the linear approximation of the damping force divided by mass, C_{1}, C_{2}, and C_{3 }are coefficients of the cubic nonlinearity. The constant f=F/m, where F is the amplitude of the excitation force and m is the mass, Ω is the angular frequency of the excitation force. The gravitation constant g is involved in the equation in a similar way as for the two already discussed models. We have also tried to find the best fitting curve in the forms based on equations (13), (18). As in the previous models, the equation for the best fitting curve was interpreted as a second order differential equation. The solution of the equation was computed and plotted in the same system of coordinates as the data curve. However, the balance between the viscosity and excitation force was too delicate, and most solution curves became stable spirals. The balance worked with an FFT threshold of 0.5, at least for some curves, based on equation (18),

z_{tt}+vz_{t}+kz^{2}_{0} (z−z_{0})−k(z−z_{0})3=f cos(Ωt − ϕ). (20)

**Example**

In **Figure 11** the red data curve is approximated by a homogeneous differential equation (13) (blue curve), and by Equation (20) containing the damping and restoring forces (green curve). The mean squared error for the homogeneous equation is 0.01. The non-homogeneous equation gives a better approximation with the mean squared error 0.002.

**Plans for future research**

**Classification of curves:** One of the ideas for future research is to compare and classify curves for different participants. The curves can have individual features, see for example **Figure 12** (FFT threshold 0.03).

Also the left and right leg strides are almost indistinguishable if we chose an FFT threshold 0.3. But for some participants, we observe a difference in the slopes for even and odd strides if the FFT threshold is smaller, suggesting sensitivity in identification of gait asymmetries. In **Figure 13** the green curve (FFT threshold=0.03) is split into two parts. Recall that different slopes correspond to different angular frequencies. The black curve shows the same data smoothed with an FFT threshold of 0.3. Bipedal models were investigated using a rather mechanical approach in several papers as for example in Hartmut Geyer, et al. and Horst-Moritz Maus, et al. ^{[17,23]}.

**Horizontal components**

In this paper we have modelled only the vertical component of the movement of the COM. We plan to also investigate the horizontal components in the future. Our hypothesis is that the horizontal component that is perpendicular to the forward movement can be approximated by a second order differential equation with negative stiffness^{ [24]}.

**Investigating possible reasons for injuries**

Our model gives an insight into the dynamics of human gait. The Duffing stiffness at a particular speed is approximately the same for each person who maintains a stable gait. If the Duffing stiffness differs, it points to some uncomfortableness in walking or running. The reason for this uncomfortableness can be investigated by looking at the curves closely. We plan to investigate possible reasons for injuries by determining how the stress is generated ^{[25,26]}.

One possible idea is to investigate why female runners have more frequent ACL (Anterior Cruciate Ligament) tears than men^{ [27]}.

We study the human gait modelled by the movement of the centre of mass of the test person. This is an example of a biological process which can be considered as a periodical dynamic system. Roughly, this movement behaves in a similar way to a vibrating mass suspended on a spring, but it is more complex. The vertical component of the movement during walking or running can be visualised as an oscillogram: a graph of the position as a function of time. We suggest a visualisation of the data in 3D space, where the coordinates describe position, velocity, and acceleration. Our new visualisation method allows us to model the movement of a person’s centre of mass by a nonlinear differential equation. The resulting curve for an ideal spring-mass movement, without viscosity or external force, is an ellipse in the suggested 3D space. The shape of the data curve shows at which position an additional force was applied, or the movement slowed down. Some deviations are common for all test persons and others are different. In the future we plan to investigate the reasons for these deviations, such as different running techniques or the presence of injuries.

We have introduced an innovative method for the investigation of human gait. This method is based on the visualisation of the vertical component of the movement of the center of mass during walking or running, in the space of the coordinates position, velocity, and acceleration of the center of mass. This method allows us to compute leg stiffness for an established harmonic oscillator model in a new way, as the slope of the best fitting plane. We determined some differences and some common features for walking and running curves. In this paper we concentrated on features that are common for walking and for running. We suggested a model by a non-linear homogeneous differential equation. For this model we need to know the average of the COM, z_{0}, during the gait motion, and a coefficient K that depends on a participant and on a walking/running speed. Assuming h=2z_{0}, we get the Duffing equation for a softening oscillator

mz_{tt}+Kz^{2}_{0} (z−z_{0})−K(z−z_{0})^{3}=0. (21)

Statistics showing the dependence of the Duffing stiffness on the speed was shown separately for walking and for running speeds, using both our collected dataset and an open dataset.

We also had a partial success in approximation of the movement by a second order non-linear non-homogeneous equation

mz_{tt}+νz_{t}+Kz^{2}_{0} (z−z_{0})−K(z−z_{0})^{3}=F cos(Ωt−ϕ). (22)

We believe that our new method is a tool that could lead to other interesting and novel results.

- Sharon R Bullimore, et al. Consequences of forward translation of the point of force application for the mechanics of running. J Theor Biol. 2006;238:211–219.
[Crossref] [Google Scholar] [PubMed]

- Susanne W Lipfert, et al. A model experiment comparison of system dynamics for human walking and running. J Theor Biol, 2012;292:11–17.
[Crossref] [Google Scholar] [PubMed]

- Michael J Brennan, et al. The Duffing equation: Nonlinear oscillators and their behaviour. John Wiley and Sons. Newyork, USA. 2011.
- Yvonne Blum, et al. Effective leg stiffness in running. J Biomech. 2009;42:2400–2405.
[Crossref] [Google Scholar] [PubMed]

- Ferris DP, et al. Running in the real world: Adjusting leg stiffness for different surfaces. Proc Biol Sci. 1998;265:989–994.
[Crossref] [Google Scholar] [PubMed]

- Nikooyan AA, et al. Mass spring damper modelling of the human body to study running and hopping an overview. Proc Inst Mech Eng H. 2011;225:1121–1135.
[Crossref] [Google Scholar] [PubMed]

- Amy Silder, et al. Running with a load increases leg stiffness. J Biomech. 2015;48:1003–1008.
[Crossref] [Google Scholar] [PubMed]

- Sharon R Bullimore, et al. Ability of the planar spring mass model to predict mechanical parameters in running humans. J Theor Biol. 2007;248:686–695.
[Crossref] [Google Scholar] [PubMed]

- Claire T, et al. Farley and Octavio Gonzalez. Leg stiffness and stride frequency in human running. J Biomech. 1996;29:181–186.
[Crossref] [Google Scholar] [PubMed]

- Arthur D Kuo, et al. Stabilization of lateral motion in passive dynamic walking. Int J Robot Res. 1999;18:917–930.
- Maxwell Donelan J, et al. Simultaneous positive and negative external and mechanical work in human walking. J Biomech. 2002;35:117–124.
[Crossref] [Google Scholar] [PubMed]

- Hyunhwa Hong, et al. Spring like gait mechanics observed during walking in both young and older adults. J Biomech. 2013;46:77–82.
[Crossref] [Google Scholar] [PubMed]

- Chang Kuen Jung, et al. Compliant bipedal model with the center of pressure excursion associated with the oscillatory behaviour of the center of mass reproduces the human gait dynamics. J Biomech. 2014;47:223-229.
[Crossref] [Google Scholar] [PubMed]

- Seyoung Kim, et al. Leg stiffness increased with speed to modulate gait frequency and propulsion energy. J Biomech. 2011;44:1253–1258.
[Crossref] [Google Scholar] [PubMed]

- Hyunggwi Song, et al. A springy pendulum could describe the swing leg kinetics of human walking. J Biomech. 2016;49:1504–1509.
[Crossref] [Google Scholar] [PubMed]

- Christian Ludwig, et al. Multiple step model experiment matching allows precise definition of dynamical leg parameters in human running. J Biomech. 2012;45:2472–2475.
[Crossref] [Google Scholar] [PubMed]

- Horst-Moritz Maus, et al. Constructing predictive models of human running. J R Soc Interface. 2015;12:20140899.
[Crossref] [Google Scholar] [PubMed]

- Robert Hook, et al. Micrographia: or, Some physiological descriptions of minute bodies made by magnifying glasses, with observations and inquiries thereupon. Printed by J. Martyn and J. Allestry. London, England. 1665: 246.
- Jeffrey M Hausdorff, et al. Gait dynamics, fractals and falls: Finding meaning in the stride to stride fluctuations of human walking. Hum Mov Sci. 2007;26:555–589.
[Crossref] [Google Scholar] [PubMed]

- Dominic Jordan, et al. Nonlinear ordinary differential equations. 4th Edition. Oxford University Press, United States. 2007.
- Claudiane A, et al. A public dataset of running biomechanics and the effects of running speed on lower extremity kinematics and kinetics. Peer J. 2017;5:3298.
[Crossref] [Google Scholar] [PubMed]

- Claudiane A, et al. A public dataset of over ground and treadmill walking kinematics and kinetics in healthy individuals. Peer J. 2018;6.
[Crossref] [Google Scholar] [PubMed]

- Hartmut Geyer, et al. Compliant leg behaviour explains basic dynamics of walking and running. Proc Biol Sci. 2006;273:2861-2867.
[Crossref] [Google Scholar] [PubMed]

- Blickhan R, et al. The spring mass model for running and hopping. J Biomech. 1989;22:1217-1227.
[Crossref] [Google Scholar] [PubMed]

- Georges Dalleau, et al. The spring mass model and the energy cost of treadmill running. Eur J Appl Physiol Occup Physiol. 1998;77:257–263.
[Crossref] [Google Scholar] [PubMed]

- Simon Mochon, et al. Ballistic walking. J Biomech. 1980;13:49–57.
[Crossref] [Google Scholar] [PubMed]

- Raymond Stefani, et al. Kinesiology analysis of athletics at the ancient Olympics and of performance differences between male and female Olympic champions at the modern games in running, swimming and rowing. Athens J Sport. 2017;4:123–138.