Application Center - Maplesoft
  Submit your work
Application Center Applications Lesson 30 -- Application: Nonlinear Springs Preview

App Preview:

Lesson 30 -- Application: Nonlinear Springs

You can switch back to the summary page by clicking here.

Learn about Maple
Download Application

 

Lesson30.mw

ORDINARY DIFFERENTIAL EQUATIONS POWERTOOL

Lesson 30 -- Application: Nonlinear Springs

Prof. Douglas B. Meade

Industrial Mathematics Institute

Department of Mathematics

University of South Carolina

Columbia, SC 29208

URL:   http://www.math.sc.edu/~meade/

E-mail: meade@math.sc.edu

Copyright  2001  by Douglas B. Meade

All rights reserved

-------------------------------------------------------------------

>

Outline of Lesson 30

30.A Model for a General Spring-Mass System with Damping

                30.A-1 Linear Spring

                30.A-2 Soft Spring

                30.A-3 Hard Spring

>

Initialization

> restart;

> with( DEtools ):

> with( plots ):

> with( student ):

> with( LinearAlgebra ):

> with( VectorCalculus ):

Warning, the name changecoords has been redefined

Warning, the names `&x`, CrossProduct and DotProduct have been rebound

Warning, the assigned names `<,>` and `<|>` now have a global binding

Warning, these protected names have been redefined and unprotected: `*`, `+`, `.`, D, Vector, diff, int, limit, series

>

30.A Model for a General Spring-Mass System with Damping

As in Lesson 28, consider a spring, not necessarily linear, suspended vertically from a support.  Instead of Hooke's linear law for the spring, assume that the force of the spring on the mass is given by a function such as S(x) .  A mass m is attached to the lower end of the spring, and the mass is allowed to sink to its equilibrium position where the origin of a coordinate system is established.  Measure displacements from this origin by x(t) taken positive upward.  Finally, add linear damping to the system so that Newton's Second Law of Motion generates the ODE

m d^2*x/(dt^2)+b dx/dt-S(x(t)) = 0

Note that the displacement must remain with certain bounds to remain physically realistic. If the displacement becomes too large, the spring will stretch or even break. The limit on negative displacement from compression is the point at which the coils of the spring bunch together.

If this second-order ODE is written as

> ode1 := m * diff( x(t), t$2 ) + b * diff( x(t), t ) - S(x(t)) = 0;

ode1 := m*(diff(x(t), `$`(t, 2)))+b*(diff(x(t), t))-S(x(t)) = 0

>

the corresponding first-order system would be

> U := [ x(t), v(t) ]:

> f := (x,v) -> v:

> g := (x,v) -> ( S(x) - b*v ) / m:

> sys1 := equate( diff( U, t ), [ f( x(t), v(t) ), g( x(t), v(t) ) ] );

sys1 := {diff(x(t), t) = v(t), diff(v(t), t) = (S(x(t))-b*v(t))/m}

>

The general condition for an equilibrium solution is

> q1 := solve( eval( sys1, { x(t)=x0, v(t)=v0 } ), { x0, v0 } );

q1 := {v0 = 0, x0 = RootOf(S(_Z))}

>

The most note-worthy property of this result is that the equilibrium solutions are independent of the damping force.

>

30.A-1 Linear Spring

Hooke's Law assumes the spring force is linear in the displacement so that it is

> SS[linear] := x -> -k*x;

SS[linear] := proc (x) options operator, arrow; VectorCalculus:-VectorCalculus(VectorCalculus:-VectorCalculus(k, x), -1) end proc

>

with k, m, b, g all positive.

There is only one equilibrium solution to the linear system, namely,

> S := SS[linear]:

> equil[linear] := q1;

equil[linear] := {v0 = 0, x0 = 0}

>

The nature of this equilibrium point is determined by the calculations

> for X0 in [equil[linear]] do

