Spring-Mass Parameter Estimation

Kurt Bryan and SIMIODE

A notebook to estimate spring and damping constants from experimental data. 


The Data: In this example the Excel file containing the data lies in "Base/spring_mass_data_clean.xls".

s = Import["Base/spring_mass_data_clean.xls"];
exceldat = s[[1]];

The data consists of n = 1460 data points in an N x 2 array, the time (seconds) and  position (meters)  of the mass every 1/50 of a second. The time starts at t = 0.82 seconds for the first position measurement.

The column headings are

exceldat[[1,1]]
exceldat[[1,2]]

Arrange the data into an n x 2 Table, with time rescaled to start at t = 0.0 (by subtracting 0.82 from time values).

n = 1460;
udat0 = Table[{exceldat[[j,1]]-0.82,exceldat[[j,2]]},{j,2,n+1}];
Tmax = udat0[[n,1]]

A plot of the position data over time.

ListPlot[udat0,AxesLabel->{"Time (seconds)","Position (meters)"},ImageSize->Scaled[0.8],Joined->True]

Recenter the data values so that the mean position is 0. First compute mean value of data over the total time interval (might be best to do this over an integer number of cycles):

uave = Sum[udat0[[j,2]],{j,1,n}]/n

Recenter data, amalgamated into an array "udat" consisting of (time, position) pairs:

udat = Table[{udat0[[j,1]],udat0[[j,2]]-uave},{j,1,n}];

Plot recentered data

plt1 = ListPlot[udat,AxesLabel->{"Time (seconds)","Position (meters)"},ImageSize->Scaled[0.8],Joined->True]

Estimate the spring and damping constants: We will fit a function y(t) of the form

y[t_] = d1*Exp[-alpha*t]*Cos[omega*t];

to this data. 

Based on the plot above we can guess that the initial amplitude d1 is about 0.05. We can estimate alpha by the rate of decay of the amplitude of the oscillations. For example, at time t = 25 the amplitude is down to about 0.035, so 0.05*exp(-alpha*25) = 0.035, which leads to alpha equal to about 0.14 (solve 0.05*exp(-alpha*25) = 0.035 for alpha).

We can estimate omega by estimating the period of the motion, e.g., count how many complete oscillations the mass undergoes during the approximate 29 second data set. We can then plot y(t) with these estimated values, as

plt2 = Plot[y[t]/.{d1->0.05,alpha->0.014,omega->7.0},{t,0,Tmax},PlotStyle->{Red},ImageSize->Scaled[0.8]];
Show[plt1,plt2]

It may help to plot on a smaller time range, at first:

Show[plt1,plt2,PlotRange->{{0,5},{-0.05,0.05}}]

Obviously some adjustment in omega and perhaps the other parameters is in order.

Adjust d1, alpha, and omega to obtain the best (visual) fit possible, then use the formulas in Modeling Exercise 6.3.5 to estimate the spring and damping constant.