Part III. First Order Ordinary Differential  Equations

Section 1. Entering, Solving, Plotting

Unevaluated derivatives are used to enter a differential equation, DSolve solves it. The Table function can generate and plot families of solutions satisfying specified initial conditions.

0. Pull down the Help menu, choose Help Browser... ,Go to Built-in Functions/Algebraic Computation/Equation Solving and read the help page for DSolve.
Open a new Mathematica worksheet and do the following

1. Enter the differential equation y' + y = sin(t). Name it DE. Obtain the  general solution and the solution satisfying y(0) = 0.

2. Continuing 1. Plot the solutions to DE satisfying this initial conditions y(0) = -2, -1, 0, 1, 2. Do it by first creating a list of solutions, then plotting them as shown below.

       solns = Table[ DSolve[ {DE, y[0]==y0}, y[t], t], {y0,-2,2}]
      Plot[ Evaluate[Table[y[t]/.solns[[k]], {k,5}]], {t,0,12}]

3. Obtain the general solution to the equation y ' = t/cos(y)by entering the equation with the name DE and using DSolve.

4. Continuing 3. Obtain the solution to DE satisfying y(0) = 1. Call it soln. Note that there are two solution formulas. Plot both solutions over the interval [-0.2,0.2]. Are they the same? Can you tell from the solution formulas that they are the same? Hint. To make the first plot use the entry

       Plot[ y[t]/.soln[[1]], {t,-0.2,0.2} ]

5. Obtain the general solution to the equation  

y ' (t) = (y^2 + 1)/(t^2 + 1)

Call is soln. Plot the solutions corresponding to C[1] = -2, -1, 0, 1, 2. Hint. Experiment with the horizontal plot range until you get nice picture displaying all 5 curves near t =  0. Use the option Framed->True.

Hint. Mathematica will not permit substitution for C[1]. Get around that by using the following initial plot entry.

      Plot[ Evaluate[Table[ y[t]/.soln/.C[1]->c, {c,-2,2}]], {t,-1,1} ]

6. Continuing 5. Obtain the solution to DE satisfying y(1) = 1. Plot it over a reasonable interval  containing t = 1. Use Table and MatrixForm to generate an array of solution values for t = 0 , 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75, and 2.0. Hint. See Exercise 8 in Part II, Section 3.

7.* Enter the differential equation y' + y = cos(t) with the name DE. Obtain the general solution.

a.   Plot some solutions starting at points evenly spaced on the y axis.   
b.  Use plot some solutions starting at points evenly spaced on the t axis.  
c.   Plot solutions starting at points evenly spaced around the unit circle in the style of the two plots on page 35 of the manual. That is, one picture runs time forward, another runs time backward.

Hint. Use Table to generate a list of solutions, then plot them. See the example below.

In[28]:=

DE = y'[t] + y[t] == Cos[t];
solna = Table[ DSolve[{DE,y[0]==y0},y[t],t], {y0,-2,2}];
Show[ ListPlot[ Table[{0,k},{k,-2,2}], PlotStyle->PointSize[0.02]],
      Plot[ Evaluate[Table[y[t]/.solna[[k]],{k,5}]], {t,0,12}]]

[Graphics:../HTMLFiles/Exercises_18.gif]

8. Look before you leap. Consider the following differential equation y ' (t) = y/(t - 1)   

a.   What is the formula for the function f  such that this equation is equivalent to  y ' (t) = f(t, y)   
b.   Based upon the statement of the Unique Solution Theorem at what points do you expect that solutions will fail to exist and/or fail to be unique?   
c.   Plot the solutions described in  parts a, b, and c of Exercise 7 except replace the unit circle with the circle of radius 0.5.

9. Obtain an informative picture of the family of solutions to  y ' (t) = y/(t^2 - t) over the interval t > 1.

Section 2. Working with Solutions: Modeling

Getting the solution formula is only the beginning of the story in most applications. Mathematical  modeling requires solution manipulation. The exercises in this section are similar to the examples in the  manual.

0. Pull down the Help menu, choose Help Browser... Choose Built-in Functions/Graphics and Sound/2D Plots. Read the information about the Plot, ListPlot, and ParametricPlot functions.
Open a new Mathematica notebook and do the following.

