# Math 2250 # Maple Lab 7: Earthquake project # S2016 # Project 7. Solve problems L7-1 to L7-5. The problem headers: # # _______ PROBLEM L7.1. EARHQUAKE MODEL FOR A BUILDING. # _______ PROBLEM L7.2. TABLE OF NATURAL FREQUENCIES AND PERIODS. # _______ PROBLEM L7.3. UNDETERMINED COEFFICIENTS STEADY-STATE SOL # _______ PROBLEM L7.4. PRACTICAL RESONANCE. # _______ PROBLEM L7.5. EARTHQUAKE DAMAGE. # # FIVE FLOOR Model. # Refer to the textbook of Edwards-Penney, section 5.4 Application (after the section 5.4 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 L7-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.4 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 L7-2. TABLE OF NATURAL FREQUENCIES AND PERIODS. # Refer to figure 5.4.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.4.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 L7-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.4. 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 L7-4. PRACTICAL RESONANCE. # Consider the forced equation x'=Ax+cos(wt)b of L7-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.4.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 L7-5. EARTHQUAKE DAMAGE. # The maximum amplitude plot of L7-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 L7-4. # # (a) Replot the amplitudes in L7-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);