Stochastic approach to modelling of near-periodic jumping loads

Mechanical Systems and Signal Processing 24 (2010) 3037–3059

Contents lists available at ScienceDirect

Mechanical Systems and Signal Processing journal homepage: www.elsevier.com/locate/jnlabr/ymssp

Stochastic approach to modelling of near-periodic jumping loads V. Racic , A. Pavic 1 Department of Civil and Structural Engineering, University of Sheffield, Sir Frederick Mappin Building, Sheffield S1 3JD, United Kingdom

a r t i c l e in fo

abstract

Article history: Received 10 June 2009 Received in revised form 9 February 2010 Accepted 27 May 2010 Available online 2 June 2010

A mathematical model has been developed to generate stochastic synthetic vertical force signals induced by a single person jumping. The model is based on a unique database of experimentally measured individual jumping loads which has the most extensive range of possible jumping frequencies. The ability to replicate many of the temporal and spectral features of real jumping loads gives this model a definite advantage over the conventional half-sine models coupled with Fourier series analysis. This includes modelling of the omnipresent lack of symmetry of individual jumping pulses and jump-by-jump variations in amplitudes and timing. The model therefore belongs to a new generation of synthetic narrow band jumping loads which simulate reality better. The proposed mathematical concept for characterisation of near-periodic jumping pulses may be utilised in vibration serviceability assessment of civil engineering assembly structures, such as grandstands, spectator galleries, footbridges and concert or gym floors, to estimate more realistically dynamic structural response due to people jumping. & 2010 Elsevier Ltd. All rights reserved.

Keywords: Vibration serviceability Human–structure dynamic interaction Jumping forces Stadia Grandstands

1. Introduction Predicting dynamic performance of civil engineering structures due to crowd-induced loading is an increasingly critical aspect of vibration serviceability design process for assembly structures, such as grandstands, spectator galleries and concert halls, which are becoming more slender and lighter than ever before [1]. Broadly speaking, the procedure for predicting dynamic response of a structure under consideration involves specifying design load and determining dynamic properties of the structure in terms of modal mass, stiffness and damping obtained from a structural model. There are degrees of uncertainty and latitude in each element and different degrees of guidance on procedures. However, of all these elements, determining the design load has the greatest uncertainty. A vast majority of relevant design guidelines around the world has recognised jumping as the most important type of crowd-induced load on an assembly structure [2,3]. This is because it is the most severe and frequent to happen in practice. There have been cases of impaired vibration serviceability under crowds jumping on footbridges (due to vandal loading as described by Zivanovic et al. [4]), grandstands [5] and concert halls [6]. Therefore, there is a need for a reliable model of jumping forces to facilitate vibration serviceability checks of these structures. It is now well established that jumping by individuals and crowds is not a deterministic and ‘perfectly’ periodic, but rather a stochastic and narrow band process [7,8]. This is because of the so called intra- and inter-subject variability between naturally imperfect humans taking part in the jumping. However, the vast majority of jumping

 Corresponding author. Tel: + 44 114 222 5727; fax: +44 114 222 5700.

E-mail addresses: v.racic@sheffield.ac.uk (V. Racic), a.pavic@sheffield.ac.uk (A. Pavic). Tel.: + 44 114 222 5721; fax: +44 114 222 5700.

1

0888-3270/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.ymssp.2010.05.019

3038

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

Nomenclature fs

Dt Fi N Ti 0 Ti Ti 0 Ti AT A0 T

s2T

ST(fm) ST(f)

sampling rate time step measured force number of jumping pulses real periods of jumping cycles synthetic periods of jumping cycles variations of real jumping periods variations of synthetic jumping periods Fourier amplitudes of T i series 0 Fourier amplitudes of T i series variance of T i the ASD of T i fit of the ASD ST(fm)

Iw,i I0 w,i Is,i

real weight normalised impulses synthetic weight normalised impulses real unity normalised impulses ai real scaling factors a0 i synthetic scaling factors r0 and r1 autoregression coefficients DIw,i autoregression error oi angular frequency of rotation y angular coordinate Zi(t) and Zi(y) shapes of unity normalised jumping pulses Wj and Air Gaussian weights cj, tir and yr Gaussian centres dj, bir and br Gaussian widths FFT fast Fourier transform Z correlation coefficient

models in the published literature and design guidelines assume that jumping forces are deterministic and periodic, presenting their modelling as sinusoids capable of exciting pure resonance of structures [9]. This assumption is often over-conservative. The resulting excessive levels of vibration predicted can therefore preclude efficient design. The same problem was observed in the past for walking as another type of human-induced dynamic force [10]. Recently, there have been serious attempts to resolve this issue by developing a new generation of near-periodic mathematical models for jumping forces induced by a single person, groups or crowds [7,8]. Sim’s modelling is the most recent and relevant step in the right direction, but has shortcomings which need addressing. To fit individual jumping pulses, it utilises a cosine squared function symmetric about a vertical axis through its peak. Therefore, it does not take into account the lack of pulse’s symmetry observed in the real force measurements [11]. Also, only a limited number of jumping frequencies are modelled, the majority of them being in the range of moderate and fast jumping rates. As will be seen later, the key reason for this is a considerably more complex shape of the pulses generated for jumping at slow motion, featuring two peaks, which cannot be modelled using a squared cosine function. Subsequent pulses for a single jumper are shifted in time by different amounts determined by autoregressive modelling based on approximately 1000 force records. The same procedure is repeated for additional persons yielding a model for a crowd dynamic loading as a simple sum of individual synthetic forces. It is clear that a key ingredient of a reliable crowd jumping model is an accurate model for a single person jumping. Therefore, considering advances made by Sim et al. [8], there is still a need to develop a good quality jumping model for a single person applicable to the whole range of possible jumping frequencies. This model has to be narrow band and random taking into account all aspects of the inter- and intra-subject variability. This paper offers a solution to this problem by proposing a novel pulse modelling technique which makes use of an extensive database of measured jumping forces gathered at the University of Sheffield in 2008 and 2009. This database has an impressive coverage of 14 jumping frequencies in the range 1.4–2.8 Hz, which are observed to be comfortable for rhythmic, repeated motion [12]. The technique proposed is motivated by an existing procedure for modelling recordings of electric waves being generated during heart activity [13], generally known as electrocardiogram (ECG) signals. These signals have similar near-periodic features as the force signals due to jumping. Therefore, the paper utilises technology for measuring jumping forces, which has been available for more than 20 years, and addresses the current lack of an appropriate mathematical modelling to simulate accurately what is being measured. Moreover, the paper presents a logical extension of the force models published recently by the authors [11,14]. Having analysed additional above mentioned jumping force records, major differences between the models can be observed in the selection of modelling parameters, the way they are mathematically described and mutually related to simulate more reliably the measured data.

2. Background reviews Two periods can be clearly identified in the history of developments in jumping force modelling. First, consisting of models developed prior to about 2004, utilised methods for modelling perfectly periodic signals. Subsequent methods recognised the need for a narrow band modelling. These are also characterised by rapid advancements in computational power (whereby complex models can be better utilised), structural analysis, and advanced measurement technology which became affordable. The next two sections present a critical overview of the force models which marked these two periods.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3039

2.1. Periodic models Dynamic forces generated by people jumping on a structure are commonly determined by direct measurement of the interface forces between the feet and the structure itself, hence they are known as ground reaction forces (GRFs). The measured force signal for individuals jumping is typically a series of distinctive pulses (Fig. 1), which are the reaction to the force the body exerts on the supporting ground during ‘contact phase’ of jumping. The pulses are separated by zero-force intervals which indicate ‘aerial phases’ of jumping when both feet leave the ground. A common practice until now is to idealise a continuously measured jumping force signal as periodic with the period being the average time between two consecutive jumps. This means that actual forces due to continuous jumping can be represented by a sequence of identical pulses on a jump-by-jump basis. A number of authors [15,16] fitted simple half-sine function to the average measured jumping pulse. Fig. 2 shows an example of directly measured pulse extracted from Fig. 1 and the associated half-sine model based on modelling parameters suggested elsewhere [16]. As illustrated in Fig. 2, the symmetric half-sine pattern cannot fit a typically asymmetric shape of the real jumping pulse with visually apparent good matching with measured data [11]. For dynamic analysis, a set of identical half-sine pulses can be represented more efficiently if expressed in terms of Fourier series with the fundamental harmonic having frequency identical to the jumping rate [9,15]. Such models can be found in the current British Standard BS 6399-1 [2] and Commentary D of the National Building Code of Canada [3]. In these, jumping loads are defined for use in the design of structures likely to be subjected to significant vertical occupant motion, such as footbridges, grandstands, and concert or gym floors, where human-induced serviceability issues may govern design. The main disadvantage of this method is that it requires many terms to describe satisfactorily the original half-sine approximation. For ease of use, only the first three harmonics (including six coefficients, i.e. three amplitudes and three phases) are typically considered [3,9,16]. However, even the sum of the first six Fourier harmonics (12 coefficients), which

