Numerical Methods In Engineering With Python 3 Solutions Guide

t_test = 2.0 velocity = central_diff(position, t_test) print(f"Velocity at t=2s (central diff): velocity:.2f m/s") distance = simpsons_rule(acceleration, 0, 5, 10) print(f"Distance (integrated): distance:.2f m") 5. Ordinary Differential Equations (ODEs) Euler, Runge–Kutta 4th Order (RK4) def euler(f, y0, t0, tf, h): t = np.arange(t0, tf + h, h) y = np.zeros(len(t)) y[0] = y0 for i in range(len(t)-1): y[i+1] = y[i] + h * f(t[i], y[i]) return t, y def rk4(f, y0, t0, tf, h): t = np.arange(t0, tf + h, h) y = np.zeros(len(t)) y[0] = y0 for i in range(len(t)-1): k1 = f(t[i], y[i]) k2 = f(t[i] + h/2, y[i] + h k1/2) k3 = f(t[i] + h/2, y[i] + h k2/2) k4 = f(t[i] + h, y[i] + h k3) y[i+1] = y[i] + h/6 * (k1 + 2 k2 + 2*k3 + k4) return t, y Example: cooling of an engine block (Newton's law of cooling) def cooling(t, T): T_env = 25 # ambient temp (°C) k = 0.05 # cooling constant return -k * (T - T_env)

slope = (n * sum_xy - sum_x * sum_y) / (n * sum_x2 - sum_x**2) intercept = (sum_y - slope * sum_x) / n return slope, intercept def poly_fit(x, y, degree): coeffs = np.polyfit(x, y, degree) return np.poly1d(coeffs) strain = np.array([0.0, 0.05, 0.10, 0.15, 0.20]) stress = np.array([0.0, 35.2, 68.4, 99.7, 128.5]) Numerical Methods In Engineering With Python 3 Solutions

This guide gives you for typical engineering numerical methods problems. Each block can be extended to full assignments or projects. t_test = 2

t_euler, T_euler = euler(cooling, 100, 0, 60, 2) t_rk4, T_rk4 = rk4(cooling, 100, 0, 60, 2) t_euler, T_euler = euler(cooling, 100, 0, 60, 2)

def beam_ode(x, y): # y = [y, dy/dx, d2y/dx2, d3y/dx3] w = 10.0 EI = 20000.0 dydx = y[1] d2ydx2 = y[2] d3ydx3 = y[3] d4ydx4 = w / EI return [dydx, d2ydx2, d3ydx3, d4ydx4] def shooting_method(): L = 5.0 # Initial conditions at x=0: y=0, d2y/dx2=0 # Guess dy/dx(0) and d3y/dx3(0) from scipy.integrate import solve_ivp # Use secant method to satisfy y(L)=0 and y''(L)=0 # Simplified: for this problem, analytical solution exists. # Numerical approach: def residual(guess): # guess = [dy/dx(0), d3y/dx3(0)] sol = solve_ivp(beam_ode, (0, L), [0, guess[0], 0, guess[1]], t_eval=[L]) return [sol.y[0, -1], sol.y[2, -1]] # y(L) and y''(L)

slope, intercept = lin_regress(strain, stress) print(f"Linear (Young's modulus): slope:.1f MPa")

We solve by converting to 1st-order system.