Cs 312 Project 1

Project 1: Scientific Computing Posted: September 4th, 2008. Due: Thursday, September 25th, 2008 Objectives: 1. This project should introduce you to the idea of Symbolic vs. Numerical solvers, and help you understand the difference between “analytic” vs. “algorithmic” solutions to problems. 2. You should understand the importance of difference and differential equations for modeling various phenomena, including the efficiency of some algorithms. 3. You should understand how “numerical methods” use difference equations to approximate the solutions of differential equations, and become familiar with the concept of stability. 4. You should be able to solve any second order linear, constant coefficient difference equation with geometric forcing terms by hand (e.g. on a test). Background: Functions that map the integers (usually the non-negative integers) into the real or complex numbers are frequently used to represent a variety of phenomenon. These phenomenon include the number of packets received by a particular router each minute, the sizes of a sequence of email messages, the closing price of various securities traded on the New York Stock Exchange, the number of points scored by a particular sports team over a sequence of games, or even the efficiency of a particular algorithm. Describing such functions can be difficult. One method would be to describe the function as a combination of a predefined library of “simple” functions (such as sin, cos, exp, log, etc.). Although such a description would be easy to reason about, only very special functions can be easily and exactly described in this way. Another way to describe the function would be to list the values of the function as a table. Nevertheless, this representation would only describe the function over part of its domain, since the table would have to be finite. Moreover, reasoning about the function, represented as a table, could be very awkward. Another way to describe a function would be to provide a constraint (i.e. an equation) that is satisfied by only the particular function of interest. Analysis tools can then be used to reason about this equation as a proxy for the function. Difference equations are a type of constraint that can be designed to be satisfied by one, and only one, function such as those described above. The difference equation can be thought of as a “rule” satisfied by a family of functions, and the particular function of interest is then uniquely identified by specifying a particular set of initial values. Difference equations can be very useful for modeling complex phenomenon, especially when characterizing functions that are indexed by a single parameter, such as time. For this reason they are sometimes called “dynamic” or “dynamic systems” and they are a natural model for iterative processes of all types. Often when modeling various phenomenon, such as the efficiency of recursive algorithms, or the change of balance in a checking account, it is easier to identify the

“rule” associated with the difference equation than directly expressing the function of interest. One then needs a procedure to “solve” the difference equation, that is, to describe the desired function given an associated difference equation with the corresponding initial conditions. “Solving” a difference equation involves finding another description of the function that satisfies the constraint characterized by the equation. When the function is described in terms of a library of predefined “simple” functions, we say the solution is in “closedform.” When the function is described as a table of values, we call it a “numerical” solution of the difference equation. “Symbolic solvers” are procedures that compute a closed-form solution to the difference equation. Since they are represented by an explicit formula, such solutions are called “analytic.” Note that most difference equations do not admit analytic solutions. Procedures that compute numerical solutions are called “numerical solvers,” or “numerical methods.” The solutions resulting from numerical procedures are called “algorithmic” because increasingly refined descriptions of the solution (e.g. extending the size of the table of values describing the solution function) are found through more iterations of the procedure. Note that numerical methods apply to a much broader class of difference equations than symbolic solvers. A certain class of difference equations, called linear, constant coefficient with geometric forcing terms, have analytic solutions. This means that we know a procedure that can be systematically executed to obtain the exact solution in terms of a library of predefined “simple” functions. This project will build such a solver for a subset of this class of difference equations that are not too large. We then explore how such methods can be extended to apply to differential equations by approximating the differential equation with an appropriate difference equation. Such procedures form the backbone of scientific computing. To Do: 1. Write the pseudo code for the algorithm for solving linear difference equations with constant coefficients. This should break down into three sequential steps. First, write the algorithm for solving a second order, linear, constant coefficient homogeneous difference equation. Next, extend the algorithm to deal with a geometric forcing term of the form brk, where b and r are real-valued constants, and k is the “time” index variable. Finally, extend the algorithm to use realvalued initial condition information to compute the values of the constants in the general solution. 2. Implement your procedure in C# using Visual Studio. Using Visual C#, write a program implementing your pseudo code of the algorithm. As input, the program will take in the parameters a1, a0, b, r, y(0), y(1) specifying equations of the form: y(k + 2) + a1 y(k + 1) + a0 y(k) = br k with initial conditions y(0), y(1). Output will be the unique solution of this equation with these initial conditions, expressed in

the form: y(k) = c1l 1k + c 2 l k2 + dl k3 or y(k) = c1l 1k + c 2 l k2 + dkl k2 or y(k) = c1l 1k + c 2 kl 1k + dl k2 or y(k) = c1l 1k + c 2 kl 1k + dk 2 l 1k . In other words, you need to compute values for the λ’s, the c’s, and d, and then express the solution in the correct functional form in terms of the symbol k. Make sure the cases of distinct roots, repeated roots, and imaginary roots are appropriately treated and tested. Report: 1. Briefly discuss the difference between analytic and algorithmic solutions to problems. 2. Briefly discuss what you have learned about symbolic vs. numerical solvers. 3. Briefly discuss how the ideas described above for difference equations apply to differential equations. Comment on the importance of difference and differential equations as a modeling tool. List and explain three phenomena that you discovered where differential or difference equations could be used to mathematically model their behavior. 4. Include your pseudo code. 5. Think about how you would analyze the run-time efficiency of your algorithm. Where would most of the computational effort be spent in solving an nth order linear constant coefficient difference equation with geometric forcing terms? 6. One issue when building a numerical solver for a differential equation is discretizing the time axis to transform the differential equation into a difference equation where computational solution techniques can be applied. Consider the first order, linear, homogeneous, constant coefficient differential equation dy(t) + ay(t) = 0 . What is the general solution to this equation (hint: think back dt to calculus!)? What is the unique solution to this equation if you also know that y(0) = 3? Now, suppose that you wanted to sample a solution, y(t), that satisfies this differential equation, at specific time intervals t =kT, k = 0, 1, 2, …, where T is a fixed sampling interval. For any of these intervals, we could think of y(KT) as the initial condition. Using your knowledge of the solution of this differential equation, what would the value of y be at the next sample point i.e. at y(KT+T)? Write a difference equation that exactly samples this differential equation. 7. Derive a difference equation that approximates the differential equation in Part 6 dy(t) using the Forward Euler method, which approximates the derivative by the dt y(KT + T) - y(KT) quantity . Compare the resulting difference equation with the T one found in Part 6 that exactly samples the differential equation, and discuss how the Forward Euler method employs a truncation of the series expansion of function e at to achieve its approximation. List two other techniques that are used by numerical methods to approximate a derivative by a difference. 8. Plot (you may use MATLAB if you wish) the solution to the differential equation in Part 6 for the case a=1 and y(0)=1. Plot the corresponding solution to the exact sampling difference equation found in Part 6. Superimpose these plots so you can check whether they coincide well or not. Now superimpose the solution

to the Forward Euler approximation for different values of the sampling interval, T. Try T=0.5, T=1.5, T=2, and T=3. What happens to the solution as an approximation to the differential equation as |1-T | becomes larger than 1? For what sample intervals is the Forward Euler system a good approximation of the original differential equation? The problem you are observing is called “instability” and it plays a major role in the design of algorithmic solutions for a variety of problems. How does stability of the difference equation relate to correctness of the algorithm for finding a solution to the specified differential equation?