Force [N]

3000

2000

1000

0 0

5

10 Time [s]

15

20

Fig. 1. Typical measured GRF signal generated by a single person jumping at 2 Hz.

3500 measured data

3000

half-sine model

Force [N]

2500 2000 1500 1000 500 0 0

0.1

0.2

0.3

0.4

0.5

Time [s] Fig. 2. Measured jumping pulse extracted from Fig. 1 vs. corresponding half-sine model based on modelling parameters suggested by Ellis and Ji [16]. The shaded areas represent the difference.

3040

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

is the maximum number reported in the literature [16,17], cannot match adequately enough the original half-sine forcing function for all contact times (Fig. 3), let alone real jumping force measurements (Fig. 1). Brownjohn et al. [10] showed that there were significant differences between the resonant responses due to real walking forces (which are narrow band and therefore ‘imperfectly’ periodic) and the equivalent periodic simulation. The effect is more pronounced for higher harmonics of walking loading where the simulated vibration response is regularly overestimated. This issue has not been researched in great detail for jumping excitation, so the available literature is limited. Nevertheless, a similar analysis can demonstrate that ‘synthetic’ periodic jumping forces generate a considerably higher resonant response compared with its more realistic narrow band counterpart. The reason for this becomes evident if, instead of an average jumping pulse, a window comprising a number of consecutive jumping pulses is used to derive their Fourier amplitude spectra. Fig. 4 shows the dominant harmonics around 2, 4 and 6 Hz of real jumping at an average rate of 2 Hz. Other spectral lines, having lower amplitudes around the dominant harmonics, are a consequence of the narrow band nature of the actual force signal. Also, the effect is frequency dependant: the higher harmonic centre frequency the greater the ‘spread’ of excitation energy. This phenomenon has been observed for walking forces elsewhere [10,18]. Hence, if resonance of a single degree of freedom system representing a mode of interest is assumed, the spread of energy into nearby frequencies results in reduced structural response for higher harmonics compared with the prediction using a ‘perfectly’ periodic model. Parametric studies can also demonstrate that the differences reduce with increase in damping of the structure, but are still significant even for damping as high as 3% found, say, in modern grandstands [19]. A mathematical characterisation of naturally irregular jumping pulses from individuals is hence one of the key problems which needs addressing by modern design guidelines if they are to represent realistically dynamic structural response due to individual people, groups and crowds jumping. If reliable models existed representing the complete Fourier spectrum of continuously measured jumping forces (see Figs. 5a and b), the reconstruction of the force in the time domain would be possible. Although a quite good model of the Fourier amplitudes (Fig. 5a) could be obtained (e.g. as suggested for the walking forces elsewhere [10,18]), variations in the Half-sine model Fourier series model

Force amplitudes/body weight

5 4 3 2 1 0 -1 0.0

0.2

0.4

0.6

0.8

1.0 1.2 Time [s]

1.4

1.6

1.8

2.0

1.4 1.2 1 0.8 0.6 0.4 0.2 0 1

1.5

2 Frequency [Hz]

2.5

3

0.6

Fourier amplitude/body weight

1.6

Fourier amplitude/body weight

Fourier amplitude/body weight

Fig. 3. Half-sine function and sum of the first six Fourier terms to represent the force history given in Fig. 1 (after Ellis and Ji [17]).

0.5 0.4 0.3 0.2 0.1 0

0.10 0.08 0.06 0.04 0.02 0

3

3.5

4 Frequency [Hz]

4.5

5

5

5.5

6

6.5

7

Frequency [Hz]

Fig. 4. Fourier amplitudes of measured jumping force (black) and corresponding periodic model (grey) due to jumping at 2 Hz in the vicinity of the dominant harmonics at a) 2 Hz, b) 4 Hz and c) 6 Hz.

Fourier amplitudes/body weight

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3041

1.5

1

0.5

0 0

2

4

6

8

10

6

8

10

Frequency [Hz]

Phase [rad]

4 2 0 -2 -4 0

2

4 Frequency [Hz]

4000

Force [N]

3000 2000 1000 0 -1000 0

5

10 Time [s]

15

20

Fig. 5. (a) Amplitude and (b) phase FFT of real-life measured jumping data, demonstrating apparently ‘random’ measured phases; (c) regenerated time history based upon Fourier amplitudes in (a) and phase lags based upon the uniform random domain [  p, p].

jumping force between subsequent jumps (Fig. 1) yielding widely varying harmonic phase lags (Fig. 5b) would be very difficult to characterise analytically. If they are, however, assumed to be uniformly distributed in the range [  p, p] (which is currently the only known modelling strategy for this phenomenon [18]), the sum of dominant Fourier series sinusoids normally does not match the real jumping force time history. This can be clearly illustrated by comparison of Figs. 1–5c. Therefore, randomising phases is not the way forward. Their variation seems to be more subtle than a simple set of uniformly distributed random numbers. Bearing all this in mind, a more advanced modelling strategy than Fourier series approach is needed to approximate reliably the narrow band nature of the actual jumping loading.

2.2. Narrow band models First attempts to account for the near-periodic nature of jumping forces can be attributed to studies by Ellis and Ji [20] and Kasperski and Agu [21]. They modelled the jumping frequency and Fourier coefficients of each jumping pulse using probability density functions. However, although the parameters generated by such models are random in nature, they are independent from values of the parameters calculated for the preceding jumps. This seems not to be the case in reality.

3042

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

With regard to this, Sim et al. [8] showed that there is some structure in the slight variation of timing between peaks of subsequent jumps, which can be predicted by the first order autoregression model. This implies that a jumper adjusts the timing of the current jump according to the timing of the previous jump. The variations in the timing were further related to variations in amplitudes of jumping pulses yielding near-periodic GRF time series. However, Sim and co-investigators did not manage to model reliably the full frequency range contained in measured jumping forces. A good compatibility between frequency spectra of measured and modelled force signals was found only for the first two dominant harmonics. The reason for this was the assumption that a very smooth, symmetric, cosine-squared fitting function could represent the irregular shape of measured pulses (Fig. 6). As such, it requires smaller number of Fourier components, thus cannot model accurately high frequency content. For vibration serviceability assessments, the number of Fourier harmonics to be taken into account depends on their contribution to the vibration response and type of structure under consideration. The latest BRE digest 426 [16] suggests that even small energy of the force signal around the sixth dominant harmonic can cause vibration serviceability issues. This might not be relevant for design of grandstands and spectator galleries, but it is relevant for vibration serviceability assessment of, say, multi-story apartment buildings including a fitness centre. The cosine-squared function could not also represent a wide variety of the pulse shapes which can be generated at different jumping rates. As noted by Sim [7], there are three characteristic pulse shapes: double peaked, merging and single peaked, as illustrated in Fig. 6. The double peaked shape is often generated by jumping at rates below 2 Hz, when the landing and launching actions are separated visibly by a long bounce. The majority of people generate merging shapes when jumping around 2 Hz. They are either due to a brief bounce or because the feet hit the ground at slightly different times. Also, a sharp peak may appear at the beginning of the landing phase if the heels are the first to strike the ground, depending on the style of jumping. At rates higher than 2 Hz, the pulse shapes are mostly smooth and single peaked. This happens when people hit the ground with the left and right toes ‘simultaneously’ and propel straight away into the air

Shape 1

Shape 2

2 1.5 1 0.5