1.  Joe's savings account contains $12,000 dollars. He does not know the interest rate (compounded continuously) but after 60 days he checks and discovers that he has $12,130.    

a.  What is his annual interest rate?  
b. How much money will be in the account one year from the day he has $12,130?   
c. Plot the graph of P(t) including the points corresponding to 60 days and 60 days + 1 year.

Hint. First solve the exponential model initial value problem P' = r P, P(0) = 12000 with variable r and  then substitute the data for t = 60/365 to determine the value of r. If you call the solution "soln" the substitution and r calculation can be made as follows.

In[1]:=

DE = P'[t] == r P[t]
soln = DSolve[ {DE, P[0]==12000}, P, t]

Define g to be the function that solves the equation.

In[3]:=

g = soln[[1,1,2]];
g[t]

Find the r value (numeric).

In[6]:=

Solve[g[60/365]==12130, r]//N

Give that value to the r variable and calculate g(1+60/365)) as follows.

In[7]:=

r = %[[1,1,2]]
g[1+60/365]

The graph should show the solution curve (and the two points).  Don't forget to clear the g and r variables if you want to use them again.

2. Continuing 1. In Exercise 1 you discovered that Joe's savings account has an annual interest rate of  6.555 percent (approximately). Suppose that on the 60th day Joe also has a credit card debt of $560 at a 9% annual interest rate. Starting on that date he decides to pay off the credit card debt continuously from  his savings account at the annual rate of $400 per year.    

a.  How many days later will the credit card debt be paid off?   
b. At the time the credit card debt is paid off, how much money will Joe have in his savings account?

3.* The following problem is adapted from Ledder, Chapter 1, Section 1.  Suppose you borrow $12,000 to buy a car. The loan is to be paid in 60 equal monthly installments at an  interest rate of 5% per year.  

a.  Assume the payments are actually made continuously at whatever rate is needed to pay off the loan in  60 months. Determine the continuous rate per month that would be required.   
b.  Compare the answer to part a to the answer if 60 equal monthly payments are made at a constant  annual interest rate of 5% applied to the outstanding balance. In other words, the first payment, due one month after the loan is made, would be

12000 (e^0.05/12 - 1) + P[1]

     where P[1] is the amount that is put towards reducing the principal (the %12,000) in the first month. The second payment is

(12000 - P[1]) (e^0.05/12 - 1) + P[2]

    where P[2] is the amount that is put towards reducing the principal in the second month. The nth payment is

(12000 - Underoverscript[∑, k = 1, arg3] P[k]) (e^0.05/12 - 1) + P[n]

Hint for b. The second payment is supposed to equal to the first. Therefore, the following equation must  be satisfied.

In[30]:=

eqn = 12000(Exp[0.05/12]-1)+P[1] == (12000-P[1])*(Exp[0.05/12]-1) + P[2]

Out[30]=

50.1043 + P[1] 0.00417536 (12000 - P[1]) + P[2]

This determines P[2] in terms of P[1].

In[38]:=

Solve[ eqn, P[2] ]//Simplify

Out[38]=

{{P[2] 0. + 1.00418 P[1]}}

Create an iterated function that calculates P[3], P[4], ... , P[60] in terms of P[1], then find P[1] using the fact that

Underoverscript[∑, k = 1, arg3] P[k]= 12000

Solution. The anwer to part a is $226.04 per month. The answer to part b is $226.51 per month or about $30 more over the five years of the loan. See the Solutions.

4. Continuing 3. Consult with a banker and/or the internet to determine the amount a bank would charge  per month for a 5% loan of $12000 over 60 months.

5. The following problem is adapted from Ledder, Chapter 1, Section 1. There is a power failure in your house at 1:00 P.M on a winter afternoon and your heating system stops  working. The temperature in your house is 68 degrees F when the power goes out. At 10:00 P.M. the  temperature in the house is down to 57 degrees F. Assume that the outside temperature is 10 degrees F.   

