# Math 2250 Maple Project 2,  March 2003. Tacoma Narrows.
# 
# NAME ________________________________________________  CLASSTIME
# ______ 
# 
# There are six (6) problems in this project. Please answer the
# questions A, B, C , ... associated with
#  each problem.  The original worksheet "project2-spring-2003.mws" is a
# template for the solution;
#  you must fill in the code and all comments. Sample code can be copied
# with the mouse. Use pencil
#  freely to annotate the worksheet and to clarify the code and figures.
# 
# The problem headers for the Spring 2003 revision of David Eyre's
# project (original was year 2000). 
#  __________ 2.1. OVERDAMPED FREE OSCILLATIONS.                     The
# first four problems study the
#  __________2.2. UNDERDAMPED FREE OSCILLATIONS.                  
# simplest linear model for x(t).
#  __________2.3. UNDAMPED FORCED OSCILLATIONS (c=0).
#  __________2.4. DAMPED FORCED OSCILLATIONS (c>0).
#  __________2.5. LARGE SUSTAINED OSCILLATIONS.                      
# This is the second, nonlinear model.
#  __________2.6. MCKENNA NON-HOOKES LAW CABLE MODEL.  This is the
# third, nonlinear model.
# 
# 
# 2.1.  PROBLEM (OVERDAMPED FREE OSCILLATIONS)
# 
# FREE OSCILLATIONS. Consider the general problem of free linear
# oscillations
# 
#            m x'' + c x' + k x=0,
#            x(0)=x0,  x'(0)=v0,
# 
#  where  m=1,  k=16 and  c  is a non-negative constant.  The symbols  
#  x0  and  v0  are the initial position and initial velocity, 
# respectively.
# 
#     A. Suggest a value for parameter c  > 0 so that the free
# oscillations are overdamped. This
#          value will be used in item B below. Check your answer by
# solving the characteristic equation
#          using Maple's "solve" command (example: solve(r^2-1=0,r);).
#     
#     B. Use x(0)=1 and x'(0)= -2 for the initial conditions and Maple's
# "dsolve" to find the explicit
#         real solution x(t). Plot the solution x(t) for t=0 to t=5
# using Maple's "plot" command.
# 
# EXAMPLE(Wrong parameters! Change it!)
# de:=3*diff(x(t),t,t)+1.5*diff(x(t),t)+4*x(t)=0:  # Define the
# differential equation
# ic:=x(0)=0,D(x)(0)= -1:                                 # Define the
# initial conditions
# dsolve({de,ic},x(t),method=laplace);              # Symbolically solve
# for x(t)
# X:=unapply(rhs(%),t):                                    # Capture the
# dsolve symbolic answer as a function X(t)
# plot(X(t),t=0..5);                                           # Plot
# the solution
> 
> #2.1-A
> # overdamped means mr^2+cr+k=0 has two real roots.
> #2.1-B         
> 
# 2.2. PROBLEM (UNDERDAMPED FREE OSCILLATIONS)
# 
# FREE OSCILLATIONS. Consider the problem of free linear oscillations
# 
#            m x'' + c x' + k x=0,
#            x(0)=0,  x'(0)=1,
# 
#  where  m=2,  c=4 and  k  is a non-negative constant. The overdamped
# case
#  is studied here, for which the solution x(t) has infinitely many
# oscillations. 
# 
#     A.  Find a  Hooke's constant k > 5  so that the solution x(t)
# changes
#          sign infinitely many times and decays to zero at t=infinity.
# Display
#          the value of  k  and the exact symbolic solution.
# 
#     B. Plot the exact symbolic solution x(t) from t=0 to t=5. Estimate
# from the 
#         graph the decimal value of the pseudoperiod. Display the
# graphical estimate
#         and also the exact pseudoperiod 2Pi/w, where  w  is the
# natural frequency
#         of the trigonometric term in the solution  x(t)  found in item
# 2.2.A.
# 
# Maple tip: Click with the mouse on the graphic to print the cursor
# location (left upper 
# corner of the maple window). The coordinates printed are of the form
# (x,y). From this
#  coordinate information, a simple subtraction estimates the period
> 
> #2.2-A Define  k, then solve and plot.
> #2.2-B 
> 
# 2.3. PROBLEM (UNDAMPED FORCED OSCILLATIONS (c=0))
# 
# FORCED LINEAR OSCILLATIONS.  Consider the forced problem
# 
#       x'' +  1.5625 x = 10 cos(wt),
#       x(0)=0,  x'(0)=0,
# 
#      where  w is a  non-negative constant.
#   
#     A.  Choose  w=3.75,  so that the forcing  frequency w is 3 times
# larger 
#           than the natural frequency w0=1.25. Solve for x(t) using
# dsolve().
#           Plot the solution x(t) on a suitable interval in order to
# show the 
#           global behavior of  the solution x(t).
# 
#     B. The solution x(t) is the sum of two functions, one of period 
# 2Pi/w  and the 
#          other of period  2Pi/w0. Approximate the period of x(t) by
# examining the 
#          graph in 2.3.A. Display the estimate for the period and also
# the exact period, 
#          as calculated from the solution formula for  x(t) -- see page
# 341 for details 
#          about how to calculate the exact period.
# 
#     C. Suggest a value for the forcing frequency  w  so that the
# oscillations exhibit
#          resonance.  Show resonant behavior on a graph.
> #   
> #2.3-A 
> #2.3-B
> #2.3-C
> 
# 2.4. PROBLEM (DAMPED FORCED OSCILLATIONS (c>0))
#    Consider the forced problem
# 
#       x'' + c x' + 21 x = 5cos(w t),
#       x(0)=0,  x'(0)=0,
# 
#     A. Consider the damping constants  c=2, c=1 and c=1/2.  Compute
# the
#         amplitude function  C(w)  [page 346] for these three
# equations, then plot
#         for w=0 to w=20 the three amplitude graphs on a single set of
# axes.
# 
#     B.  For each case  c=2, c=1, c=1/2, print the values  w*, C* 
# where 
#          C*=C(w*)=max {C(w) : 0 <= w <= 20}. The three data pairs
# should 
#          show that C* becomes larger as c tends to zero. 
#          SAVE YOUR MAPLE FILE FREQUENTLY
# 
#          Maple Hint: Use Maple's mouse interface on the graphic of
# Part C. 
#          Specifically, click on a possible maximum (horizontal
# tangent) in the 
#          graph to display the values  w*, C* on the screen. Copy the
# values on paper.
# 
# EXAMPLE(Beware! Wrong values!)
# F:=15: m:=1: k:=25: c:='c': w:='w':
# C:=(w,c)->F/sqrt((k-m*w*w)^2+(c*w)^2):
# plot({C(w,4),C(w,3),C(w,2)},w=0..15,color=black);
> 
> #2.4-A Plot C(w), three graphics on one set of axes
> #2.4-B Table of six data values for w*, C*
> 
# 2.5. PROBLEM (LARGE SUSTAINED OSCILLATIONS)
# 
# NONLINEAR MODEL WITH GEOMETRY INCLUDED.
# Consider the nonlinear, forced, damped oscillator equation for
# torsional motion, 
# with bridge geometry included,
# 
#     x'' + 0.05 x' + 2.4 sin(x)cos(x) = 0.06 cos (12 t/10) ,
#     x(0) = x0,  x'(0) = v0
# 
# and its corresponding linearized equation
# 
#    x'' + 0.05 x' + 2.4 x = 0.06 cos (12 t/10) ,  
#    x(0) = x0,  x'(0) = v0.
# 
# The spring-mass system parameters are m=1, c = 0.05,  k = 2.4,  w =
# 1.2 ,  F = 0.06.
# Maple code used to solve and plot the solutions appears below.
# 
#    # WARNING: set the parameters on the second line! Use "copy as
# maple text" for maple 6+.
#    m:=1:  F := 0.06:  w := 1.2: m:=1: c:= 0.05: k:= 2.4:  
#    x0:=0: v0:=0: a:=0: b:=50:
#    deNonLinear:= m*diff(x(t),t,t) + c*diff(x(t),t) +
# k*sin(x(t))*cos(x(t)) = F*cos(w*t):
#    deLinear:= m*diff(x(t),t,t) + c*diff(x(t),t) + k*x(t) = F*cos(w*t):
#    with(DEtools):  opts:=stepsize=0.1:
#   
# DEplot(deNonLinear,x(t),t=a..b,[[x(0)=x0,D(x)(0)=v0]],opts,title='NonL
# inear');
#   
# DEplot(deLinear,x(t),t=a..b,[[x(0)=x0,D(x)(0)=v0]],opts,title='Linear'
# );
# 
# 
#     A. Let  x0=0,  v0=0.  Plot the solutions of the linear and
# nonlinear equations from t=160 to t=260.
#          These plots represent the steady state solutions of the two
# equations.
# 
#    B. Let  x0=1.2,  v0=0. Plot the solutions of the linear and
# nonlinear equations from t=220 to t=320.
#        These plots represent the steady state solutions of the two
# equation, with new starting value x0=1.2.
# 
#     C.  Argue in a sentence why the two linear plots have to be
# identical, based upon the 
#           superposition formula x(t)=xh(t)+xss(t), even though the
# homogeneous solution 
#           xh(t) is different for the two plots. Please include a
# discussion of the amplitude of  xh(t) 
#           on the corresponding  t-interval. 
# 
#     D.  Determine the ratio of the apparent amplitudes (a number > 1)
# for the nonlinear plots.
#          Explain why "large sustained oscillations" is an appropriate
# description of the
#          nonlinear steady-state behavior.
# 
> 
> #2.5-A
> #2.5-B
> #2.5-C
> # The plots are identical because ....
> #2.5-D
> 
# 2.6.  PROBLEM. ( MCKENNA'S NON-HOOKE'S LAW CABLE MODEL)
# 
#  MCKENNA'S NON-HOOKE'S LAW CABLE MODEL FOR THE TACOMA NARROWS BRIDGE
# 
# The model of McKenna studies the bridge with a nonlinear, forced,
# damped 
# oscillator equation for torsional motion that accounts for the
# non-Hooke's law
# cables coupled to the equations for vertical motion.  The equations in
# this
# case couple the torsional motion with the vertical motion.  The
# equations are:
# 
#     x'' + c x' - k G(x,y) = F sin wt,   x(0) = x0,  x'(0) = x1,
#     y'' + c y' + (k/3) H(x,y) = g ,     y(0) = y0,  y'(0) = y1,
# 
# where x(t) is the torsional motion and y(t) is the vertical motion.
# The functions
# G(x,y) and H(x,y) are the models of the force generated by the cable
# when 
# it is contracted and stretched. Below is sample code for writing the
# differential
# equations and for plotting the solutions. It is ready to copy with the
# mouse.
# 
# with(DEtools):                                                        
#  
# w := 1.3:  F := 0.05:  f(t) := F*sin(w*t):                            
#  
# c := 0.01:  k1 := 0.2:  k2 := 0.4:  g := 9.8:  L := 6:                
#  
# STEP:=x->piecewise(x<0,0,1):
# fp(t) := y(t)+(L*sin(x(t))):
# fm(t) := y(t)-(L*sin(x(t))):
# Sm(t) := STEP(fm(t))*fm(t):
# Sp(t) := STEP(fp(t))*fp(t):
# sys := {             
#      diff(x(t),t,t) + c*diff(x(t),t) -
# k1*cos(x(t))*(Sm(t)-Sp(t))=f(t),
#      diff(y(t),t,t) + c*diff(y(t),t) + k2*(Sm(t)+Sp(t)) = g}:
# ic := [[x(0)=0, D(x)(0)=0, y(0)=27.25, D(y)(0)=0]]:
# vars:=[x(t),y(t)]:   
# opts:=stepsize=0.1:  
# DEplot(sys,vars,t=0..300,ic,opts,scene=[t,x]);
# 
# The amazing thing that happens in this simulation is that the large
# vertical  oscillations take 
# all the tension out of the springs and they induce large torsional
# oscillations.
# 
#  A. TORSIONAL OSCILLATION PLOT. Get the sample code above to produce
# the plot of x(t) 
#       [that's what scene=[t,x] means]. 
# 
#  B. Estimate the number of degrees the roadway oscillates based on the
# plot; recall that x in the
#      plot is reported in radians. Comment on the agreement of this
# result with historical data.
# 
#      Tip: Average the five largest amplitudes in the plot to find an
# average maximum amplitude
#      for t=0 to t=300. Convert to degrees using  Pi  radians = 180
# degrees.
# 
#  C. VERTICAL OSCILLATION PLOT. Modify the DEplot code to scene=[t,y]
# and plot the
#       oscillation y(t) on t=0 to t=300. The plot is supposed to show
# 30-foot vertical oscillations
#       that dampen to 7-foot vertical oscillations after 300
# seconds.Comment on the agreement
#       between these oscillation results and the historical data for
# Tacoma Narrows, especially the
#       visual data present in the film clip of the bridge disaster.
#       
> 
> #2.6-A Torsional plot t-versus-x
> #2.6-B  Roadway oscillation estimate in degrees + comments.
> #2.6-C Vertical plot t-versus-y + comments.
>