Force/Body Weight [-]

2.5

0

2 1.5 1 0.5 0

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.1

2 1.5 1 0.5 0

0.4

0.5

0.6

0.3

1 0.5 0

0.7

0.1

0.2

0.4

2

1 0.5

0.5

0

0.1

0.2

0.3

0.4

0.5

0

1.5 1 0.5

3.5 3 2.5 2 1.5 1 0.5 0

0

0.5

0.6

0.7

0.1

0.2

0.3

0.4

0.5

Time [s]

Force/Body Weight [-]

Force/Body Weight [-]

2

0.4

5 4.5 4 3.5 3 2.5 2 1.5 1 0.5 0

Time [s] 4

2.5

0.3

Time [s]

0 0.2

3 Force/Body Weight [-]

0.3

1.5

Time [s]

2.5 Hz

0.2

Force/Body Weight [-]

3 2.5

Force/Body Weight [-]

Force/Body Weight [-]

3

0.1

2 1.5

Time [s]

2.5

0

3 2.5

0 0

Time [s]

2.0 Hz

Shape 3 3.5

2.5 Force/Body Weight [-]

1.5 Hz

Force/Body Weight [-]

3

4.5 4 3.5 3 2.5 2 1.5 1 0.5 0

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

Time [s]

Time [s]

Time [s]

Fig. 6. Examples of different shapes of measured force traces due to jumping at 1.5, 2 and 2.5 Hz.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3043

2000 I i-1

Ii

I i+1

I i+2

Force [N]

1500 P i-1

Pi

P i+1

P i+2

1000 T c,i

T c,i-1

T c,i+2

T c,i+1

500 T i-1 0 4.4

T P,i-1 4.6

4.8

T i+1

Ti T P,i 5.0

5.2

T P,i+1 5.4 Time [s]

5.6

5.8

T P,i+2 6.0

6.2

6.4

Fig. 7. A portion lasting 2 s of measured force signal due to a single person jumping at 2 Hz. The complete signal lasts 25 s.

having no time for the heel contact. Bearing all this in mind, the quality of Sim’s model depended on jumping rate. At higher rates, the cosine-squared function could approximate the shape of actual pulses better [7,8]. Hence, the best fit was for jumping when people generate landing and launching impulses which are almost the same, yielding fairly symmetrical single peaked shapes. However, because of symmetry and smoothness of the fitting function, double peaked pulses were not considered in the modelling at all. More recently, Racic and Pavic [11,14] made a step forward and used a sum of two Gaussian functions to account for the omnipresent lack of symmetry of the single peaked shapes. Although by varying the overlapping between the Gaussians their sum could fit reasonably well main features of all three pulse shapes (such as the number of peaks and the general form), modelling local irregularities yielding higher Fourier harmonics (e.g. the sharp peak in Shape 3 at 2 Hz in Fig. 6) still remains a point of concern. A way to overcome this problem is to increase the number of Gaussians in the sum, as it will be shown in Section 3.5. Not only can the pulses change shape when a person is jumping at different rates, but successive pulses can also take different shapes for jumping at a single frequency. This typically happens for moderate rates around 2 Hz when the successive pulses switch the shape randomly between merging and single peaked profiles (see Fig. 7). To the best knowledge of the authors, there is no model available in the literature which takes into account changing the shape between successive jumping pulses. Sim developed her stochastic model using the most comprehensive database of human jumping loads available worldwide, which comprises approximately 1000 force histories from about 100 individuals jumping alone on a rigid force plate [22]. Such a large number of the GRF time series provided a statistically reliable platform for the study of the inter-subject variability. However, the poor resolution of the measured jumping rates, which included only four tempos at 1.5, 2.0, 2.67 and 3.5 Hz, resulted in a lack of statistical rigour when studying certain aspects of the intra-subject variations. For example, a study of changes of the force patterns for different jumping rates would be very rudimentary and fragmented if based on such a coarse dataset. This means that there remains a requirement to establish a sufficiently large database of GRF time series which will also have a fine resolution of measured jumping rates to provide statistical reliability in the study of both inter- and intra-subject variability. Establishment of such a database is the key aspect outlined in the rest of this paper, together with the utilisation of the database as a solid foundation for developing a more realistic mathematical model of jumping forces which can be used reliably in vibration serviceability assessment of assembly structures.

3. Concept: from measured to synthesised force The purpose of this section is to describe the concept of the new model development, so that the reader can follow the rationale for the remaining parts of the paper. This will be done by demonstrating, step-by-step, how a single measured force trace can be utilised to generate its synthetic counterpart. Two components are needed: a good quality measured force trace from a single test subject and an appropriate mathematical model of this trace. Therefore, the method presented in this section accounts only for the intra-subject variability. It will be extended later, in Section 5, to account for inter-subject variability as well. This will be done based on utilisation of a database of measured jumping force histories generated by a diverse human population, as presented in Section 4.

3044

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3.1. Step 1: Measured force records Fig. 7 illustrates a portion of a continuously measured force time history generated by a single person jumping for 25 s in response to a regular metronome beat at 2 Hz. The force signal was recorded at a sampling frequency of 200 Hz. By looking at this signal, it is apparent that there is some irregularity on a jump-by-jump basis in terms of variation of shape of jumping pulses, duration and force amplitudes. Therefore, a mathematical model for characterisation of irregular jumping pulses must include modelling of the signal parameters which can represent this natural variability best. Selection of these parameters is the key aspect outlined in the next section. 3.2. Step 2: Basic processing of measured force time history parameters

1950

0.34

1900

0.33

Contact time Tc,i [s]

326 324 322 320 318 316 314 312 310 308 306 0.47 0.48 0.49 0.50 0.51 0.52 0.53 0.54

Peak amplitude Pi [N]

Impulse Ii [Ns]

A jumping cycle is the period of time between any two nominally identical events in the jumping process. In the context of this paper, the instant at which the feet hit the ground (also known as ‘initial contact’) yielding a new pulse was selected as starting (and completing) event (Fig. 7). From the 25 s long force signal illustrated in Fig. 7 yielding about 50 jumping pulses, a window comprising 42 successive jumping cycles was selected for further analysis. A total of eight cycles were discarded from the start and end of this time history. The force threshold marking the start of the pulse was set to 15 N. From each of the 42 cycles, contact time Tc,i, period Ti (i.e. duration of a jumping cycle), peak timing TP,i, peak amplitude Pi and impulse Ii are extracted on a cycle-by-cycle basis (Fig. 7). In the past, statistical models of these parameters and their mutual relationships were used to describe intra-subject variability of jumping forces. For example, Sim et al. [8] derived their model based on the theoretically derived linear relationship between the impulse size Ii and timing T P,i ¼ ðTP,i þ 1 TP,i1 Þ=2. However, using the measured values extracted from the 42 cycles, Fig. 8a illustrates a lack of linear correlation between these two parameters. This is probably because the model proposed by Sim et al. [8] was based on an assumption that the impulse is of very short duration i.e. instantaneous, whereas in reality it lasts as long as the

1850 1800 1750 1700 1650 0.48

0.49

0.50

0.51

0.52

0.32 0.31 0.30 0.29 0.28 0.48

0.49

0.50

0.51

0.52

Period of jumping cycle Ti [s]

Period of jumping cycle Ti [s]

Timing TP,i [s]

Fig. 8. Examples of poor correlation between (a) timing T P,i and impulse size Ii (the correlation coefficient Z = 0.104), (b) periods Ti and peak amplitudes Pi (Z = 0.235)and (c) periods Ti and contact times Tc,i (Z = 0.418).

Weight normalised impulses Iw, i [s]

0.52 measured data linear fit 0.51

0.50

0.49

0.48 0.48

0.49

0.50

0.51

0.52

Period of jumping cycle Ti [s] Fig. 9. Correlation between periods Ti and weight normalised impulses Iw,i (Z =0.801).

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3045

