State-Space model of spring-mass system

[Number2Pencil]

Homework Statement

Derive the state space model of a spring-mass-damper system. Neglect friction, wind resistance, etc. Neglect the force of gravity. Assume the spring's equilibrium point is at y=0 and there is an arbitrary

Homework Equations

[X'] = AX + BU
y = CX + DU

The Attempt at a Solution

Newton's law:
$$F_{NET}=ma = m\ddot{y}$$

The sum is an arbitrary input force, u, the recoil force that works in the direction towards the spring's equilibrium point, and the dampening coefficient which works against velocity:
$$u-c\dot{y}-ky=m\ddot{y}$$

Rearranging terms:
$$\ddot{y} = \frac{1}{m}u - \frac{c}{m}\dot{y}-\frac{k}{m}y$$

Setting up state-space integral form:
$$y = \int( -\frac{c}{m}y + \int (-\frac{k}{m}y + \frac{1}{m}u) dt) dt'$$I call the entire outer integral X1 and the inner integral X2, using this nomenclature:

$$y = X_1$$

$$\dot{X_1} = -\frac{c}{m}X_1 + X_2$$

$$\dot{X_2} = -\frac{k}{m}X1 + \frac{1}{m}u$$

$$\left[ \begin{array}{cc} \dot{X1}\\ \dot{X2}\\ \end{array} \right] = \left[ \begin{array}{cc} -\frac{c}{m}&1\\ -\frac{k}{m}&0\\ \end{array} \right] \left[ \begin{array}{cc} X1\\ X2\\ \end{array} \right] + \left[ \begin{array}{cc} 0\\ \frac{1}{m} \end{array} \right]u$$$$y = \left[ \begin{array}{cc} 1&0 \end{array} \right]\left[ \begin{array}{cc} X1\\ X2\\ \end{array} \right] + 0u$$So now I wanted to test it using a sine-wave input force to the system. Intuitively, I would expect a sine-wave output as the mass should be bobbing up and down on the spring. What I actually get is at the beginning, there is a much slower frequency component than my input sine-wave, as well as the input-frequency component, almost like the sum of two different sine waves. This slower-component then dies off and I am left with something that matches more closely to my input sine-wave. What is strange is that during that initial slow-frequency response, the amplitude gets quite large, much larger than when it is finally settled.

My question is: What is going on? Is my intuition just wrong and truly this mass has a large initial "push"? Is there an error in my state-space model? Or something I am forgetting to take into account. I tested that having a zero-force input, the spring remains at zero.

Edit: Also, my test constant values were m= 0.72kg, k = 0.5, c =1

 

[milesyoung]

Seems fine to me, you're solving it numerically I guess? MATLAB/Simulink?

Maybe you could post the code you're running etc.

I'm off to bed..

 

[Number2Pencil]

Sure. Here is the MATLAB code I used to produce the time-domain plots. I had trouble with simulink...so I really can't verify if this is implementing state-space correctly.

%Physical Inputs
mass = 0.72;
springK = 1;
damperC = 0.5;

%Input Parameters
sineWaveAmp = 1;
sineWaveFreq = 0.5;

%Timing Parameters
timeStart=0;
timeStop=50;
dt = 0.001;

%Create State-Space Matrices
A=[-damperC/mass, 1; -springK/mass, 0];
B=[0; 1/mass];
C=[1,0];
D=[0];
stateSpace = [0;0];

%Loop Variables
numSamps = (timeStop-timeStart)/dt;
currTime = timeStart;

%Create Storage Elements
posVals(1:numSamps) = 0;
timeVals(1:numSamps) = 0;
u(1:numSamps) = 0;


%Step 1, create sine wave input values
sineWavePeriod = 1/sineWaveFreq;
for n=1:numSamps
	cycleTime = mod(currTime,sineWavePeriod);
	u(n) = sineWaveAmp*sin(2*pi*sineWaveFreq*cycleTime);
	currTime = currTime + dt;
end

%Step 2, simulate state space
currTime = timeStart;
for n = 1:numSamps

	% Y = Cx' + Du
	posVals(n) = (C*stateSpace)+(D*u(n));
	integratorInput = (A*stateSpace)+(B*u(n));
	
	%Integrate. x' = int(A*x + B*u)
	stateSpace = stateSpace+(integratorInput*dt);

	%Prepare next iteration
	timeVals(n) = currTime;
	currTime = currTime + dt;
end


%Step 3, plot
figure(1);
plot(timeVals,posVals);
title('Mass-Spring Position over time');
ylabel('Meters');
xlabel('Seconds');

figure(2);
plot(timeVals,u);
title('Mass-Spring External Input over time');
ylabel('N');
xlabel('Seconds');

 

[milesyoung]

I did a quick comparison of your solution with that of the Dormand-Prince solver in the MATLAB Control Systems Toolbox and everything checks out, so good job

The transient you see at the start is the sum of the decaying natural response and the forced response.

Is there a reason for using that particular state vector? I'd have thought that $\begin{bmatrix}y & \dot{y}\end{bmatrix}^T$ would be more convenient.

 

[Number2Pencil]

Alright, thanks miles. I guess my intuition was just off :p

For choosing that state vector, it's simply a product of me being very methodical on solving these kinds of things.

I'll check out this Dormand-Prince toolbox you mentioned, thank again