>  bb := eval( [ f(x0,v0), g(x0,v0) ], X0 );

>  A := eval( Jacobian( [f(x0,v0),g(x0,v0)], [x0,v0] ), X0 );

>  linsys1 := equate( diff( U, t ), evalm( A . U + bb ) );

>  print( `-------------------------------------------------------------` );

>  print( `Linearization at [x0,v0]`=eval([x0,v0],X0) );

>  print( 'A' = A );

>  print( lambda = Eigenvalues( A ) );

>  print( `-------------------------------------------------------------` );

> end do:

`-------------------------------------------------------------`

`Linearization at [x0,v0]` = [0, 0]

A = Matrix([[0, 1], [-k/m, -b/m]])

lambda = Vector[column]([[-1/2*(b-(b^2-4*m*k)^(1/2))/m], [-1/2*(b+(b^2-4*m*k)^(1/2))/m]])

`-------------------------------------------------------------`

>

It is not difficult to see that there are three cases to consider for these eigenvalues.  These cases correspond to underdamping, critical damping, and overdamping studied in Lesson 28.

Case 1: 0 < b <  2*sqrt(m*k)

The eigenvalues are complex conjugates with negative real part. Hence the equilibrium is a stable spiral.

This is the underdamped case considered in Lesson 28, Section B-1.

>

Case 2: b = 2*sqrt(m*k)

With a repeated negative eigenvalue, the equilibrium is a stable node (sink).

This is the critically damped case discussed in Lesson 28, Section B-2.

>

Case 3: b > 2*sqrt(m*k)

Both eigenvalues are negative so the equilibrium is a stable node (sink).

This is the overdamped case discussed in greater detail in Lesson 28, Section B-3.

>

Thus, in all cases, the linear spring-mass system returns to its equilibrium position.

>

This linear spring model, i.e., Hooke's Law, is valid only for "small" displacements. Beyond these limits, the linear model ceases to provide a good model of the spring's response.

Nonlinear spring forces can have wider applicability. For a nonlinear spring with spring force S(x) , the stiffness of the spring is -d/dx S(x) , a definition that is consistent with the linear spring where S(x) = -k*x .

If the stiffness decreases as a function of displacement, its derivative would be negative, and the spring is said to be soft. Given the definition of stiffness, the soft spring is characterized by -d^2/(dx^2) S(x) < 0 , or `S''`(x)*`>`*0 .

If the stiffness increases as a function of displacement, its derivative would be positive, and the spring is said to be hard.  Given the definition of stiffness, the hard spring is characterized by -d^2/(dx^2) S(x)*`>`*0 , or `S''`(x) < 0 .

With both k and alpha positive, simple models for soft and hard springs are

> SS[soft] := x -> -k*x + alpha*x^3;
SS[hard] := x -> -k*x - alpha*x^3;

SS[soft] := proc (x) options operator, arrow; VectorCalculus:-VectorCalculus(VectorCalculus:-VectorCalculus(VectorCalculus:-VectorCalculus(k, x), -1), VectorCalculus:-VectorCalculus(alpha, x^3)) end p...