contact period. Also, Sim et al. [8] assumed that a launching impulse which projects the jumper into the air and the subsequent landing impulse are the same. However, this is only true for the frictionless impact between two infinitely stiff i.e. rigid bodies. In reality, it is possible for the stiffness of the legs to be adjusted between successive jumps. For example, Farley et al. [23] showed that when humans hop in place, the neuromuscular system can significantly alter the stiffness of the leg to accommodate variations in hopping height at a given frequency. Differences in launching and landing pulses might also be the most likely explanation for the asymmetry of jumping pulses. More examples of the lack of correlations between other jumping force parameters are given in Figs. 8b and c. However, certain correlations between the jumping parameters seemingly do exist. In the study presented here, it was found that there is some relationship between period Ti and weight normalised impulses Iw,i calculated over integration time Tc,i as Iw,i ¼

n 1 X F Dt, mg i ¼ 1 i

Dt ¼

1 fs

ð1Þ

where m is the body mass, g is gravity (g =9.81 m/s2), Fi is the force, fs is the sampling rate and n is the total number of samples of Tc,i. This relationship shows a linear trend (Fig. 9), hence provides an opportunity to describe weight normalised impulses as a function of the period. This is an important feature for the proposed methodology, as will be shown later in Section 3.4. 3.3. Step 3: Analysis of periods of jumping pulses The problem now is to model slight variations of period Ti (i= 1, y, 42) between each jump. This effect can be represented by a sequence of numbers T i calculated as Ti ¼

Ti mT

mT mT ¼ meanðTi Þ

ð2Þ

Two methods were tried for modelling T i data. It was first assumed that the current T i value is a linear combination of previous k values T i1 ,. . .,T ik , which could be described by an autoregressive model [24]. This method proved unreliable since weak correlation was found between the previous and subsequent T i values for increasing order of regression from one to four (k= 1, y, 4). Therefore, the rest of this section will be focused on the alternative and more successful approach based on utilisation of the auto spectral density (ASD) of T i . Given the single-sided spectral density Sx(f) of the real random process x(t), the variance of x(t), s2x , can be computed by using the relation [25]: Z 1 s2x ¼ Sx ðf Þ df ð3Þ 0

0

The point here is to use the ASD (as illustrated in Fig. 10) to artificially generate synthetic T i series having the same standard deviation as the actual T i series. When doing this, the newly generated set of numbers is assumed to have the

x10-3 3 measured Gaussian fit

2.5

ASD [1/Hz]

2 1.5 1 0.5 0

c5 0

0.1

c10 0.2 0.3 Quasi-frequency

c20

c15 0.4

Fig. 10. Single-sided ASD of discrete data T i .

0.5

3046

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

same ASD as the measured set of numbers. This applies regardless of the number of data points in the old and new sets, which will be demonstrated later in this section. This assumption also means that the standard variation of T i does not change for the given jumper, jumping rate and jumping interval. Moreover, the ASD preserves a ‘frequency structure’ between T i values (if there is any), which the autoregressive models could not represent. The ASD of T i can be calculated as A2 ðfm Þ

m T , m ¼ 0,. . .,20 ð4Þ , fm ¼ 42 2Df where AT ðfm Þ is a single-sided discrete Fourier amplitude spectrum having spectral line spacing Df= 1/42. The ASD ordinates do not depend on the number of discrete data points T i used for the calculation of AT ðfm Þ but it is coarse due to the limited number of points (21 single-sided FFT). More points might reveal a much richer structure but this requires a longer jumping force record. However, continuous jumping demands significant effort of a test subject, so the duration is limited to avoid causing fatigue which can influence the force records. This will be discussed further in Section 4.2.The ASD ST(fm) can be analytically described by a series of Gaussian functions (Fig. 10): ST ðfm Þ ¼

S0T ðf Þ ¼

21 X

2

Wj eððf cj Þ

2

=2dj Þ

ð5Þ

j¼1

Here, parameter Wj is the height of the jth Gaussian peak, cj is the position of the centre of the peak, and dj controls the width (i.e. time duration) of the corresponding bell-shaped curve. The Gaussian centres cj, j =1, y, 21, are placed in each sample on the quasi-frequency axis in order to fit exactly the measured ASD (Fig. 10). For such fixed positions of Gaussian centres cj and predefined widths dj = Df, Gaussian heights Wj (also called weights) can be computed using the non-linear least-square method [26]. The results of the fitting are given in Table 1. Representation of the discrete ASD ST(fm) in the form of the continuous function S0 T(f) enables calculation of the ASD ordinates for an arbitrary spectral line spacing Df. As the ASD ordinates do not depend on the number of discrete data 0

points used for its calculation via FFT, the continuous function S0 T(f) can be used to generate a set of data points T k (k= 1, y, N) which would have the same ASD as the measured data. This is done by calculating the new discrete spectral Table 1 Parameters of function S0 T. j

cj (dimensionless)

Wj (dimensionless)

j

cj (dimensionless)

Wj (dimensionless)

1 2 3 4 5 6 7 8 9 10 11

0.0000 0.0238 0.0476 0.0714 0.0952 0.1190 0.1429 0.1667 0.1905 0.2143 0.2381

 0.0040 0.0027  0.0001 0.0008 0.0020 0.0008  0.0001 0.0009 0.0002 0.0004  0.0001

12 13 14 15 16 17 18 19 20 21

0.2619 0.2857 0.3095 0.3333 0.3571 0.3810 0.4048 0.4286 0.4524 0.4762

0.0003  0.0001 0.0002 0.0001 0.0000 0.0001 0.0001 0.0001  0.0001 0.0000

Variation of jumping periods Ti [-]

0.03 measured simulated

0.02 0.01 0 -0.01 -0.02 -0.03 1

5

10

15

20

25

30

35

40 42

Jumping cycle [-] Fig. 11. Variations of peak-to-peak intervals on a cycle-by-cycle basis.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3047

line spacing Df 0 ¼ 1=N. For a sequence of discretely spaced frequency points fn =n Df0 (where n ¼ 0,. . .,N=21) corresponding to the new set of data points to be generated, a new set of ASD amplitudes is calculated using S0 ðfn Þ. T qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi This is then used in Eq. (4) to generate a new set of Fourier amplitudes A0 ðfn Þ ¼ 2Df 0 S0 ðfn Þ. Finally, assuming a set of T

T

0

randomly distributed phases in the range [ p, p], a new set of N variations T k can be generated by inverse Fourier transform making use of A0 ðfn Þ and the phases. Different realisations of the random phases may be specified by varying the seed of the T

0

random number generator, hence many different series T k may be generated with the same spectral properties (Fig. 11). 0 Because of the relation defined by Eq. (3), all series T k generated also have the same variance s2T regardless of their 0 length N. Bearing all this in mind, it is possible to generate T k series of arbitrary length (e.g. N 542), which will have statistically the same properties of variations on the sample-by-sample basis as the measured set of 42 actual T i data 0 points. According to Eq. (2), scaling T k by mT and adding the offset value mT calculated from the measured data, results in a 0 series of synthesised jumping periods T k , as would be generated by the test subject during nominally identical jumping exercises. Empirical evidence for this is presented in Section 3.6. 3.4. Step 4: Analysis of impulses of jumping pulses As previously demonstrated (Fig. 9), the relationship between the weight normalised impulses Iw,i and the durations of the jumping cycles Ti is approximately a linear function. Therefore, this relationship can be used to generate a new set of synthetic impulses I0 w,k based on the generated set of synthetic periods T0 k from the previous section. The measured weight normalised impulse values Iw,i can be expressed as a function of period Ti (Fig. 9) using the following linear regression model [24]: Iw,i ¼ r1 Ti þ r0 þ DIw,i ð6Þ Here, r1 = 0.832 and r0 = 0.082 are regression coefficients and DIw,i is the subsequent error (also known as a disturbance term), which is a random variable. Given the regression coefficients r1 and r0, DIw,i can be calculated from the real data as

DIw,i ¼ Iw,i r1 Ti r0

ð7Þ

It is common to model DIw,i as Gaussian noise [24] having a probability density function jDI(x) given by [27]

jDI ðxÞ ¼

2 1 2 pffiffiffiffiffiffi eððxmI Þ =2sI Þ ,

sI 2p

x2R

ð8Þ