a.    Estimate the temperature in your house at 7:00 A.M. the next morning. Should you worry about your  water pipes freezing?   
b.   Suppose the power goes on at 8:00 A.M. providing heat flow into the house that would increase the  temperature at the rate of 10 degrees per hour if there were no heat loss. At what time will the  temperature in the house be 68 degrees F again assuming the outside temperature stays at 10 degrees F  throughout the day?  
c.   What is the answer to the question in part b if the heating system provides heat sinusoidally with a period of 3 hours and a maximum temperature increase of 10 degrees:
H(t) = 10 - 10 cos((2π t)/3) ?

6. Verify that there is only one value of k that satisfies the equation on page 41 of the manual.
Hint. Plot, as functons of k, the constant function 50 and the expression h(2).

Section 3. Slope Fields

The PlotVectorField function plots vector fields. Normalize the vectors to obtain direction fields (or slope fields). It must be loaded from the Graphics package.

In[9]:=

<<Graphics`PlotField`

Open a Mathematica notebook and do the following.

1. Consider the first order differential equation (t + 1) y ' (t) = 1 - y(t)  taken from Example1 in Section 2.3 of Ledder.

a.   Use DSolve to obtain the general solution and also the solution satisfying y(0) = 0. Plot the second  solution using in the window -2 ≤ t ≤ 2, -2 ≤ y ≤ 2 and compare the picture to the solution curve shown in Figure 2.3.4 in  Ledder and reproduced below.

[Graphics:../HTMLFiles/Exercises_30.gif]

b.   Comment on the existence and uniqueness of solutions using the statement of the Unique Solution  Theorem. Your comments should be based upon the properties of the function  f (t,y) where  y' = f (t,y).

2.* Use PlotVectorField and to make a nice-looking slope field for the autonomous equation y' = sin(y). (Use the  window -6 ≤ t ≤ 6, -6 ≤ y ≤ 6).    

a.   Put some solution curves into the plot using DSolve and Table. Comment  on the relationship between one curve and the next.   
b.   Obtain the general solution formula.
c.   Use
DSolve to obtain the solution satisfying the initial condition y(0) = 1. What is the value of this  solution when t = 1? Get the exact value and an approximation.  
d.   Comment on the long-term behavior of solutions to this differential equation
.

3. Repeat Exercise 3 using the differential equation y ' (t) = 1/sin(y(t)) . Use the same plot window.  Note that the symbolic solutions are much simpler than the ones in Exercise 3. Explain why.

Section 4. Approximate Solutions

Iterated functions and user defined procedures are featured in this set of exercises. 0.
Pull down the Help menu, choose Help Browser... . Then  choose The Mathematica Book/Principals of Mathematica/Modularity and the Naming of Things. Read the pages explaining Modules and local variables.
Open a new Mathematica notebook and do the following.

1. Create the modular function called Euler defined on in Part 3, Section 4 of the manual. Test it on the differential equation y' = cos(t) y with the initial condition y(0) = 1 (as in the manual).

2. Use Euler to obtain a tabular display of Euler approximations to the initial value problem (IVP)

y ' (t) = (8e^(-t))/(3 + y(t)) , y(0) = 1

for t = 0, 0.2, 0.4, 0.6, 0.8, 1.0. (I.e. h = 0.2). Plot these points and the line segments connecting them  along with the actual solution. Name the plot P1. Note. You will have to create the solution using DSolve.

3.  (Continuing 2) Make a similar plot named P2, displaying the solution and the approximation for h = 0.1.  

a.   Use the Show function to display plots P1 and P2 together. Comment on the error displayed in the picture.   
b.    If you have the Ledder text, compare the plot displayed in part b to the one displayed in Figure 2.5.7.  

4.  (Continuing 3) Use your calculator to check that the Error has been cut in half (approximately) when going from h = 0.2 to h = 0.1.

5. Create the modular function called ModEuler and test it on  y' = cos(t) y with the initial condition y(0) = 1 (as  in the manual).

6. Repeat 2 - 5 using ModEuler in place of Euler. Comment on the graphs and the reduction of error in  the Matrix. Hint. Use Copy and Paste, then make minor changes in the code.

7.* Modify the ModEuler procedure to make a modular function that implements the classical Runge-Kutta  algorithm (see Ledder, Section 2.6). Test in by applying it to the problems described in Exercises 2 - 6.  (Copy and Paste). Is the error cut in half in problem 4?


Created by Mathematica  (December 6, 2004)