\documentclass{article}
\marginparwidth 0pt \oddsidemargin  0pt \evensidemargin  0pt 
\marginparsep 0pt
\topmargin   -1in \textwidth   6.5in \textheight  8.5 in
\begin{document}
\small\rm
\begin{verbatim}
 Math 2280
 Maple Lab 8: Earthquake project
 S2015
 
 Name _____________________________________   Class Time __________
 
 Project 8.  Solve problems L8-1 to L8-5. The problem headers:
 
 _______     PROBLEM L8.1. EARHQUAKE MODEL FOR A BUILDING.
 _______     PROBLEM L8.2. TABLE OF NATURAL FREQUENCIES AND PERIODS.
 _______     PROBLEM L8.3. UNDETERMINED COEFFICIENTS STEADY-STATE SOL
 _______     PROBLEM L8.4. PRACTICAL RESONANCE.
 _______     PROBLEM L8.5. EARTHQUAKE DAMAGE.
 
 FIVE FLOOR Model.
 Refer to the textbook of Edwards-Penney, section 5.3 Application (after
the section 5.3 exercises).
 Consider a building with five floors each weighing 50 tons. Each floor
 corresponds to a restoring Hooke's force with constant k=5 tons/foot.
 Assume that ground vibrations from the earthquake are modeled by
 (1/4)cos(wt) with period T=2*Pi/w.
 
 
 PROBLEM L8-1. BUILDING MODEL FOR AN EARTHQUAKE.
 Model the 5-floor problem in Maple. 
 Define the 5 by 5 mass matrix M and Hooke's matrix K for this system
 and convert Mx''=Kx into the system x''=Ax where A is defined by
 textbook equation (1), section 5.3 Application.
 Sanity check:  Mass m=3125, and the 5x5 matrix contains fraction 16/5.
 
 Then find the eigenvalues of the matrix A to six digits, using the 
 Maple command "linalg[eigenvals](A)."
 Sanity check: All six eigenvalues should be negative.
 
 # Sample Maple code for a model with 4 floors.
 # Use maple help to learn about evalf and eigenvals.
 # A:=matrix([ [-20,10,0,0], [10,-20,10,0],
 [0,10,-20,10],[0,0,10,-10]]);
 # with(linalg): evalf(eigenvals(A));
 
 
 # Problem L8.1
 # Define k, m and the 5x5 matrix A.
 # with(linalg): evalf(eigenvals(A));
 
 PROBLEM L8-2. TABLE OF NATURAL FREQUENCIES AND PERIODS.
 Refer to figure 5.3.17 in Edwards-Penney.
 
 Find the natural angular frequencies  omega=sqrt(-lambda) for the
 five story building and also the corresponding periods
 2PI/omega, accurate to six digits.  Display the answers in a table .
 Compare with answers in Figure 5.3.17 (actually a table), for the 7-story case.
 
 # Sample code for a 4x3 table, 4-story building.
 # Use maple help to learn about nops and printf.
 # ev:=[-10,-1.206147582,-35.32088886,-23.47296354]: n:=nops(ev):
 # Omega:=lambda -> sqrt(-lambda):
 # format:="%10.6f  %10.6f  %10.6f\n": 
 # seq(printf(format,ev[i],Omega(ev[i]),2*evalf(Pi)/Omega(ev[i])),
 i=1..n);
 
  # Problem L8-2
 # ev:=[fill this in]: n:=nops(ev):
 # Omega:=lambda -> sqrt(-lambda): format:="%10.6f  %10.6f  %10.6f\n":
 # seq(printf(format,ev[i],Omega(ev[i]),2*evalf(Pi)/Omega(ev[i])),i=1..n);
 
 PROBLEM  L8-3. UNDETERMINED COEFFICIENTS 
               STEADY-STATE PERIODIC SOLUTION.
 Consider the forced equation x'=Ax+cos(wt)b where b is a constant
 vector. The earthquake's ground vibration is accounted for by the extra
 term cos(wt)b, which has period T=2Pi/w. The solution x(t) is the
 5-vector of excursions from equilibrium of the corresponding 5 floors.
 Sought here is not the general solution, which certainly contains
 transient terms, but rather the steady-state periodic solution, which
 is known from the theory to have the form x(t)=cos(wt)c for some vector
 c that depends only on A and b.
 
 Define  b:=0.25*w*w*vector([1,1,1,1,1]): in Maple and find the vector c
 in the undetermined coefficients solution x(t)=cos(wt)c.  Vector c
 depends on w. As outlined in the textbook, vector c can be found by
 solving the linear algebra problem -w^2 c = Ac + b; see equation (32),
 section 5.3. Don't print c, as it is too complex; instead, print c[1]
 as an illustration.
 
 # Sample code for defining b and A, then solve for c, 4-floor case.
 # See maple help to learn about vector and linsolve.
 #  w:='w': u:=w*w: b:=0.25*u*vector([1,1,1,1]):
 #  A:=matrix([ [-20,10,0,0], [10,-20,10,0],
 #            [0,10,-20,10],[0,0,10,-10]]);
 # Au:=evalm(A+u*diag(1,1,1,1));
 #  c:=linsolve(Au,-b): 
 #  evalf(c[1],2);
 
 # PROBLEM L8-3
 # Define w, u, b, A, Au, c  
 # evalf(c[1],2);
 
 PROBLEM L8-4. PRACTICAL RESONANCE.
 Consider the forced equation  x'=Ax+cos(wt)b  of L8-3 above with
 b:=0.25*w*w*vector([1,1,1,1,1]).
 Practical resonance can occur if a component of  x(t)  has large
 amplitude compared to the vector norm of b. For example, an earthquake
 might cause a small 3-inch excursion on level ground, but the
 building's floors might have 50-inch excursions, enough to destroy
 the building.
 
 Let Max(c) denote the maximum modulus of the components of vector c.
 Plot  g(T)=Max(c(w)) with w=(2*Pi)/T for periods T=1 to T=5, ordinates
 Max=0 to Max=10, the vector c(w) being the answer produced in L8.3
 above.  Compare your figure to the textbook Figure 5.3.18.
 
 # Sample maple code to define the function Max(c), 4-floor building.
 # Use maple help to learn about norm, vector, subs and linsolve.
 #  with(linalg):
 #  w:='w': Max:= c -> norm(c,infinity); u:=w*w:
 #  b:=0.25*w*w*vector([1,1,1,1]):
 # A:=matrix([ [-20,10,0,0], [10,-20,10,0], [0,10,-20,10],
 [0,0,10,-10]]);
 #  Au:=evalm(A+u*diag(1,1,1,1));
 #  C:=ww -> subs(w=ww,linsolve(Au,-b)):
 #  plot(Max(C(2*Pi/r)),r=1..5,0..10,numpoints=150);
 
 # PROBLEM L8.4. WARNING: Save your file often!!!
 # w:='w': Max:= c -> norm(c,infinity): u:=w*w:   
 # Define b
 # Define A
 # Define Au
 # Define C  
 # plot(Max(C(2*Pi/r)),r=1..5,0..10,numpoints=150);
 
 PROBLEM L8-5. EARTHQUAKE DAMAGE.
 The maximum amplitude plot of L8-4 can be used to detect the of
 earthquake damage for a given ground vibration of period T.  A ground
 vibration (1/4)cos(wt), T=2*Pi/w, will be assumed, as in L8-4.
 
 (a) Replot the amplitudes in L8-4 for periods 1.5 to 5.5 and amplitudes
     5 to 10. There will be several spikes.
 (b) Create several zoom-in plots, one for each spike, choosing a
     T-interval that shows the full spike. 
 (c) Determine from the several zoom-in plots approximate intervals for
     the period T such that some floor in the building will undergo
     excursions from equilibrium in excess of 5 feet.
 
 # Example: Zoom-in on a spike for amplitudes 5 feet to 10 feet,
 # periods 1.97 to 2.01. This example for the 4-floor problem.
 # with(linalg): w:='w': Max:= c -> norm(c,infinity); u:=w*w:
 # Au:=matrix([[-20+u,10,0,0],[10,-20+u,10,0],[0,10,-20+u,10],[0,0,10,-10+u]]);
 # b:=0.25*w*w*vector([1,1,1,1]):
 # C:=ww -> subs(w=ww,linsolve(Au,-b)):
 # plot(Max(C(2*Pi/r)),r=1.97..2,01,5..10,numpoints=150);
 
 # PROBLEM L8-5. WARNING: Save your file often!!
 #(a) Re-plot the five spikes.
 # plot(Max(C(2*Pi/r)),r=1.5..5.5,5..10,numpoints=150);
 #(b) Plot five zoom-in graphs.
 # one:=1.79..1.83:plot(Max(C(2*Pi/r)),r=one,5..10,numpoints=150);
 # two:=???:plot(Max(C(2*Pi/r)),r=two,5..10,numpoints=150);
 # three:=???:plot(Max(C(2*Pi/r)),r=three,5..10,numpoints=150);
 # four:=???:plot(Max(C(2*Pi/r)),r=four,5..10,numpoints=150);
 # five:=???:plot(Max(C(2*Pi/r)),r=five,5..10,numpoints=150);
 #(c) Print period ranges.
 # PeriodRanges:=[one,two,three,four,five];
\end{verbatim}
\end{document}
