Thank you for visiting our site! You landed on this page because you entered a search term similar to this: free answers for comparing linear combination methods for solving systems of equation in a chart. We have an extensive database of resources on free answers for comparing linear combination methods for solving systems of equation in a chart. Below is one of them. If you need further help, please take a look at our software "Algebrator", a software program that can solve any algebra problem you enter!
|
Hey! Here's our bungee jumper! I think she's smiling! A spring-mass system with a constant downward force (gravity) as well. |
In this chapter we return to our mechanical systems and start solving our equations of motion to see how the motion of the system evolves over time. We have introduced all the essential elements of a mechanical system: dampers, springs and masses. We have discussed how to combine their constitutive laws, equilibrium equations and continuity equations to come up with state equations. If we have more than one spring or mass in a system, we will generally have more than one equation of motion and these equations will be coupled together: variables from one equation appear in another. These are therefore "coupled differential equations" that must be solved together simultaneously. Here we will revisit the Euler method to solve differential equations and modify it to handle several equations simultaneously.
State variables vs. parameters
Our larger purpose is to be able to compute how systems move. By "move", we mean that there are some variables that describe the moment-to-moment behavior of the system. The value of these variables is a function of time. These variables are called state variables. If we know the value of these variables at a given moment we know what the system is doing right now, and also, we can tell what it will do one moment into the future.
For mechanical systems, the term "move" is very literal. Among the state variables of a mechanical system are the displacements of components of the system - the length x of a spring for instance. When this spring gets stretched, x changes. It's a function of time: x(t).
For electrical systems, the state variables are also said to "move" because they change with time, even though there is no literal "motion" to be seen. We might better say that the state variables "evolve" with time.
There are lots of other variables that describe a system, besides the state variables. Let's name some of these so that we can keep clear which variables are state variables.
Consider a spring and a mass. The spring has spring constant k, and we've learned that a simple constitutive law states that the spring will exert a force f which is equal to kx, where x is its displacement from its rest length. I like to think of the constitutive law as being a spring's whole purpose in life: to continually look at its own length, and decide what force to exert as a consequence.
The mass, which has mass m kilograms, obeys Newton's second law which says: look at the total force acting on you (f), and decide what rate to change your velocity as a consequence: v' = (1/m) f
Parameters vs. State Varables: m and k are not state variables. They are parameters, simply tokens that we use instead of numbers such as 4.3 kg or 13 N/m. Parameters describe what the system is, while state variables describe what it is doing.
Initial conditions
Consider the spring-mass system above. Suppose that to begin with, the spring is stretched an amount xs=2, and the mass has some velocity, vm =3. By this we mean that at some particular point in time, which we might as well call t=0, we know the values of xs and vm. Recall that the state variables are functions of time. What we know is
xs(t=0) = 2 vm(t=0) = 3
These are called initial conditions. We don't have to be told the values of xs and vm at time zero; any particular instant of time would do. But there is no way for us to predict what will happen next if we don't know what the system is doing right now.
Think about it: if you know that the spring is stretched an amount xs =2, but you don't know the velocity vm of the mass, you can't even tell which direction the mass is moving. It could have positive velocity and be moving to the right, or negative velocity and be moving to the left. A moment into the future, will the spring be more or less stretched? You can't tell.
This sheds some light on state variables. In order to tell you the whole present state of the system, I need to tell you the value of all its state variables at a particular instant. If I give you the value of all the state variables, you can figure out how these values will evolve into the future using the state equations. If I don't give you the values of all the state variables, you can't.
So, the state variables are a minimum set of variables needed to fully describe what the system is doing right now, and with this information we can predict what will happen next using the state equations.
System evolution for Spring-Mass System
In the figure below, on the left side, I've carefully labeled the positive directions of the force acting on the spring. (Just a reminder -- this means the arrow points in the direction we will call positive. It does not mean that at any particular moment the actual force will be in the direction of the arrow. fspring can be negative, after all.)
On the right side, I've drawn the positive direction of velocity Vmand acceleration am for the mass. Remember again that the velocity of the mass is an absolute velocity, while that of the spring (or a damper) is a relative velocity.
Geometric continuity. If the spring is getting longer, the mass moves to the right:
Vm = vs Force balance. Draw the FBD of the mass and impose Newton's laws:
Using the constitutive law of the spring, let's take the numerical example that k=3.5 and m=4
At this point we could identify our two state variables as xs and vm, and go right ahead and reduce the above equations (2-4) into a couple of state equations, which would come out like this
xs' = vm vm' = -k/m xs These are a couple of coupled differential equations, because the equation for xs (5a) involves vm , and the equation for vm (5b) involves xs . They cannot be solved individually. (If you can see exactly why that is, in terms of the Euler algorithm, then you understand diffeqs pretty well.) Note that the two equations (5a) and (5b) are almost exactly what we obtained for our bungee jumper problem, except that in that problem we had an external force, gravity, acting which showed up as a constant in (5b). In that section, we plotted the solution to these coupled differential equations, promising to come back to how to solve them. That's what we're going to do now: we will extend our Forward Euler's Method of solution to a set of coupled differential equations. In the earlier section, we were only solving systems with one state equation.
We could compose the coupled diffeqs into a matrix form by defining a state vector x = (xs, vm)T, and writing the state equations as
where A is a 2x2 matrix (you figure out the details). We won't do that now. Instead, we'll work with the original equations 2-5 numerically for a while. We'll return to matrix form in a bit, since it will be much easier for systems with more than 2 equations.
Forward Euler's Method
The first moment
At t=0 our initial condition (1) tells us xs = 2, so the constitutive law for the spring (4) dictates fs = 7. The spring is in tension with a force applied to it of 7.
Equation (3) tells us that am = -fs/m = (-7)/4 = -1.75.
Finally, equation (2) tells us that since at t=0 Vm = 3 (our initial condition), the spring must be lengthening at a rate vs= 3.
What we just did was apply equations (2)-(4) individually to obtain the time derivatives of our state variables at time=0. We could obtain the same result faster by just applying our state equations (5) directly to the ICs.
A moment later
Just as we did in the unit on differential equations, we are going to solve these equations incrementally and we are going to do so by assuming that both velocity and acceleration are constant in each tiny time step.
Let's look ahead a moment of time dt. In principle this should be an infinitesimal moment. For our numerical example here we'll use a finite moment dt=0.1
If that acceleration of the mass (-1.75) persists for the whole 0.1 seconds, the velocity of the mass will change according to:
Of course the mass did not suddenly get this new velocity at t=0.1. Nor did it have this velocity immediately after t=0. In reality its velocity changed smoothly. Also, we assumed that the acceleration of -1.75 persisted for the whole 0.1 seconds, which isn't quite right either. All of these are approximations, which would be more accurate if we had used a shorter "moment" dt.
If the mass had a velocity of 3 for the whole 0.1 seconds, it will have moved 0.3 to the right. Again, we know it did not have this velocity for the whole 0.1 seconds, but let's leave that aside. If it moved 0.3 to the right, the spring will no longer be stretched xs =2, it will be stretched xs =2.3
So now, at t=0.1, our state variables have evolved from their initial values
xs(0) = 2 vm(0) = 3
to new values (even though these are somewhat approximate)
xs(0.1) = 2.3 vm(0.1) = 2.825
Generalizing algebraically
Let's encode the above repetitive process as algebra, then as differential equations, and also as a Matlab program. We start with initial conditions for the state variables
xs(0) = 2 vm(0) = 3
Now apply our state equations (5)
xs' = vm vm' = -k/m xs to obtain values for the derivatives of the state variables at time equal zero: xs' (0) and vm' (0).
By assuming that that velocity vs = xs' and acceleration am = v'm persist virtually unchanged for a short time step of Dt, using the definitions of velocity and acceleration
v = Dx/Dt and a = Dv/Dt we approximate
xs(t+Dt) = xs(t) + Dt * vs(t) vm (t+Dt) = vm(t) + Dt * am(t) or perhaps, better, rewrite using the state variable derivative directly:
xs(t+Dt) = xs(t) + Dt * xs'(t) vm (t+Dt) = vm(t) + Dt * v 'm(t) This brings us up to the next moment of time, and we can repeat the above process indefinitely
- solve (5) for (xs'(t), v 'm (t))
- solve (8) for (xs(t+Dt)+, vm (t+Dt))
Note that the state equations are specific to our particular system: a mass m and a spring with spring constant k in a particular arrangement. The "update" equations (7) & (8), on the other hand, serve to get us from one moment to the next. They don't depend on anything about the particular devices in our system.
Packing state equations into matrix form
Solving a two-state diffeq numerically
Now let's move to solving a set of coupled diffeq's numerically. First we will put them in matrix form and then solve numerically.
Let's look at another specific example with yet one more state equation. Here we will apply the procedure we just discussed to set up our system of coupled differential equations and update equations for the Forward Euler method numerical solution. Since we will now have three state equations, it will be much nicer to put everything in matrix form. Let's try it!
Consider this system of two springs, a damper and a mass.
- What are the state variables?: x1, x4, v3. So we will have three state equations. To form these, it might help you to go back and check our list of things to do to get our state equations.
For practice, you should go through the details of deriving the state equations. We'll summarize here. Draw FBD's and write force balance:
f1 = f2 
Geometric continuity provides:
The constitutive laws are:
Combine these to obtain our state equations - remember derivatives of the state variables on the left hand side and on the right hand side only state variables, their first derivatives and constants:
Define a 3-vector, named x, which we will call the state-vector. It contains three components, which are (not surprisingly) the three state variables. By define I mean decide in which order the three state variables are to be stored in the state vector x. Write your decision down clearly for the benefit of your readers (or the grader!):
Here, x' is the vector with the derivatives of our state variables. Note that x and x' contain the state variables in the same sequence. Note that I can write these column vectors as a row vector by writing them as transpose:
Note that Matlab uses apostrophe ' for transpose. Don't confuse this with the common mathematical usage (and our usage here) which is ' for time-derivative, T for transpose.
Next, we wish to encode the state equations (1) as a much simpler looking single equation:
where A is a matrix. Here we will do it for this system:
Note that the vectors consist of names of dynamic variables, and the matrix contains nothing but numbers. Some of the numbers are zero, and the others are made up of products or ratios or sums of parameters like k, b, etc, but: parameters are numbers. There can be no dynamic variables in A.
Recalling our numerical procedure, we need the state equations, which we now have in nice form in (3), and the update equations. Note I can write my update equation from before directly in matrix form like this:
This works because the update equation is simply taking an approximation in the definition of the derivative relationship between velocity and displacement and between acceleration and velocity. This same procedure will be used in the electrical and hydraulic domains: our state variables will have special time derivative relationships to one another. We can thus generally write our update equation for any system in vector form as
Thus, numerically we need to first use our initial conditions to find x(0). Then we use (3) to find x'(0). Then we use (4) to find x(0+Dt). So, we repeat:
- use state equations, (3), to find x'(t)
- and update equations, (4), to find x(t+Dt)
Solving our coupled diffeq's numerically
The Matlab program twostate.m is an implementation of an Euler-type diffeq solver, which uses a vector form of the state equations and of the update equation. It works for the two-state spring-mass system that we worked with before:
Study the program, see how (5) is part of it, run it. Notice that the array for state variables is a 2x500 matrix. Each column contains the values of xs and vm at increasing time steps. Conceptually the array for x looks like this:
and at each step in the "for loop" another column of the matrix is filled with values. You will get a chance to modify this program for homework. You should see that if we have a system with three state equations, our big matrix for x in the matlab code would have 3 rows, one for each state variable.
Here is our computational solution for this system with parameters and IC's as before: xs(t=0) = 2, vm(t=0) = 3, k=3.5 and m=4. Here I have taken Dt=0.05:
You can see that the mass and the spring oscillate, which seems intuitively right. You might also notice something strange: the oscillations grow bigger with time. Why is this? Here are the standard four possibilities when you get a surprising numerical result:
- Spring-mass systems really do oscillate bigger and bigger with time.
- Our model of a spring-mass system is in some way wrong
- Our reduction of the model to differential equations is wrong.
- Our numerical integration of the differential equations is wrong.
Real-life experience and conservation of energy rules out #1. #2 is really subtle, and if we can't find the problem elsewhere it is the place we will look. #3 would be just a blunder - we can check that we got the equations correctly. #4 is easy to test, and is a frequent problem with computational solutions. The best way to check our numeric solution is if we have an exact, analytic solution and we will do this in the next section. But here, suppose we do not have an exact solution, so let's check the key assumption in our numerical approximation: Dt was supposed to be infinitesimal, but we used Dt = .05. If our time step being too large is the cause of the behavior we observe, changing Dt will change the amount of that behavior. Below we try Dt = 05:
Much better! So we weren't taking a small enough time step! But it seems that Dt of 0.05 was pretty small already. It turns out that the Forward Euler method is not the most accurate algorithm for solving these differential equations numerically. In the next chapter, we will present a few better algorithms. You will revisit all these in greater depth in EA4.
Finally, note that also a part of problem #4 is that perhaps we just coded something wrong -- our algorithm is robust and correct, but something in the coding is wrong. This is also a frequent problem with numerical solutions and you will need to become proficient at trouble-shooting. Often it helps to step through the first couple iterations of your program by hand and make sure it is doing what you think you told it to do!
Another way to check our numerical solution is to see if we can obtain an analytical solution.
Analytic solution to a coupled differential equation
Coupled differential equations:
dxs/dt = vmdvm/dt = -k/m xs
Let's return to our simple spring-mass system again. We just solved this numerically and initially we saw a solution growing with time that didn't make sense. One thing we should always do if we can to check numerical solutions is find the analytic solution. If we take the first differential equation and plug it into the second, we obtain
This is a second order differential equation. "Order" refers to the highest derivative: differential equations that involve second derivatives are called second order. When taking two first order diffeq's and replacing with a single second order diffeq, note that part of the point is to simplify our life: we want to replace our two equations with two state variables with a single equation and only one state variable. Here, we have chosen xs to be our single state variable and vm no longer appears in the equation. Note that in the end, once we have a solution for xs, we can always solve for vm by differentiation (by one of our first order diffeqs in (1)). Be sure when you move to "second order form" for other systems, that the state variable(s) you eliminate do not appear in your second order diffeq(s). We will not always move to second order form, but it is useful now and will be again in the normal modes section.
As before, equation (2) is a linear differential equation with constant coefficients. It tells us that the second derivative of xs is related to xs itself by a constant. Remember what we said earlier about this class of diffeqs in which the derivatives of the functions are related to the functions themselves: we will always be able to write the solution as a form of exponential. In some cases, the solution can also be written as sines and cosines -- but we will see here in just a minute that these are special forms of exponentials!
So, let's try a solution of the form:
The first and second derivatives of this are
Plugging these into our diffeq (2):
Hmmm - this means that r2 is a negative number, so that r is the square root of a negative number. Here are those complex numbers! To review a few basic concepts on complex numbers, check here.
OK, now let's continue, armed with those facts:
And our solution looks like
Using the Euler formula for the exponential of a complex number, we can rewrite this as
Now, we still have the imaginary j in this solution, which is a little unsatisfying. We just need to know a little bit more about solutions to differential equations to resolve this dilema. Recall that we learned what defines a solution of a differential equation. Now, we need to know about superposition of solutions: a linear combination of solutions to a differential equation is also a solution to the differential equation. We can use this fact to eliminate the j from our solution above. First write the +/- versions of the solution separately and add:
Above, the A1 is simply a new constant (in this case equal to 2A). So, one solution to our differential equation (2) is simply cosine! Similarly, we can multiply (5) through by first +j and then by -j and take + for the +/- term in both cases and then add to obtain:
So, we see that our complex exponential form (4) is equivalent to a general solution in real form like:
Such a solution makes sense: oscillatory motion for our mass/spring system. Note that there is a single frequency for this solution, determined by the values of the mass and the spring constant. This is the natural frequency for the system. Regardless of how the system is set into motion (initial conditions), the response is described by this single natural frequency!
Note that indeed (6) is our general solution: we have satisfied the diffeq, but haven't applied those IC's yet. So (6) is the solution for the spring-mass system subjected to any initial conditions. To find the particular solution for our problem we apply the two initial conditions we were given before: xs(t=0) = 2, vm(t=0) = 3:
So our particular solution is
Let's plot this analytic solution along with our numerical solution from Euler's method:
So you see we could have checked our numeric solution with the analytic solution to see that the oscillations growing in time was incorrect. If I compare this to my solution using the smaller time step, dt=05, we obtain:
Note that I have extended the time scale to show that even though the numerical solution looks pretty good for time up to 20, it is still slowly growing. So, if we were calculating the solution for a long time we would obtain the wrong answer even with this very small time step of 05 (note that we took 10000 time steps just to obtain the plot above!). This reinforces our need for some more accurate numerical algorithms which we cover in the very next chapter.
Before we go there, note a few important things:
- for a single-spring, single-mass system we have found one frequency that describes the oscillatory motion, which is independent of initial conditions. This is the natural frequency of the system.
- we have solved analytically for xs here. We can find vm simply by differentiation, using our first state equation in (1)
- later (in the modes chapter) we will solve systems of many differential equations, leaving each equation individual. That is, we will solve equations like (1) without first substituting one equation into the other.
- If we can obtain an analytic solution, we should do so to check our numerical solution.
- Often an analytic solution to a messy set of diffeq's is not available or would be very difficult to obtain. In many cases, however, you can still find an analytic solution to a simplification or special case of your problem -- and this is still a powerful aid. You can check your code against one or more "limit cases" to verify that it is doing things right in those situations. Given how easy it is to make a slight programming error, remember this tool and use it! Not just in this class, but throughout your engineering career.
Finally, to re-emphasize when we can solve diffeqs analytically with the method we know so far, look at this free-fall example. Some details are left as possible homework problems.
Complex Numbers
General Issues
Complex numbers are numbers that involve the square root of -1, which we will write as j:
Often the square root of -1 is written as i, but since we will use i for current shortly, we will use j for sqrt(-1).
The rules for dealing with complex numbers are simple. Remember these three things:
- Treat j like a variable in an algebraic equation to combine like terms etc, but use the fact that j2=-1 to simplify terms
Here are a few examples:
simple multiplication:
dividing by j: multiply top and bottom by j to get rid of j in the denominator:
In general for complex terms in the denominator, multiply top and bottom by the complex conjugate, again to get rid of j in the denominator:
Euler Formula
When we solve differential equations in general exponential form, we will sometimes obtain solutions that look like:
But what is the exponential of a complex number???
Back in calculus, you have probably seen that the exponential can be expressed as a Taylor series
Since this expansion consists of algebraic terms, we can substitute a complex number in for z to figure out what is going on:
using j2=-1 and combining real and imaginary components:
As indicated, the terms in the brackets can be recognized as series expansions for sine and cosine! Thus, we can write generally that:
This expression for the exponential of a complex number is known as the Euler Formula.
Superposition of Solutions to Linear, Homogeneous Differential Equations
We have been encountering many linear, homogeneous differential equations in this course.
Linear means that the differential equation is linear in terms of the dependent variable and its derivative. For example, 10x'' + sint*x' + exp(t)*x = t2 is a linear equation. x'' + x2=0 is not linear. Also x'' =x/x' is not linear. Homogeneous means the right hand side is zero, after collecting all dependent variable terms (and their derivatives) on the left. For example, x'' + 5x' + 2x = 0 is homogeneous. Also x'' = 12x is homogenous (can rewrite as x'' - 12x =0). However x'' + 2x =t2 is not homogeneous. Likewise x'' + x2 = 5 is not homogeneous.
One extremely important property of linear, homogeneous equations is that of superposition. Before we look at superposition, let's first examine a simple multiplicative constant:
Given a differential equation of the form
(1)
where A is some constant or even function of time, suppose I have a solution x1which satisfies my differential equation (1) and initial conditions x(0)=a, x'(0)=b. Thus I know that :
(2)
Consider any constant times that solution:
Substitute x2 into (1) to see if x2 also satisfies the differential equation:
substituting in for x1 from (2):
Yes, it works: x2 satisfies the differential equation! Note that x2 will not necessarily satisfy the initial conditions, so we will have to find the particular solution for each case after superposition of general solutions. Therefore:
A constant times a general solution to a linear differential equation is also a general solution to the differential equation!
Now, let's look at superposition.
Suppose I have TWO solutions to my differential equation, x1 and x2. Thus I know that
(3)
We can also easily show that the sum of two solutions is a solution. Consider the sum
Substituting into the differential equation (1),
Therefore,
Linear combinations of general solutions to a linear differential equation are also general solutions to the differential equation!
Finally, even though this was derived here for a single linear differential equation of a particular form, we can easily extend the above proofs to any linear differential equation. And indeed to systems of linear equations. For example, consider x"=[A]x , where x is a vector and [A] is a square matrix: if x1 and x2 are known solutions to that differential equation, then x3 = ax1+bx2 is also a solution, where a and b are scalar constants.
Superpostion of solutions is used pervasively in engineering to solve problems in all fields of study. Remember it!
Free-fall with air resistance
|
|
This example has some components that may be solved further as homework. An object falls under gravity. The law for gravity (within a few thousand miles of the Earth) is pretty accurately modeled as
(0)
The constant is mg, where m is the mass of the object (kg) and g is the earth's surface gravity, 9.8 m/sec2.
Newton's second law applied to the mass is v'm = (1/m) SFm
Air resistance is quite complicated, but we'll use a simple model in which the force due to air resistance is proportional to the velocity of the object squared: Fair= A v2. Here v is the speed of the object that is moving through the air, so the velocity that counts (as far as air resistance is concerned) is vm.
Intuitively, what happens:
Starting from rest (this is our initial condition) gravity exerts a (down) force on the object, giving it an increasingly negative velocity. At first its velocity is small, so air resistance does not exert much force in the opposite (up, positive) direction. Therefore the object's velocity becomes increasingly negative and it goes faster and faster. As its velocity increases, the rate at which its velocity is increasing slows, because Fair cancels most of Fgravity. Eventually these two balance, the net force on the mass Fm is zero, so velocity does not change any more. This is called "terminal velocity".
Free body diagram.
The figure above is a free body diagram already. Choose upward to be positive for velocity, force and acceleration. Be consistent!
Finding equations of motion (state equations)
1. What are the states? One x state for each spring and one v state for each mass, so here:
vm
2. What will the state equation(s) look like? There's just one, and it's like this:
v 'm = some function of vm
3. What are the constitutive laws?
Fgravity = mg
Fair= A vm2
4. Force balance equations? There's only one, and we use the FBD of the mass to find it:
S Fm = Fair - Fgravity = m v 'm
5. Geometric continuity equations ?
None!
If I had been really pedantic and said that the constitutive law for air resistance were Fair= A vair2 instead of simply Fair= A vm2, in other words that air speed determines air-resistance, then I would need a geometric continuity equation vair = vm to make the point that the speed of the air rushing by the object is the same as the speed of the mass rushing through the air...
6. Find the state equation
v 'm = (1/m) ( Fair - Fgravity )
v 'm = (1/m) ( A vm2 - mg )
v 'm = (A/m) vm2 - g
And we stop here because this is the proper form for a state equation:
[time-derivative of a state-variable] equals [some function of state variables only].The remainder of this example is left for homework, with some clarifications here.
Computational solution
- Write a Matlab program that follows the evolution of a free-falling object, starting with an initial velocity of zero. g is 9.8 m/sec2. You choose m and A.
Steady state
Note the terminal velocity on the graph that results. We could have found the terminal velocity analytically. Terminal velocity is steady state, meaning that the state variables are not changing anymore.
- If nothing is changing, v 'm = 0. Use this fact and the state equation to find the terminal velocity analytically, for the same parameter values (m, A, g) that you used above. Compare to your computational result.
Analytic solution
Look at the graph. Try to guess a general form for the shape you see. It looks like an exponential decay, v = Ce-Bt, but that would decay to zero unfortunately. So let's guess a spectrum of solutions
v = Ce-Bt - D