## Calculus Without Tears## Synopsis of Mathematical Modeling and Computational Calculus I## IntroductionThis book will take you from not being able to spell ‘calculus’ to doing calculus just the way I did it for twenty years as an engineer at high tech firms like Lockheed and Stanford Telecom. You will learn how physical processes are modeled using mathematics and analyzed using computational calculus. Systems studied include satellite orbits, the orbits of the earth and moon, rocket trajectories, the Apollo mission trajectory, the Juno space probe, electrical circuits, oscillators, filters, tennis serves, springs, friction, automobile suspension systems, lift and drag, airplane dynamics, and robot control. And not a single theorem in sight. This book focuses on differential equation models because they are what scientists and engineers use to model processes involving change. Historically, this has presented a problem for science education because while creating a model is usually straightforward, solving the model differential equations is usually impossible, that is, most differential equations do not have analytic solutions. This has been a tremendous roadblock to science and science education. But, that was before computers; now, with computers, solutions to differential equations can be computed directly, without the need for an analytic solution, and without advanced mathematics, using the formula ## Constant Velocity MotionWe start from ground zero and cover calculus as applies to constant velocity motion. Everything there is to know about constant velocity motion is in the formula Constant velocity motion is represented mathematically by a function of form The derivative of a function is the rate of change of the function, the derivative of Suppose we have a position function where we've written the function name in bold letters to remind us that its definition is unknown. An equation containing the derivative of an unknown function is called a differential equation. To solve this differential equation we need to find a definition of p(t) that makes the differential equation true. Can you solve this differential equation? Hint: see the previous paragraph.## Mathematical Models of a Falling ObjectAccording to the story he told himself, Newton discovered the Law of Universal Gravity after watching an apple fall from a tree on his country estate. Gravity is the force that pulls the apple toward the earth and causes it to fall. It also acts on the moon to hold it in its orbit. The formula, it could have been a lucky guess, is: F = G * M _{earth} * M_{apple} / (R * R), where
G is the gravitational constant M _{earth} is the mass of the earth
M _{apple} is the mass of the apple
R is the distance between the centers of earth and apple Having discovered gravity, Newton understood that there was a force acting on the apple; now he wanted to understand the apple’s trajectory as it falls, that is, he wanted to know the position function for a falling apple. First, he needed to know how the force of gravity affects the apple's trajectory, and to this end he discovered the Second Law of Motion (after discovering the First Law of Motion), which is force (F) equals mass (M) times acceleration (A), or
Putting the two formulas together gives us Newton's model for a falling apple, that is,
v(t) = -10 * t + V_{0}. That is, if we drop the ball at time t = 0 with initial velocity 0, then after 1 second its velocity is -10 m/sec, after 2 seconds its velocity is -20 m/sec, and so on. ## Graphing Solutions to Differential EquationsHow can we plot the trajectory of a falling apple, using our simple model _{0}, with V_{0} = 0 ? (Note: computing the solution for Newton's model involves just one more computation, but we're going for total simplicity here.) As the apple falls its downward velocity increases, so we're not dealing with our favorite type of motion, constant velocity motion. Can we convert the problem into a constant velocity problem? We can. Suppose we want to compute the trajectory for 10 seconds, what we can do is divide the time interval into shorter subintervals, and assume that the velocity is constant on each subinterval ! This is how we do it.We start by dividing the 10 second interval into two 5 second subintervals. We need a velocity for the 1 ^{nd} subinterval? The apple's velocity at the start of the subinterval is ’(5) = -10 * 5 = -50, so that is the velocity we use for the second subinterval. So, we have the approximation to the trajectory of a falling apple shown. Not so great, but, we're just getting started.
pThe problem with our first approximation is that the subinterval size was too large. If we divide the 10 seconds into 10 one second intervals, and follow the method used above to produce a linear approximation for each subinterval, we get the graph to the left. This is a pretty good approximation to the true trajectory of a falling apple.
So, let's apply it to this spring-box assembly sliding on a frictionless table. If the spring is at rest when (t) / m, and we have our model.pIt's a second order model which we can rewrite as two first order models, (t) and v'(t) = -k*v(t) / m and we're ready to compute a
solution.
pOops, because we didn't want to do too many computations by hand we chose a large subinterval size, and it was too big, and the box has gotten away from us ! We should have picked a smaller subinterval size. ## An Ounce of TrigonometryWe will be modeling 2-dimensional processes, e.g. orbits (a falling object in 2-d), and that requires a little bit of trigonometry to resolve vectors into their x and y components, essentially the trig necessary to convert from rectangular to polar coordinates and back, that is (x, y) = (R, ## Introduction to MATLAB/OCTAVE/FREEMATMATLAB is a user friendly computer programming language designed for ease of use and applicability to engineering problems. We do not use any advanced features of MATLAB, only the commands for adding, subtracting, etc., and saving the results, an occasional if statement, for loops, and commands for plotting. Given a model, we might do the necessary calculations 3 times by hand so that they are understood completely, and then use a MATLAB program to iterate the calculations 10,000 times. The program only performs operations that we have already done by hand, plus plotting. The program the computes the 1-second subinterval trajectory above can be seen here (you thought I really did it by hand?) The downside of MATLAB is that it is expensive. Fortunately two completely free and legal clones are available, FREEMAT and OCTAVE. The clones execute exactly the same commands as MATLAB and are for our purposes equivalent ## Orbits, Satellites, RocketsThe moon, and satellites around the earth, are falling bodies and hence the mathematical model for their trajectories is the same as the model for a falling apple, with one difference, that is the simplified model won't work, it takes Newton's model to get accurate elliptical orbits. Using that model with the proper initial conditions produces the graph of the moon's orbit shown here: Note: the illustration is not to scale, the graph is. The program the computes the moon's orbit (minus the outline of earth) can be seen here A geostationary satellite orbit has a period of 24 hours so that the satellite in orbit above the equator rotates at the same rate as the earth and appears fixed in the sky when viewed from the earth. The same model works, initialized to the initial state of a geosynchronous satellite orbit: Note: the illustration is not to scale, the graph is.
Adding a rocket to the model of the moon's orbit, we can launch (in this case shoot) the rocket from the earth and aim for the moon. A hit ! The Apollo mission did not hit the moon but rather just missed, and the lunar module went into orbit around the moon. We hit the moon in the last project, we can also miss, but no matter how closely we come to the moon, we can't get into orbit around the moon without a guidance boost, so we program a directed boost just as the rocket passes the moon, et viola'. The plot shows the Apollo spacecraft going into orbit around the moon and making several orbits, the outline of the moon is drawn at its position at the end of the program run. The US launched a probe to Jupiter on Aug. 5, 2011. The trip to Jupiter will take five years, and utilizes a
slingshot, or gravity assist, trajectory. The boost rocket put a Centaur upper stage and the Jupiter probe into a low
orbit around the earth. The Centaur then launched the probe into a deep space orbit around the sun that will return to the earth in two years! The gravitational field of the earth will then accelerate the probe and send it on its way to Jupiter. The flight of the probe is un-powered except for two deep space maneuvers a year into the flight. A diagram of the actual planned trajectory is shown in the figure.
Once the probe detaches from the Centaur booster it is acted on only by the gravitational fields of the earth and the sun (minus the deep space maneuvers). The model is the falling object model times three, the earth is falling towards the sun, the rocket is falling towards the sun, and the rocket is also falling towards the earth. We program these three models and add a programmed deep space guidance boost to get the right trajectory, and viola'. The graph shows the 3 year trajectories of the earth around the sun, which traces out a circle 3 times, and the rocket, the big ellipse followed by the slingshot trajectory to the outer planets.
## Electrical CircuitsThe models for the primary circuit elements are We can use graph the capacitor voltage when the switch closes as shown
An RLC oscillator circuit is shown, it includes a voltage generator V _{S}(t) - (t)*R - iv_{C}(t))/L
and the capacitor model gives us the differential equation for v _{C}(t),
v_{C}'(t) = (t)/C
iWe can 'program' the voltage source by first initializing it to zero VS = zeros(1,N), and then setting non-zero values as we like. Suppose the subinterval length is 0.01 seconds and we want a unit (integral) impulse to occur at 1 second into the simulation, then VS(100) = 100 does it. The graph shows the V _{S} voltage.
This graph shows the response of the system to a unit impulse at 1 second, see the program here. ## DynamicsAn actual 3-d object not only moves from place to place, it rotates, think of a kick serve in tennis, the ball spins as it travels through air. The motion of a 2-d object (the bookkeeping for a 3-d object is considerably more complicated) can be decomposed into the motion of the center of mass, and the rotation of the object about the center of mass. Newton's 2 We revisit the spring-box assembly, this time with an added hydraulic damper, similar to an automobile shock absorber. The photo shows a spring and damper combination. The force exerted by the damper is proportional to the velocity of the plunger in the damper and is C*p'(t) where c is the
damping coefficient. The model of the spring-damper assembly is (t) - c*p'(t). Let's jump ahead a bit and put the spring-damper in a tennis ball that bounces and skids.pThen let's have Maria Sharapova give the ball a good whack, now the ball is traveling
through the air and spinning. When the ball hits the ground it will skid and bounce; the skid is modeled as sliding friction, also called Couloumb friction, the model is simply F = -u*N*sign(v) where u is the coefficient of friction, N is the force (usually the weight) of the object on the surface, and sign(v) is 1 for v positive, 0 for v = 0, and -1 for v negative. Thus the magnitude of sliding friction is constant and it is in the direction opposite the apparent velocity of the ball on the surface. The skid generates a force which accelerates or decelerates the ball and which creates a torque which changes the spin rate. Putting all the models together we can plot the trajectories for kick (topspin), flat, and slice (backspin) serves. This is all spelled out in detail in the book. The spring-damper combination is used in motorcycle (shown) and auto suspension systems. We have one axel (shown) and two axel models. For the one axel model the passenger just goes up and down. For the two axel model there is one spring-damper assembly of the back axel and one on the front,
and the auto chassis is a rigid object acted on by gravity plus the forces of the two spring/shock absorbers, so it rocks. We can plot the up/down motion of the car when it hits a bump in the road. We program the bump just the way we programmed the impulse to the high-pass filter circuit above. The graph shows the bump in the road and the cars response to it. See the program here. The forces acting on an airplane in flight are thrust, lift, drag, and gravity. In addition the pilot can activate control surfaces that generate force on the airplane. In order to study aircraft flight we need models these forces. Gravity we know, we will program thrust directly, and we will program one control torque, the torque resulting when an tail elevator is activated, directly. That leaves lift and drag. Lift and drag of an airplane depend on the airplane's angle of attack (AoA), that is the angle between the centerline of the airplane and the airplane's velocity vector, as shown. The relationship between the lift coefficient (CL) and drag coefficient (CD) for a typical aircraft are shown in the figure. We model these functions with ad hoc piecewise linear functions and lift is determined using the formula, and a similar formula used for drag,
We initialize the mass, moment of inertia, and wing surface area from the data for a single-engine airplane, and we're ready to fly. Once we find the thrust that keeps us in level flight we have only one control to activate, so we raise the tail elevator causing the tail to go down, the plane to go up, and then lose velocity and lift, and then go down as a result, as we see in the graph ## Appendices## App. I - Installing Octave/Freemat## App. II - Computational Calculus vs. Analytic CalculusClick on 'New Computatinal Calculus Vs. Old Analytic Calculus' to the left.## App. III - Other Approximation MethodsThe Runge-Kutta methods## App. IV - Convergence of Euler's MethodHey, it works !## App. V - Mathematical Models and Physical LawsSummary of models and laws used, click here to see App. V## App. VI - Debugging## MMCC I Answer Key - here |