where x is a synthetic random variable which corresponds to measured DIw,i, mI = 0 is the mean value and sI is the standard deviation of real DIw,i. Given a set of artificially generated periods of jumping pulses T0 k, as explained in Section 3.3, a series of corresponding synthetic weight normalised impulses I0 w,k can be therefore calculated using Eqs. (6) and (8). This will be demonstrated in Section 3.6. 3.5. Step 5: Analysis of shape of jumping pulses By looking at the force signal given in Fig. 7, each jumping pulse can be observed as a function of time having distinctive size and shape. Therefore, the aim is to extract as many as possible different jumping pulses from a continuously measured jumping force time history and mathematically describe them. The models will be then used to artificially generate a force signal which includes an arbitrary number of jumping pulses. This assumes that the extracted pulses can represent the individual’s long term performance for the nominally identical jumping exercises. A question arises here about the minimum number of successive jumping pulses needed to ensure reliable representation. At given jumping rate, the total number of measured pulses depends on the duration of the force record. Jumping is a demanding activity requiring a lot of effort, so the duration should be just long enough to avoid causing fatigue, physical distress or ethical issues. This will be discussed further in Section 4.2. Published research on the subject is very rare and limited. Parkhouse and Ewins [22] claimed that a minimum of 30 consecutive jumping cycles should constitute the shortest duration of the force signal to calculate reliably the corresponding Fourier coefficients, provided a person does not tire visibly or obviously change their jumping pattern. However, this suggestion was made without giving much justification, hence it requires verification. Using a statistical technique called sequential estimation analysis, Rodano and Squadrone [28] found that at least 12 single and nonconsecutive jumping pulses (i.e. with pauses between jumps) were needed to obtain a stable mean of joint kinetics, such as hip, knee and ankle internal forces and moments, derived from vertical jumping forces. Assuming that at least 12 jumping cycles are also necessary to reach stability of the mean of continuously measured jumping pulses, the set of 42 successive, weight normalised pulses are used further in the analysis of the corresponding pulse shapes. Each weight normalised pulse was individually extracted into a half-second segment and scaled by the corresponding peak amplitude ai yielding unity normalised pulses (Fig. 12), so that differences in pulse shape could be investigated independently of amplitude. These scaling factors can be expressed as the ratio:

ai ¼

Iw,i , Is,i

i ¼ 1,. . .,42

ð9Þ

3048

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

1.2 Gaussian functions Gaussian fit Gaussian centres measured signal

Force/body weight [-]

1 0.8 0.6 0.4 0.2 0

0

0.1

0.2

0.3

0.4

0.5

Time [s] Fig. 12. Example of unity normalised jumping pulses extracted from Fig. 7.

Gaussian functions Gaussian fit Gaussian centres

1 0.8

z

0.6 0.4 0.2 0

r -1

-0.5

0

0.5

1

θ(t)

1 0.5

-1

-0.5

0

Fig. 13. Trajectory Zi(t) in a three-dimensional (3D) space. Figs. 13 and 12 represent the same data.

where Is,i are impulses of unity normalised pulses. The ratio ai will be used later in Section 3.6 to describe the smooth modulation of measured amplitudes of subsequent pulses in a synthetic jumping force signal. For fs = 200 Hz sampling rate, each unity normalised pulse consists of 100 samples (Fig. 12). Therefore, the underlying pattern of the ith pulse can be modelled mathematically as a sum of 100 Gaussian functions Zi(t): Zi ðtÞ ¼

100 X

2

Air eððttir Þ

=2b2ir Þ

,

t 2 ½0,0:5,

i ¼ 1,. . .,42

ð10Þ

r¼1

where the parameter Air is the height of the rth Gaussian peak, tir is the position of the centre of the peak, and bir controls the width of the corresponding Gaussians (Fig. 12). The Gaussian centres tir ¼ tr ¼ n Dt, n =0, y, 99, Dt = 0.005 s are placed in each sample on the time axis to fit exactly the measured pulse amplitudes, thus to reflect completely the corresponding Fourier amplitude spectrum. A less dense distribution of Gaussian exponentials (e.g. centres are placed in every second or third sample) will make the fit smoother causing the high frequency components to vanish. For such fixed positions of Gaussian centres tir and predefined widths bir ¼ br ¼ Dt, Gaussian heights Air can be optimised using non-linear least-square curve fit [26]. Analytical functions Zi(t) can be used to replicate unity scaled jumping pulses. Furthermore, they can be individually scaled in terms of amplitude (vertically) and time (horizontally) to model realistic impulses and periods of jumping cycles, respectively. This will be demonstrated in the next section. 3.6. Step 6: Dynamic model Periodic behaviour of jumping could be better visualised when individual pulses are ‘wrapped’ around the surface of a cylinder (Fig. 13). Now the first and the last sample of each pulse overlap, thus they become closed orbits that encircle the cylinder.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3049

The 42 consecutive jumping pulses can be generated by retracing the Gaussian fits Zi(t), given by Eq. (10), around a circle of unit radius in the (r, y) plane (Fig. 13). The time required to complete one revolution on this circle is equal to the period of corresponding jumping cycle and lasts Ti seconds. The angular frequency of rotation oi can be calculated as

oi ¼

 2p  rad=s Ti

ð11Þ

The equations of motion of a point moving around the circle are therefore given by a set of two equations [14]: rðtÞ ¼ 1 yðtÞ ¼ oi t,

y 2 ½0,2p

ð12Þ

Further, Eq. (10) can be rewritten as a function of the angular coordinate y instead of time t: Zi ðyÞ ¼

100 X

2

ðyy2r Þ

Air e

2br

,

y 2 ½0,2p, i ¼ 1,. . .,42

ð13Þ

r¼1

Unity-normalised pulses [-]

where yr ¼ oi tr and br ¼ oi br are in radians. The time positions of the Gaussian peaks tr now correspond to fixed angles yr along the unit circle, as illustrated in Fig. 13. For pairs ðoi ,Zi ðyÞÞ sorted in numerical order i=1, y, 42, coupled system of Eqs. (11)–(13) through 42 iterations generates a synthetic signal which is identical to the real signal comprising the 42 consecutive unity normalised jumping pulses. However, the aim is to artificially generate a force record of arbitrary duration, i.e. which includes an arbitrary number N of jumping pulses (e.g. N 542), the test subject could generate during nominally identical jumping exercises. Let T0 k (k =1, y, N) be a series of periods of jumping cycles computed as explained in Section 3.3. The corresponding angular frequency o0 k of circular motion is then given by Eq. (11). It can be further assumed that the duration of the jumping cycle does not influence the general shape of unity normalised pulses at a given jumping rate. Under this assumption, any of Zi(y) can be assigned randomly and equally likely to each o0 k yielding pairs ðo0k ,Zk ðyÞÞ. Even when the

1.2 1 0.8 0.6 0.4 0.2 0 0

5

10

15

20

25

15

20

25

Unity-normalised pulses [-]

Time [s]

1.2 1 0.8 0.6 0.4 0.2 0 0

5

10

Unity-normalised pulses [-]

Time [s]

1.2 1 0.8 0.6 0.4 0.2 0 0

5

10

15

20

25 Time [s]

30

35

40

45

50

Fig. 14. (a) Measured, and examples of synthetic unity-normalised signals when (b) N =50 and (c) N = 100. Jumping rate is 2 Hz.

3050

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

same Zi(y) is assigned to different o0 k, the resulting jumping pulses generated by the coupled system of Eqs. (11)–(13) will be slightly different. This is because higher angular frequencies o0 k generate the pulses faster, hence they compress them resulting in shorter intervals of jumping cycles and vice versa. Fig. 14 illustrates examples of 25 and 50 s long synthetic unity normalised signals generated by Eqs. (11)–(13) when the total number of pairs ðo0k ,Zk ðyÞÞ is N =50 and 100, respectively. To reflect the changes in both the amplitude and timing of the real jumping force from one pulse to another, amplitudes of each synthetic unity normalised pulse of the kind shown in Fig. 14 need to be multiplied by the corresponding element in series a0 k defined by Eq. (9). Unity normalised impulses Is,k in the denominator of ratio a0 k are fixed values for each unity normalised pulse Zk(y) and can be calculated as a definite integral of the corresponding function Zk(t) between 0 and 0.5 s: Is,k ¼

Z 0

0:5

Zk ðtÞ dt 

99 X

