From 380f7244e82cc5798519956714b9bd1720cdbf6c Mon Sep 17 00:00:00 2001 From: Joe Anderson Date: Mon, 31 Oct 2011 07:59:21 -0500 Subject: Angular Momentum graphs --- throw.py | 51 +++++++++++++++++++++++++++++++++++++-------------- 1 file changed, 37 insertions(+), 14 deletions(-) diff --git a/throw.py b/throw.py index deaea5b..3fcfdd8 100644 --- a/throw.py +++ b/throw.py @@ -66,13 +66,13 @@ sin = taylor_sin ### radii: global r0 -r0 = 1. #meters -- the radius of the first full-circle rotation +r0 = 1. #meters -- the length of cable in the first full-circle rotation r_anchor = .2 #meters -- the radius of the anchor point r_winch = .2 #meters -- the max radius of the roll of cable length = 1000 #meters -- the length of the cable phi_max = 2 * (length - r0) / r_winch #rad -#r = lambda phi: r0 + r_winch*phi # the radius of the satellite -r = lambda phi: r0 + r_winch*(phi - phi**2/(2*phi_max)) +r = lambda phi: r0 + r_winch*(1 - phi/(2*phi_max))*phi # cable let out +r_sat = lambda phi: numpy.sqrt(r(phi)**2 + r_anchor) # radius of satellite ### gravitational acceleration: geopotential = 3.987 * 10**14 #cubic meters per second^2 @@ -93,10 +93,12 @@ quad_co = kappa * density * Area global m, Mass m = 39.7 #kilograms -- mass of the satellite cage (overestimate) l = 0.16 #kilograms/meter -- density of Super Max cable, 16mm: 269.8kN -max_tension = 269.8 #kN, given by choice of cable above. +l = 0.41 #kilograms/meter -- density of Super Max cable, 26mm: 647.4kN +max_tension = 647.4 #kN, given by choice of cable above. Mass = 670 #kilograms -- mass of /everything/ M = lambda phi: Mass - m - l*r(phi) -v_cm_co = lambda phi: (.5 * l*r(phi) + M(phi)) / (m + l*r(phi) + M(phi)) +r_cm = lambda phi: (m*r(phi) + .5*l*r(phi)**2) / (m + l*r(phi) + M(phi)) +v_cm_co = lambda phi: 1 - r_cm(phi)/r(phi) tan = lambda phi: r_anchor / r(phi) # tan(theta) def throw(phi, time): @@ -124,21 +126,25 @@ def interpret(soln, max_tension=max_tension): phi = soln[:, 0] phi1 = soln[:, 1] v_tan = v_cm_co(phi) * phi1 * r(phi) / 1000 #kilometers/second - cos_t = lambda phi: numpy.cos(numpy.arctan(tan(phi))) # cos(theta) + cos_t = lambda phi: cos(numpy.arctan(tan(phi))) # cos(theta) tension = (m/cos_t(phi))*(phi1**2 * r(phi) - g(tau/4))/1000 #kiloNewtons + cable_moment = l*((r(phi) - r_cm(phi))**3 - r_cm(phi)**3)/3 #kg m^2 + moment = m*(r(phi) - r_cm(phi))**2 + cable_moment - M(phi)*r_cm(phi)**2 + momentum = moment*phi1**2 #kg m^2/s^2 try: maxindex = numpy.argmax(tension >= max_tension) max_arg = numpy.argmax(v_tan[:maxindex]) except ValueError: max_arg = numpy.argmax(v_tan) - return phi, phi1, v_tan, cos_t, tension, max_arg + return phi, phi1, v_tan, cos_t, tension, moment, momentum, max_arg -def plot(time, soln, graphs=['all'], i=5): - phi, phi1, v_tan, cos_t, tension, max_arg = interpret(soln) - phi_max, phi1_max, v_tan_max, tension_max = \ +def plot(time, soln, graphs=['all'], i=7): + phi, phi1, v_tan, cos_t, tension, moment, momentum, max_arg = \ + interpret(soln) + phi_max, phi1_max, v_tan_max, tension_max, moment_max, momentum_max = \ map(lambda x: numpy.ndarray.__getitem__(x, max_arg), - [phi, phi1, v_tan, tension]) + [phi, phi1, v_tan, tension, moment, momentum]) def test(*args): b = 0 @@ -197,6 +203,22 @@ def plot(time, soln, graphs=['all'], i=5): pyplot.axvline(x=phi_max/tau) + i -= 1 + if test('moment'): + pyplot.figure(i) + pyplot.plot(phi / tau, moment) + pyplot.grid(True) + pyplot.ylabel("Moment of Inertia (kg m^2)") + pyplot.xlabel("phi (cycles)") + + i -= 1 + if test('momentum'): + pyplot.figure(i) + pyplot.plot(phi / tau, momentum) + pyplot.grid(True) + pyplot.ylabel("Angular Momentum (kg m^2 s^-2)") + pyplot.xlabel("phi (cycles)") + #import pdb; pdb.set_trace() def opt_fun(guesses, f, func, init_vars): @@ -206,10 +228,11 @@ def opt_fun(guesses, f, func, init_vars): time, soln = solve_fun(func, init_vars) - phi, phi1, v_tan, cos_t, tension, max_arg = interpret(soln) - phi_max, phi1_max, v_tan_max, tension_max = \ + phi, phi1, v_tan, cos_t, tension, moment, momentum, max_arg = \ + interpret(soln) + phi_max, phi1_max, v_tan_max, tension_max, moment_max, momentum_max = \ map(lambda x: numpy.ndarray.__getitem__(x, max_arg), - [phi, phi1, v_tan, tension]) + [phi, phi1, v_tan, tension, moment, momentum]) f.write("%s, " % r0) f.write("%s, " % m) -- cgit v1.2.3-54-g00ecf