tutor-notes-2005-sheet-5

tutor-notes-2005-sheet-5 - Nafz 5 ‘ICV‘ fifbrs The...

Info icon This preview shows pages 1–12. Sign up to view the full content.

View Full Document Right Arrow Icon
Image of page 1

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 2
Image of page 3

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 4
Image of page 5

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 6
Image of page 7

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 8
Image of page 9

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 10
Image of page 11

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full Document Right Arrow Icon
Image of page 12
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: . . . Nafz 5 ‘ICV‘ fifbrs The UmverSIty of Queensland School of Engineering MECH2700 Engineering Analysis I (2005) Exercise Sheet 5 1. As part of the Thermodynamics course, there is an experiment on “The Performance of a Diesel Engine”. Pressure and crank angle are recorded for a single cycle of the engine and are written to a file with three values to a line. The first value is the sample number and is unimportant. The second value is the crank-angle signal and the third value is the voltage output from the pressure transducer. Write a Python program that takes the LabVieW output file (outputlab) containing sampled values for the crank angle and pressure transducer signalsand computes (i) the crank angle, 6, in radians, (ii) the cylinder volume, V, in m3 and (iii) the pressure, P, in kPa. Write the transformed data to a file called “pv.dat” and plot both the pressure versus crank angle and a P-V diagram. A second phase of your program should integrate the PV data to obtain the work done by the gas during the compression stroke (7r < 0 g 27r) and the expansion (power) stroke (271' < 0 g 37r). To integrate the sampled data, use W = f P dV a 2 Eu), + P--1)(V; — V_1) Where subscript 2' refers to the current sample point and i —- 1 refers to the previous sample point. k: Top Dead Centre M C I F: Bottom Dead Centre I I ‘ \ h2 , r T \ i l \ | I ‘ -i- ’ The pressure transducer responds linearly with a sensitivity of 1185 kPa/ Volt. The crank-angle signal maps linearly from the range ——5 S s S 15 to crank angle in radians ~7r g 0 S 371‘. The volume inside the cylinder is given as VZl/1+A[T+O—h1—h2] Where T is crank angle h1 = 1" cos T h2 2 V02 — (rsinT)2 L is the stroke, 0.075 In D is the cylinder bore, 0.080 In C is the connecting—rod length, 0.133 m A is the area of the cylinder cross-section, 7rD2/4 r is the crank radius, L/2 V2 —— V1 is the swept volume, A x L Vfi / V1 is the compression ratio, 21.2:1 . Recall the first minor assignment on the flight of a projectile. Rewrite the equations of motion for the projectile as a set of first-order ordinary differential equations 5": f(2', t) using the state vector ’03! Leave the drag as the symbol CD in these equations. . Implement an ODE integrator in Python using the Euler method (Schilling & Harris Section 8.2) and another using the fourth-order Runge—Kutta method (Schilling 85 Harris Section 8.3). Use both of these functions to integrate the projectile equations of motion for minor assignment 1. Use CD = 1.0 and select a suitable time step to get an accurate estimation of the trajectory. Compare the number of steps needed for each method to provide a value for the distance of flight to an accuracy of 0.1%. # file: diesell.py “"“Program to convert from the Labview format to pv quantities. PJ, 18—Nov—2003 "'1" print “Begin diesel...“ from math import pi, sqrt, sin, cos 0.075; # stroke in metres L / 2; # crank radius, m 0.080; # piston diameter, m 0.133; # con rod length, m pi * D * D / 4; # piston area sweptVolume = A * L v1 = sweptVolume / 20.2; psensitivity = 1185.0; # kPa/volt W<1tiw H H "I!” " print “Read raw data“ fin = open("output.lab“, "r") labview_content = fin.readlines() fin.close() print “lines read=“, 1en(labview_content) print "phase 1: select and rescale the data" fout = open(“pv.dat", "w") print >> fout, “# theta (rad), volume (m**3), pressure (kPa)" tlist = [1; vlist = ; plist for line in labview_content: strl, strz, str3 = line.split() cranksignal = float(str2) voltage = float(str3) if 5.0 <= cranksignal <= 15.0: thet = cranksignal / 5.0 * pi t1 = R * sin(theta) t2 sqrt(C * C — t1 * t1) t3 R + C — R * cos(theta) — t2 volume = v1 + t3 * A pressure = voltage * psensitivity print >> fout, "%10.4e %10.4e %10.4e“ % (theta, volume, pressure) tlist.append(theta); vlist.append(volume); plist.append(pressure) fout.close() print “Points saved=“, len(tlist) I! I] print "phase 2: compute the work done in each stroke" compress_count 0 expand;count = 0 compress_work = 0.0 expand_work = 0.0 for i in range(len(tlist)): t_i = tlist[i]; v_i = vlist[i]; p_i = plist[i] if pi <= t_i <= 2.0*pi: compress_count += 1 if compress_count >= 2: compress_work += 0.5 * (p_i + p_im1) * (v‘i — v_im1) if 2.0*pi <= t_i <= 3.0*pi: expand~count += 1 if expand_count >= 2: expand_work += 0.5 * (p_i + p_im1) * (v_i — v_im1) v_im1 = v_i; p_im1 = p_i print "Compression stroke: points=", compress_count, “ Work=", compress_work, "kJ" print "Expansion stroke: points=“, expandgcount, "Work=“, expand_work, “kJ“ Begin diesel . . . Read raw data lines read= 727 phase 1: select and rescale the data Points saved= 364 phase 2: compute the work done in each stroke Compression stroke: points= 165 Work= —O.103513217933 kJ Expansion stroke: points= 199 Work: 0.72127774209 kJ Pressure—Volume Diagram from Diesel Experiment 7000 6000 5000 4000 3000 Pressure (kPa) 2000 1000 0 0.0001 0.0002 0.0003 Volume (m**3) 0.0004 # file: diesel.gnuplot set term postscript eps 20 set output “diesel_raw_data.eps“ set title "Raw Data from Diesel Experiment" set xlabel “Angle signal“ set ylabel "Pressure signal (volts)" plot “output.lab" using 2:3 notitle set term postscript eps 20 set output “diesel_pcrank.eps“ set title "Pressure—Crank-Angle Diagram from Diesel Experiment“ set xlabel “Crank angle (radians)" set ylabel “Pressure (kPa)“ plot "pv.dat“ using 1:3 notitle with lines set term postscript eps 20 set output “diesel_pv.eps“ set title “Pressure—Volume Diagram from Diesel Experiment“ set xlabel "Volume (m**3)“ set ylabel “Pressure (kPa)“ set xtic 0.0001 plot "pv.dat" using 2:3 notitle with lines # be careful after this point because xtic has been set Pressure signal (volts) Pressure (kPa) Raw Data from Diesel Experiment 1- ++ + + + + + + + + + + + + + E 4+++++ ++W+ —6—4—2 0 2 4 6 810121416 Anglesignal Pressure-Crank-Angle Diagram from Diesel Experiment 7000 6000 5000 4000 3000 2000 1 000 3 4 5 6 7 8 9 10 Crank angle (radians) mam 6W ’5 CD2 ‘GoHLchZL brajzdong. ' ]/ Afiow t Trail/4 ,‘ A? ‘ ’U’x": 4‘ 093 A _ ‘ A 2 ‘6 Pd 2PM A6“ (7,5,) I mayhde GIVE/dun GCCeQWw; ax :- ~ 123) PM W m aw ”(92.) LE1 — 9 10] m SW 2: I; mvodeswcwm g2 dL/de _ 0x \adw 49} t dfi/oxe v, 0“; 0”th 0w (”9/016 ‘9 # file: golf_hall_pde.py """Compute and plot the trajectory of a golf ball. This is a no—frills reference solution for ndnor assignment 1 for mech2700, 2004. P. Jacobs 29—Sep—2003 10—May—2004 Adapted to Sheet 5 Qn,3 import math from ode_basic import rk4_update g 9.81 # acceleration at ground level m. 0.046 # mass of ball, kg d = 0.043 # diameter of ball, m Afr nt = math.pi * 0.25 * d * d rho = 1.2 # air density, kg/mA3 def ODEs(t, 2): """Given t and state 2, returns derivative dzdt.""" global Cd # Unpack state vector XIYIVXIVY=Z vsqr = vx * vx + vy * vy v = math.sqrt(vsqr) Fd = 0.5 * rho * vsqr * Afront * Cd ax = —(vx / v) * Fd / m ay = -(vy / v) * Fd / m — g return [vx, vy, ax, ay] if __name__ == “__main__": print "Begin golf—ball simulation..." 0.0 0.0 v0 — 48.0 # m/s alph = 21.0 / 180.0 * math.pi va = math.cos(alpha) * v0 vyO = math.sin(alpha) * V0 to = 0.0 - e 0 ill" print "Initial state:", x0, y0, va, vyO # Let’s start to record the motion in some lists, # one for each state variable. ts = [t0]; xs = [x0]; ys = [yO]; vxs = [va]; vys = [vyO] dt = float( raw_input("Enter dt: ') ) nstep = int( raw_input("Maximnm number of steps: ") ) Cd = float( raw_input("Cd: “) ) print "dt=“, dt, " nstep=", nstep, " Cd=“, Cd # This loop controls the time—stepping. for i in range(nstep): (x, y, vx, vy) = rk4_update(0DEs, t0, [x0, y0, va, vyo], dt) t = to + dt print "t=", t, x, y, vx, vy if y <= 0.0: # Stop the loop if the projectile has returned to the ground. # Linearly interpolate to get the final state. fraction = (0.0 — yO)/(y — yO) tf — to + fraction * (t — t0) xf = x0 + fraction * (x — x0) yf = yo + fraction * (y — YO) vxf = va + fraction * (vx — VXO) vyf = vyO + fraction * (vy — vyO) # Record the final state ts.append(tf) xs.append(xf); ys.append(yf); vxs.append(vxf); vys.append(vyf) print “Final t=", tf, xf, yf, vxf, vyf break # We haven't stopped; record the new state. ts.append(t) xs.append(x); ys.append(y); vxs.append(vx); vys.append(vy) # Save the new state as the old state for the next step. to, x0, yo, va, vyO = t, x, y, vx, vy # Present the trajectory graphically. import Gnuplot plt = Gnuplot.Gnuplot() d1 = Gnuplot.Data(xs, ys, with="linespoints“) plt.title("Golf—ball trajectory with Cd=%g, dt=%g" % (Cd, dt)) plt.xlabel("x, metres") plt.ylabel("y. metres") plt.plot(dl) junk = raw_input(“Press return..,") plt.hardcopy("golfvball—ode—trajectory.eps", terminal="postscript“ mode="eps", fontsize=20) junk = raw_input("Press return again...") print “Done.“ I # ode_basic.py # Example solution to sheet 5, qn 3. # P.J. 10—May—04 from.copy import copy def euler_update(f, tk, xk, h): def "“"State vector update as per equation 8.2.4 -_ Euler method. Input: b‘géfi‘l—h nu" .- step size in t Returns a list of new values for the state vector x_kpl. assert callable(f) assert type(xk) == list assert len(xk) > O dzdt = f(tk, xk) x_kp1 = copY(xk) for i in range(len(xk)): function that accepts the current state and returns its derivative. independent variable current state (list of state~vector elements) x‘kpl[i] = xk[i] + h * dzdt[i] return x‘kpl Other documentation as per euler_update(). "I‘ll assert callable(f) assert type(xk) == list assert len(xk) > O indx = range(len(xk)) ql = f(tk, xk) xkl = copy(xk) for i in indx: xklIi] += h * q1[i] / 2-0 q2 = f(tk+h/2.0, xkl) xk2 = copY(Xk) for i in indx: xk2[i] += h * q2[i] / 2.0 q3 = f(tk+h/2.0, xk2) xk3 = copy(xk) rki4_update(f, tk, xk, h): """State vector update as per equation 8.3.4 —— 4th—order Runge—Kutta. # we will reuse this list of indices for i in indx: Xk3[i] += h * q3[i] q4 = f(tk+h, xk3) x_kpl = copy(xk) for i in indx: xhkpl[i] += h /6.0 * return x_kp1 def prey(t, x): """Predator~prey model as per Schiller and Harris ex 8.3.2. Returns a vector of derivative values given time t and state vector x."““ assert len(x == 2 dedt dxldt return [dedt, dxldt] __name__ == "__main ": 20 = [0.5, 0.5] t=030 h = 0.04 N = 250 # methodName = ’euler' methodName = ’rk4’ xOlist = [1; xllist =; for k in rangetN): (ql[i] + 2 * q2[i] + 2 * q3[i] + q4[i]) (1.0 — x[l]) * x[0] (—1.0 + x[0]) * x[l] tlist K euler stepper, h= 0.040 : predator~prey model Example 8.3.2 2.4 if methodName == ’euler’: z_new = euler_ppdate(prey, t, 20, h) else: z_new = rkf4_update(prey, t, 20, h) t=t+h tlist.append(t) XOList.append(z_new[O]) xllist.append(z_new[l]) 20 = copy(z_new) import Gnuplot g1 = Gnuplot.Gnuplot() dl = Gnuplot.Data(XOlist, xllist, with=’lines’) gl.xlabel('x0’) @1.ylabel(’xl’) g1.title('Example 8.3.2: predator-prey model: %3 stepper, (methodName, h)) gl.plot(d1) gl.hardcopy("example832.eps", terminal="postscript", \ mode=‘eps“, fontsize=20) “'2 CH O V‘ 9. C3 H N I Q‘ m C1 8‘ E73 ‘¢ Ln '5 '— '5 13 (3 c; E >5 9 9- ,_ h *3 523 CL ’01 09 U3 . C5 GD 32 CL E (U X ”J (3 CH 1— v- T— r- C: CD CD CD C“ T. LX =96.3f’ 16 1.4 1.2 c T) Lx 0.8 0.6 0.4 0.2 2.2 1.2 1.4 1.6 1.8 2 x0 1 0.2 0.4 0.6 0.8 0.3 60.148 025 $0.210 0.20 €0.20?) TM Ah ”a“; OK 0.15 Lo‘zja 040 (90.280 0‘05 (90.288 001 @0231 «maul; Hwy: «:04 2 0.001 (9029’) {(90231 —— 60.351] Golf-ball trajectory with Cd=1, dt=0.2 y. metres [0 (D -f-‘- 01 03 \I co to 0—4. X, metres ...
View Full Document

{[ snackBarMessage ]}

What students are saying

  • Left Quote Icon

    As a current student on this bumpy collegiate pathway, I stumbled upon Course Hero, where I can find study resources for nearly all my courses, get online help from tutors 24/7, and even share my old projects, papers, and lecture notes with other students.

    Student Picture

    Kiran Temple University Fox School of Business ‘17, Course Hero Intern

  • Left Quote Icon

    I cannot even describe how much Course Hero helped me this summer. It’s truly become something I can always rely on and help me. In the end, I was not only able to survive summer classes, but I was able to thrive thanks to Course Hero.

    Student Picture

    Dana University of Pennsylvania ‘17, Course Hero Intern

  • Left Quote Icon

    The ability to access any university’s resources through Course Hero proved invaluable in my case. I was behind on Tulane coursework and actually used UCLA’s materials to help me move forward and get everything together on time.

    Student Picture

    Jill Tulane University ‘16, Course Hero Intern