Zk ðn DtÞ Dt

ð14Þ

n¼0

As demonstrated in Section 3.4, for a generated set of synthetic T0 k, the corresponding set of the body weight normalised impulses I0 w,k in the nominator of ratio a0 k can be found from Eq. (6). After being scaled by series a0 k, signals given in Fig. 14 are still dimensionless and feature variations of the jumping force parameters on a jump-by-jump basis (shape of the pulses, periods, impulses and peak amplitudes), as the test subject could generate in reality during nominally identical jumping exercises. These signals become equivalent jumping force time histories when their amplitudes are additionally multiplied by the body weight of the test subject (Fig. 15). The similarity between the real measured and synthetic near-periodic vertical jumping force signals may be seen by comparison of Fig. 15a–c. This comparison looks much better than the one given in Section 2.1 between Figs. 1 and 5c featuring standard Fourier transform approach with randomly generated phases used to recreate the forcing function. The standard Fourier amplitude spectra are compared in Fig. 16. For the first four dominant harmonics the relative errors are within the range of 73%. Moreover, relative error in the area under the graph of the spectra (i.e. overall energy of the

Force [N]

2000 1500 1000 500 0 0

5

10

15

20

25

15

20

25

Time [s]

Force [N]

2000 1500 1000 500 0 0

5

10

Time [s]

Force [N]

2000 1500 1000 500 0 0

5

10

15

20

25 Time [s]

30

35

40

45

50

Fig. 15. (a) Measured and examples of synthetic force signals when (b) N =50 and (c) N = 100. Jumping rate is 2 Hz.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3051

Fourier amplitudes/ body weight [-]

1.5 1 0.5 0 0

2

4

6

0

6

0

6

0

Frequency [Hz]

Fourier amplitudes/ body weight [-]

1.5 1 0.5 0 0

2

4 Frequency [Hz]

Fourier amplitudes/ body weight [-]

1.5 1 0.5 0 0

2

4 Frequency [Hz]

Fig. 16. Fourier amplitude spectra of (a) measured and synthetic force signals when (b) N =50 and (c) N =100. Jumping rate is 2 Hz.

signals) is less than 5%. This indicates a good match in the frequency content between the measured and synthesised GRF signals. Perfectly identical signals in both time and frequency domains can be generated only by chance. However, the following properties are identical between measured and synthesised forces:

1) Shapes of the jumping pulses are drawn from the same source (i.e. functions Zi(t)), where each shape has the same probability of occurrence. 0 2) Because of the common ASD, the quality and quantity of variations of jumping periods T k are the same for all synthetic 0 signals. Quality means that the relationship between successive T k data points follow a nominally identical pattern for 0 all signals, whereas the quantity means that the standard deviation of T k is a fixed value. 0 3) The statistical equivalence between T k values reflects directly equivalence between the corresponding weightnormalised impulses I0 w,k according to Eq. (6). This implies that total energies of generated signals are the same.

The modelling procedure presented in this section has been applied only to a single force record. The step forward is to apply the same procedure to a sufficiently large database of force records, to extract the modelling parameters from each record and to use them to develop a stochastic model of individual jumping loads. Establishment of such a database is the key aspect outlined in the next section.

3052

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

4. Experimental data acquisition A database comprising many high-quality jumping force records is an essential component for the development of a stochastic model of jumping loads. Section 4.1 describes only facts related to the experimental setup used, so that it is clear how the data was collected. Section 4.2 explains test protocol including all details which need considering to gather a quality database. 4.1. Experimental setup Testing was carried out in the Light Structures Laboratory in the University of Sheffield. The GRFs of one person jumping at a time were recorded by a single AMTI BP-400600 force plate [29] rigidly fixed to the laboratory floor, as illustrated in Fig. 17. All forces were sampled at 200 Hz. As well as safety, the platform built around the force plate (Fig. 17) gave the impression of a bigger jumping space to avoid test subjects deliberately targeting a relatively small plate surface area (0.6 m  0.4 m), which might influence the natural variability of the GRFs. 4.2. Acquisition of quality experimental data Fifty-five volunteers were drawn from students, academics and technical staff of the University of Sheffield. Participants were adult people from diverse ethnic groups, different genders (38 males and 17 females), body size and shape (body mass 73.2720.6 kg, height 1.7270.12 m) and varying age categories (33.176.6 years). The stochastic approach to modelling jumping loads, presented in Section 5, is based on an assumption that these fifty-five persons can represent general human population. Each participant was asked to perform 15 jumping tests. Each test followed the same pattern: 1) 2) 3) 4) 5)

the participant stood on the force plate, and then was given a constant metronome beat, was allowed a brief practice for a few seconds prior to the data collection, was asked to jump for 30 s following the metronome beat (middle 25 s were recorded), was asked to leave the force plate and rest.

The constant metronome beat was chosen in a quasi-random order from 15 different rates in the range 1.4–2.8 Hz having fine resolution of 0.1 Hz. The range included slow and fast jumping frequencies which were suggested in the past as being comfortable for individuals [12] and at which synchronised jumping of groups and crowds can occur [30]. Therefore, each

force plate

safety platform

Fig. 17. Experimental setup.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3053

participant generated 15 force records yielding a database of 825 force time histories in total for all participants. Typical measured force records are given in Fig. 18. The closely spaced jumping rates gave the so established database a definite advantage over the other similar datasets collected worldwide so far, such as those published by Parkhouse and Ewins [22,31]. This brings more statistical reliability into the modelling intra- and inter-subject variations of jumping force patterns due to changes in the jumping rate, which models published so far could not represent to such an extent [7,8,20,21]. The longer a jumping test lasts, the better insight it provides into the natural variability of the measured jumping force history. However, there is a limit on ability of people to keep jumping for a long time without tiring or changing their motions. As mentioned in Section 3.5, jumping is a physically intensive activity, so ethical concerns limited duration of the tests to 30 s. In feedback from the participants, this duration was commonly considered optimal. The participants were not given any explicit instructions about their jumping technique, but they were encouraged to move as if they were enjoying a lively concert or an aerobic exercise. Prior to the experiment, the test protocol (approved by the Research Ethics Committee of the University of Sheffield) required that each participant should complete a Physical Activity Readiness Questionnaire and pass a preliminary fitness test (by satisfying predefined criteria for blood pressure and resting heart rate) to check whether they were suited to the kind of physical activity required during the measurements. All participants wore comfortable sportswear and trainers.

Force [N]

1500

1000

500

0 0

5

10

15

20

25

15

20

25

15

20

25

Time [s]

Force [N]

1500

1000

500

0 0

5

10 Time [s]

4000

Force [N]

3000 2000 1000 0 0

5

10 Time [s]

Fig. 18. Examples of measured force signals generated by three different persons jumping in response to a regular metronome beat at (a) 1.5 Hz, (b) 2 Hz and (c) 2.5 Hz.

3054

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

5. Development of stochastic jumping model After a few preparatory steps outlined in Sections 5.1 and 5.2 combines the modelling strategy developed in Section 3 and the database of measured GRF signals from Section 4 to create a stochastic model to generate synthetic human jumping loading. 5.1. Utilisation of existing database

3.6 3.4 3.2 3 2.8 2.6 2.4 2.2 2 1.8 1.6

4.5

Average α i [-]

4

Average α i [-]

Average α i [-]

The 825 force signals measured in Section 4 were classified into 15 categories (clusters) with respect to the actual jumping rate (1.4–2.8 Hz with a step of 0.1 Hz). For example, all force records with the actual rate in the range 1.950–2.049 Hz are gathered into a cluster at 2 Hz. In fact, giving the participants constant metronome beats does not explicitly mean that the beat frequency was followed. This is because not all individuals are able to synchronise their movements to the beats [21]. If they were able to do so, the clusters would comprise 55 force records each, which is the

3.5 3 2.5 2

40

50

60

70

80

90

100 110

40

Body mass [kg]

50

60

70

80

90

4.2 4 3.8 3.6 3.4 3.2 3 2.8 2.6 2.4 2.2

100 110

40

50

Body mass [kg]

60

70

80

90

100 110

Body mass [kg]