SS[hard] := proc (x) options operator, arrow; VectorCalculus:-VectorCalculus(VectorCalculus:-VectorCalculus(VectorCalculus:-VectorCalculus(k, x), -1), VectorCalculus:-VectorCalculus(VectorCalculus:-Ve...

>

30.A-2 Soft Spring

Because of the complexities of the general solution to a cubic equation, this discussion is restricted to a specific example for which the parameter values are

> param1 := { m=1, b=2/10, k=10, alpha=2/10 };

param1 := {k = 10, m = 1, alpha = 1/5, b = 1/5}

>

The damped oscillator with the soft spring is then governed by the first-order system

> S := SS[soft]:

> sys[soft] := eval( sys1, param1 );

sys[soft] := {diff(x(t), t) = v(t), diff(v(t), t) = -10*x(t)+1/5*x(t)^3-1/5*v(t)}

>

The equilibrium solutions of this system are

> equil[soft] := allvalues( eval( q1, param1 ) );

equil[soft] := {v0 = 0, x0 = 0}, {x0 = 5*2^(1/2), v0 = 0}, {x0 = -5*2^(1/2), v0 = 0}

>

or, in floating-point form,

> evalf( equil[soft] );

{x0 = 0., v0 = 0.}, {v0 = 0., x0 = 7.071067810}, {v0 = 0., x0 = -7.071067810}

>

Thus, there are three equilibria, each classified by the linearization process computed below.

> for X0 in [equil[soft]] do

>  bb := eval( [ f(x0,v0), g(x0,v0) ], X0 union param1 );

>  A := eval( Jacobian( [f(x0,v0),g(x0,v0)], [x0,v0] ), X0 union param1 );

>  linsys1 := equate( diff( U, t ), evalm( A . U + bb ) );

>  print( `-------------------------------------------------------------` );

>  print( `Linearization at [x0,v0]`=eval([x0,v0],X0 union param1) );

>  print( 'A' = A );

>  print( lambda = Eigenvalues( A ) );

>  print( `-------------------------------------------------------------` );

> end do:

`-------------------------------------------------------------`

`Linearization at [x0,v0]` = [0, 0]

A = Matrix([[0, 1], [-10, (-1)/5]])

lambda = Vector[column]([[-1/10+3/10*I*111^(1/2)], [-1/10-3/10*I*111^(1/2)]])

`-------------------------------------------------------------`

`-------------------------------------------------------------`

`Linearization at [x0,v0]` = [5*2^(1/2), 0]

A = Matrix([[0, 1], [20, (-1)/5]])

lambda = Vector[column]([[-1/10+1/10*2001^(1/2)], [-1/10-1/10*2001^(1/2)]])

`-------------------------------------------------------------`

`-------------------------------------------------------------`

`Linearization at [x0,v0]` = [-5*2^(1/2), 0]

A = Matrix([[0, 1], [20, (-1)/5]])

lambda = Vector[column]([[-1/10+1/10*2001^(1/2)], [-1/10-1/10*2001^(1/2)]])

`-------------------------------------------------------------`

>

Based on the eigenvalues, ``(0, 0) is a stable spiral and ``(5*sqrt(2), 0) and ``(-5*sqrt(2), 0) are saddle points.

Figure 30.1 containing the equilibrium solutions, representative trajectories, and the direction field for this example concludes the discussion of soft springs.

> Pequil := pointplot( [seq( eval( [x0,v0], E ), E=[equil[soft]] )],

>                     color=BLUE, symbolsize=18 ):

> icA := [ [ x(0)=-5, v(0)=-4 ], [ x(0)=-8, v(0)=8] ]:

> PsolA := DEplot( sys[soft], [ x(t), v(t) ], t=-0.1..10, icA, linecolor=green,

>                 stepsize=0.05, arrows=NONE ):

> icB := [ [ x(0)=-6, v(0)=-4], [ x(0)=-9, v(0)=8] ]:

> PsolB := DEplot( sys[soft], [ x(t), v(t) ], t=-0.2..1, icB, linecolor=green,

>                 stepsize=0.05, arrows=NONE ):

> icC := [ [ x(0)=-10, v(0)=25] ]:

> PsolC := DEplot( sys[soft], [ x(t), v(t) ], t=-0.1..2.5, icC, linecolor=green,

>                 stepsize=0.1, x=-15..15, v=-40..40, arrows=NONE ):

> icD := [ [ x(0)=-10, v(0)=30] ]:

> PsolD := DEplot( sys[soft], [ x(t), v(t) ], t=-0.1..1.1, icD, linecolor=green,

>                 stepsize=0.1, arrows=NONE ):

> icE := [ [ x(0)=10, v(0)=-10] ]:

> PsolE := DEplot( sys[soft], [ x(t), v(t) ], t=-0.1..0.4, icE, linecolor=green,

>                 stepsize=0.05, arrows=NONE ):

> icF := [ [ x(0)=10, v(0)=-15] ]:

> PsolF := DEplot( sys[soft], [ x(t), v(t) ], t=-0.2..1.8, icF, linecolor=green,

>                 stepsize=0.1, x=-15..15, v=-40..40, arrows=NONE ):

> Pdfield := DEplot( sys[soft], [ x(t), v(t) ], t=0..1,

>                   x=-15..15, v=-40..40, dirgrid=[20,20], arrows=medium ):

> display( Pequil, PsolA, PsolB, PsolC, PsolD, PsolE, PsolF, Pdfield, view=[-15..15, -40..40 ], title="Figure 30.1" );

[Plot]

>

Note that the solution curves that become unbounded cease to be physically realistic when they exceed the spring's ability to compress or extend.

>

30.A-3 Hard Spring

Repeating the analysis for a hard spring, the first-order system is

> S := SS[hard]:

> sys[hard] := eval( sys1, param1 );

sys[hard] := {diff(x(t), t) = v(t), diff(v(t), t) = -10*x(t)-1/5*x(t)^3-1/5*v(t)}

>

and its equilibrium solutions are

> q2 := evalf( allvalues( eval( q1, param1 ) ) );

q2 := {x0 = 0., v0 = 0.}, {x0 = 7.071067810*I, v0 = 0.}, {v0 = 0., x0 = -7.071067810*I}

>

Note that only one of these solutions is real; the complex-valued solutions are not physical.  This real equilibrium is

> equil[hard] := op(remove( has, [q2], I ));

equil[hard] := {x0 = 0., v0 = 0.}

>

Linearization about this equilibrium yields

> for X0 in [equil[hard]] do

>  bb := eval( [ f(x0,v0), g(x0,v0) ], X0 union param1 );

>  A := eval( Jacobian( [f(x0,v0),g(x0,v0)], [x0,v0] ), X0 union param1 );

>  linsys1 := equate( diff( U, t ), evalm( A . U + bb ) );

>  print( `-------------------------------------------------------------` );

>  print( `Linearization at [x0,v0]`=eval([x0,v0],X0 union param1) );

>  print( 'A' = A );

>  print( lambda = Eigenvalues( A ) );

>  print( `-------------------------------------------------------------` );

> end do:

`-------------------------------------------------------------`

`Linearization at [x0,v0]` = [0., 0.]

A = Matrix([[0, 1], [-10., (-1)/5]])

lambda = Vector[column]([[-0.99999999999999990e-1+3.16069612585582106*I], [-0.99999999999999990e-1-3.16069612585582106*I]])

`-------------------------------------------------------------`

>

Therefore, the sole equilibrium is a stable spiral.

Figure 30.2 containing the equilibrium solutions, a representative trajectory, and the direction field for this example concludes the discussion of soft springs.

> Pequil := pointplot( [seq( eval( [x0,v0], E ), E=[equil[hard]] )],

>                     color=BLUE, symbolsize=18 ):

> icA := [[ x(0)=-9, v(0)=0]]:

> PsolA := DEplot( sys[hard], [ x(t), v(t) ], t=0..10, icA, linecolor=green,

>                 stepsize=0.025, arrows=NONE ):

> Pdfield := DEplot( sys[hard], [ x(t), v(t) ], t=0..1,

>                   x=-15..15, v=-40..40, dirgrid=[20,20], arrows=medium ):

> display( Pequil, PsolA, Pdfield, view=[-15..15, -40..40 ], title="Figure 30.2" );

[Plot]

>

Now, all solutions appear to remain bounded and ultimately are attracted to the equilibrium solution. This illustrates a key qualitative difference between solft and hard springs.

>

[Back to ODE Powertool Table of Contents]

>