Fig. 19. Average ai vs body mass for the clusters (a) 1.5 Hz, (b) 2 Hz and (c) 2.5 Hz representing slow, moderate and fast jumping rates, respectively.

input data: jumping rate and duration of jumping number of jumps N

from the corresponding cluster get a set of the force parameters by chance random number seed

pulse shapes

ASD of variations of jumping periods S (f)

phases [-π,π]

random number seed generation of variations of jumping periods T unity normalised impulses I

assign one pulse shape to each jump generation of jumping periods T and the corresponding angular frequencies ω

coupled system of equations pdf of body mass

unity normalised GRFs

random number seed

dynamic impact factors α

body weight

synthetic GRF signal Fig. 20. Algorithm describing the procedure for generating synthetic GRF signals.

weight normalised impulses I

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3055

total number of participants. In the present study, clusters for the rates between 1.7 and 2.5 Hz include more than 55 measurements (maximum 67 for 2 Hz), whereas those in the range 1.4–1.6 Hz and 2.6–2.8 Hz comprise less than 55 measurements (minimum 42 for 2.8 Hz). This might also be an indicator that the rates in the range 1.7–2.5 Hz were the most comfortable jumping frequencies for the majority of participants in this study as they ‘attracted’ jumping rates from other clusters. Each force history within a cluster was processed using concepts described in Section 3. Multiple sets of information (Gaussian weights Wj and Air, Gaussian centres cj and yr, Gaussian widths dj and br, autoregression coefficients r1 and r0, and standard deviation sI) were developed and stored within the cluster as Matlab structural files [32]. Finally, it needs to be shown that the force is independent from the body weight (mass), so that weight normalised synthetic signals generated by the dynamic model (Section 3.6) can be scaled by the weight of any person drawn by chance from the world’s population. By doing this, randomisation of the modelling parameters will be extended to the maximum. The evidence to support this hypothesis is given in Fig. 19, where the average ratio ai for each GRF signal in clusters 1.5, 2 and 2.5 Hz was plotted against the corresponding body mass. The scattered patterns for all three jumping rates indicated no correlation between the two parameters, thus they can be treated as independent variables in the modelling process (Section 3.6). As a further random parameter, body mass can be modelled using a probability density function, as illustrated by Hermanussen et al. [33] for German, Austrian and Norwegian citizens.

5.2. Procedure for generating synthetic forces The flow chart in Fig. 20 illustrates the complete process (algorithm) of creating synthetic GRF signals. For a specified jumping rate and jumping period, the algorithm first estimates the total number N of jumping cycles included in the synthetic force signal. From the corresponding frequency cluster, a set of multiple modelling parameters (such as continuous ASD S0 T(f), regression coefficients r1 and r0 and shapes of jumping pulses) is selected randomly. At this point, the algorithm splits into two parallel actions: creating the shapes and durations of the N jumping cycles. Shapes of jumping pulses Zi(t) are randomly assigned to each jumping cycle, so they become Zi(t), k= 1, y, N. The corresponding unity normalised impulses Is,k will be used later as denominators of ratios a0 k (Section 3.5). Having the spectral characteristics of variations of jumping periods S0 T(f)and the total number of cycles N, synthetic periods T0 k (Section 3.3) and the corresponding angular frequencies o0 k (Section 3.6) can be calculated on the cycle-by-cycle basis. The set of synthetic periods T0 k are also used to calculate weight normalised impulses I0 w,k using autoregressive model given by Eq. (9). These will later be numerators of ratios a0 k. The next step integrates everything generated so far to run the dynamic model (Section 3.6) and therefore to generate the unity normalised jumping pulses of a kind shown in Fig. 14. Finally, these will be transformed into a synthetic force

2500 Force [N]

2000 1500 1000 500 0 0

5

10

15 Time [s]

20

25

30

0

5

10

15 Time [s]

20

25

30

2500 Force [N]

2000 1500 1000 500 0

Fig. 21. Examples of two synthetic force time histories generated by the model for jumping at 2 Hz. The signals correspond to two different individual jumpers.

3056

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

1.6

Fourier amplitudes/ body weight [-]

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

2

4

6

0

6

0

Frequency [Hz]

Fourier amplitudes/ body weight [-]

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

2

4 Frequency [Hz]

Fourier amplitude spectra / body weight [-]

Fig. 22. Fourier amplitude spectra of the synthetic GRFs given in Fig. 21.

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

6

7

8

9

10

6

7

8

9

10

Fourier amplitude spectra/ body weight [-]

Frequency [Hz]

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

Frequency [Hz] Fig. 23. Cluster at 1.5 Hz—Fourier amplitude spectra of all (a) measured forces and (b) their synthetic counterparts.

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3057

Fourier amplitude spectra/ body weight [-]

0 time history after they have been scaled by ratios a0k ¼ Iw,k =Is,k and body weight. As a random parameter, the body weight can be generated using statistical models available elsewhere [33]. Fig. 21 shows examples of the signals generated when the model is run twice in a row for the jumping rate 2 Hz lasting 30 s. Each of the signals corresponds to a unique individual jumper due to the inherent randomness of the modelling parameters. A visual comparison of the two signals provides convincing evidence that the model can account for the intersubject variability in the pulse shapes and force amplitudes. On the other hand, the ability of the model to generate

1.8 1.6 1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

6

7

8

9

10

6

7

8

9

10

Fourier amplitude spectra/ body weight [-]

Frequency [Hz]

1.8 1.6 1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

Frequency [Hz]

Fourier amplitude spectra/ body weight [-]

Fig. 24. Cluster at 2 Hz—Fourier amplitude spectra of all (a) measured forces and (b) their synthetic counterparts.

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

6

7

8

9

10

6

7

8

9

10

Fourier amplitude spectra/ body weight [-]

Frequency [Hz]

1.4 1.2 1 0.8 0.6 0.4 0.2 0 0

1

2

3

4

5

Frequency [Hz] Fig. 25. Cluster at 2.5 Hz—Fourier amplitude spectra of all (a) measured forces and (b) their synthetic counterparts.

3058

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

different degrees of variability in the periods of successive jumping cycles for different persons becomes more obvious from comparison between the corresponding frequency spectra (Fig. 22). The broader spread of energy around dominant harmonics (i.e. integer multiples of 2 Hz) in Fig. 22b relative to Fig. 22a indicates that the first ‘virtual’ jumper (Fig. 22a) varied the periods less than the second one (Fig. 22b). This was also observed in the actual test data as demonstrated in the next section.

6. Model verification The modelling strategy proposed in this paper was validated for each of the 15 frequency clusters. This was done by comparing between standard Fourier amplitude spectra of all measured forces in a cluster and their synthetic counterparts. Figs. 23–25 illustrate results for three clusters at 1.5, 2 and 2.5 Hz, representing slow, moderate and fast tempos, respectively. For the first four dominant harmonics the relative errors between the average measured and average synthetic spectra for each cluster are within the range of 73%. Moreover, the relative error in the area under the graph of the average spectra is less than 7%. All this indicates good match in the frequency content between the measured and synthesised GRF signals. Therefore, synthetic forces generated by the model can be utilised in vibration serviceability assessment of civil engineering assembly structures, such as grandstands, spectator galleries, footbridges and concert or gym floors, to estimate realistically vibration response due to people jumping.

7. Conclusions This paper presents a new mathematical model used to generate near-periodic synthetic jumping force signals with specified jumping rate and morphology. Similar to modelling the near-periodic human heart beats, the near-periodic nature of the jumping force is modelled using closed-loop trajectories throughout 3D space (r, y, z) around a circle of unit radius in (r, y) plane. Each revolution on this circle corresponds to the period of one jumping cycle. The trajectory replicates the size and shape of the measured jumping pulses via a sum of Gaussian exponentials. This modelling strategy can represent temporal and spectral features of the real human vertical jumping loading more effectively than the conventional half-sine pulses and Fourier series approach yielding more reliable predictions of dynamic structural response due to people jumping. The proposed Gaussian fit, coupled with equations of circular motion, has the following considerable advantages: 1) A set of Gaussian bell functions in which centres are placed in each sample of measured pulses can fit exactly any pulse shape. This includes a lack of symmetry, double peak patterns and local irregularities yielding high frequency components, as opposed to the symmetric and smooth half-sine and cosine-squared pulses which can reflect only low frequency content in the corresponding Fourier amplitude spectra. 2) Variations of the jump-by-jump intervals can be included by varying the angular frequency for consecutive revolutions around the unit circle. For each revolution, it is also possible to change the pulse shape more effectively than using the conventional Fourier series approach. As a result, the amplitude Fourier spectrum of corresponding synthetic jumping signal becomes a narrow band random phenomenon showing the leaking of energy in the vicinity of dominant Fourier harmonics. 3) Impulses and amplitudes of the synthetic jumping force signal can be changed on a jump-by-jump basis in a manner which allows the model to simulate ‘smooth’ energy transfer between consecutive jumps, as it is measured in reality. Numerous jumping force records generated by different individuals under a range of jumping frequencies resulted in a comprehensive database of jumping forces for the general human population. The database was used to develop and calibrate a new generation of stochastic models of jumping loading for individuals. This framework can be extended further to stochastic jumping loads due to groups and crowds. At the moment, individual forces can be summed with random phase lags as suggested elsewhere [34–36]. However, there are indications that this is not what is happening in reality and more research into synchronisation between people jumping is needed. This presents an opportunity to enhance the vibration serviceability assessment of civil engineering structures occupied and dynamically excited by humans such as grandstands, footbridges, floors and staircases.

Acknowledgements The authors would like to acknowledge the financial support provided by the UK Engineering and Physical Sciences Research Council (EPSRC) for grant reference EP/E018734/1 (‘Human Walking and Running Forces: Novel Experimental Characterisation and Application in Civil Engineering Dynamics’) and GR/T03017/01 (‘Stochastic Approach to Human–Structure Dynamic Interaction’).

V. Racic, A. Pavic / Mechanical Systems and Signal Processing 24 (2010) 3037–3059

3059

References [1] SCOSS. (2001) Structural Safety 2000–02: Thirteenth Report of the Standing Committee on Structural Safety. The Standing Committee on Structural Safety, The Institution of Structural Engineers, London, UK. [2] BS 6399-1: 1996 Incorporating Amendment No. 1. Code of Practice for Dead and Imposed Loads. British Standards Institution, London, UK, 2002. [3] Canadian Commission on Building and Fire Codes, in: User’s Guide—NBC 2005: Structural Commentaries (Part 4 of Division B), National Research Council Canada, Institute for Research in Construction, Ottawa, Canada, 2006. [4] S. Zivanovic, A. Pavic, P. Reynolds, Vibration serviceability of footbridges under human-induced excitation: a literature review, Journal of Sound and Vibration 279 (2005) 1–74. [5] D. Rogers, Two more ‘wobbly’ stands. Construction News, August 17, 2000. [6] D. Parker, Rock fans uncover town hall floor faults. New Civil Engineer 20 November, 2003. [7] J.H.H. Sim, Human–structure interaction in cantilever stands, Ph.D. Thesis, Faculty of Engineering, The University of Oxford, UK, 2006. [8] J.H.H. Sim, A. Blakeborough, M. Williams, Statistical model of crowd jumping loads, ASCE Journal of Structural Engineering 134 (12) (2008) 1852–1861. [9] H. Bachmann, A. Petlove, H. Rainer, Dynamic forces from rhythmical human body motions, Vibration Problems in Structures: Practical Guidelines, Appendix G, Birkhauser, Basel, Switzerland, 1995. [10] J.M.W. Brownjohn, A. Pavic, P. Omenzetter, A spectral density approach for modelling continuous vertical forces on pedestrian structures due to walking, Canadian Journal of Civil Engineering 31 (2004) 65–77. [11] V. Racic, A. Pavic, Mathematical model to generate asymmetric pulses due to human jumping, ASCE Journal of Engineering Mechanics 135 (10) (2009) 1206–1211. [12] D. Ginty, J.M. Dervent, T. Ji, The frequency ranges of dance-type loads, The Structural Engineer 79 (6) (2001) 27–31. [13] P.E. McSharry, G.D. Clifford, L. Tarassenko, L.A. Smith, A dynamical model for generating synthetic electrocardiogram signals, IEEE Transactions of Biomedical Engineering 50 (2003) 289–294. [14] V. Racic, A. Pavic, Mathematical model to generate near-periodic human jumping force signals, Mechanical Systems and Signal Processing 24 (2010) 138–152, doi:10.1016/j.ymssp.2009.07.001. ¨ [15] H. Bachmann, W. Ammann, in: Vibrations in Structures Induced by Man and Machines, IABSE–AIPC–IVBH, Zurich, Switzerland, 1987. [16] B.R. Ellis, T. Ji, in: The Response of Structures to Dynamic Crowd Loads, Digest 426, 2004 ed, British Standards Institution, London, UK, 2004. [17] B.R. Ellis, T. Ji, in: The Response of Structures to Dynamic Crowd Loads, Digest 426, British Standards Institution, London, UK, 1997. [18] S. Zivanovic, A. Pavic, P. Reynolds, Probability-based prediction of multi-mode vibration response to walking excitation, Engineering Structures 29 (6) (2007) 942–954. [19] IStructE/DCLG/DCMS Working Group, in: Dynamic Performance Requirements for Permanent Grandstands Subject to Crowd Action: Recommendations for Management, Design and Assessment, The Institution of Structural Engineers, The Department for Communities and Local Government and The Department for Culture Media and Sport, London, 2008. [20] B.R. Ellis, T. Ji, Loads generated by jumping crowds: numerical modelling, Structural Engineer 82 (17) (2004) 35–40. [21] M. Kasperski and E. Agu, Prediction of crowd-induced vibrations via simulation, in: Proceedings of 23rd International Modal Analysis Conference (IMAC XXIII), 31 January–3 February 2005, Orlando, FL, USA. [22] J.G. Parkhouse, D.J. Ewins, Crowd-induced rhythmic loading, Structures and Buildings 159 (2006) 247–259. [23] C.T. Farley, R. Blickhan, J. Saito, C.R. Taylor, Hopping frequency in humans: a test of how springs set stride frequency in bouncing gaits, Journal of Applied Physiology 71 (1991) 2127–2132. [24] S.M. Pandit, S.-M. Wu, in: Time Series and System Analysis with Applications, Wiley, New York, USA, 1983. [25] D.E. Newland, in: An Introduction to Random Vibrations, Spectral and Wavelet Analysis, 3rd ed, Pearson Education Limited, Harlow, Essex, UK, 1993. [26] D.M. Bates, D.G. Watts, in: Nonlinear Regression and its Applications, Wiley, New York, 1998. [27] D.C. Montgomery, G.C. Runger, in: Applied Statistics and Probability for Engineers, 2nd ed, John Wiley & Sons, Inc., New York, NY, 1999. [28] R. Rodano, R. Squadrone, Stability of selected lower limb joint kinetic parameters during vertical jump, Journal of Applied Biomechanics 18 (2002) 83–89. [29] AMTI Product Manuals, Advanced Mechanical Technology, Watertown, USA, 2007. [30] J.D. Littler, Frequencies of synchronised human loading from jumping and stamping, The Structural Engineer (2003) 27–35. [31] J.G. Parkhouse and D.J. Ewins, Vertical dynamic loading produced by people moving to a beat, in: Proceedings of ISMA2004, Leuven, May 28–31, 2004. [32] MathWorks (2009) /www.mathworks.comS. [33] M. Hermanussen, H. Danker-Hopfe, G.W. Weber, Body weight and the shape of the natural distribution of weight, in very large samples of German, Austrian and Norwegian conscripts, International Journal of Obesity 25 (2001) 1550–1553. [34] A. Ebrahimpour, R.L. Sack, Modelling dynamic occupant loads, Journal of Structural Engineering 115 (6) (1989) 1476–1496. [35] T. Ji and B.R. Ellis, Evaluation of dynamic crowd effects for dance loads, International Association for Bridge and Structural Engineering (IABSE) Symposium: Structural Serviceability of Buildings, Goteburg, Sweden, 1993. [36] M. Willford, An investigation into crowd-induced vertical dynamic loads using available measurements, Structural Engineer 79 (12) (2001) 21–25.