Programming Exercises

The following exercises are meant to be answered by either a script or a function Mfile (you decide). The descriptions are complete but the nature of the input/output and display options is left to you. Also, some problems are much simpler than others. Suggestion: start simple and add complexity.


[Back to main page] [Back to exercises page]
1. Compute the value of pi using the following series
How many terms are needed to obtain an accuracy of 1e-12? How accurate is the sum of 100 terms of this series? 2. The Fibonacci numbers are comuted according to the following relation: Fn = Fn-1 + Fn-2 with F0 = F1 = 1. a. Compute the first 10 Fibonacci numbers. b. For the first 50 Fibonacci numbers, compute the ratio Fn / Fn-1 It is claimed that this ratio approaches the value of the golden mean ( (1 + sqrt(5))/2 ). What do your results show? 3. The Legendre polynomials (Pn(x)) are defined by the following recurrance relation (n+1) Pn+1(x) - (2n+1)x Pn(x) + n Pn-1(x) = 0 with P0(x) = 1, P1(x) = x and P2(x) = (3x2 - 1)/2. Compute the next three Legendre polynomials and plot all 6 over the interval [-1,1]. For more information on Legendre polynomials, see the Mathworld site. 4. The present value of an annuity (a yearly sum of money) may be computed from the formula P = (A/i) [(1 + i)n - 1]/(1 + i)n where A is the annuity (in $/year), i is the nominal yearly interest rate (in decimal form), n is the number of years over which the annuity is paid and P is the present value ($). Example computation: If i = 0.15 (15%), A = $100/year and n = 10 years then P = $501.88. If you won the $1,000,000 State Lottery and the Lottery offered you the choice of $500,000 today or $50,000/year for 20 years, which would you take? You can assume an inflation (interest) rate of 5%. 5. I found the following algorithm on a web page* for computing pi: 1. Set a = 1, b = 1/sqrt(2), t = 1/4 and x = 1 2. Repeat the following commands until the difference between a and b is within some desired accuracy: y = a a = (a + b)/2 b = sqrt(b*y) t = t - x*(y - a)^2 x = 2*x 3. From the resulting values of a, b and t, an estimate of pi is Pi_est = ((a + b)^2)/(4*t) How many repeats are needed to estimate pi to an accuracy of 1e-8? 1e-12? Compare the performance of this algorithm to the result from Exercise 1, above. */Gauss_L.html 6. Write a script that asks for an integer (n) and then computes the following based on the value of the integer: While the value of n is greater than 1, replace the integer with half of its value (n/2) if the integer is even. Otherwise, replace the integer with three times its value, plus 1 (3*n + 1). Make provision to count the number of values in (or the length of) the sequence that results. Example calculation: If n = 10, the sequence of integers is 5, 16, 8, 4, 2, 1 and so the length is 6. Make a plot of the length of the sequence that occurs as a function of the integers from 2 to 30. For example, when n = 10, the length is 6 while for n = 15, the length is 17. Is there any pattern? Try larger numbers to see if any pattern occurs. Is there any integer for which the sequence does not terminate? 7. Write a script/function that converts a Roman numeral to its decimal equivalent. There are two distinct situations that you might design your program for: a. The "old" style where the order of the symbols does not matter. In this case, IX and XI both mean 10 + 1 or 11. You should be able to handle the following conversion table: Roman Decimal I 1 V 5 X 10 L 50 C 100 D 500 M 1000 b. The "new" style where the order of the symbols does matter. For example, IX is 9 (10 - 1), XC is 90 (100 - 10). The conversion table given above still holds and you may assume for this case that the only instances of "order" you will encounter are IV (4), IX (9), XL (40), XC (90), CD (400) and CM (900) The function input will be useful here. The format >> str = input('Roman numeral: ','s') will provide a way to get the Roman number into your program as a string. It would be a good idea to try case a. first. 8. Write a function that will do the inverse of the previous problem - convert a decimal number into a Roman number. 9. Compute and plot the path(s) of a set of random walkers which are confined by a pair of barriers at +B units and -B units from the origin (where the walkers all start from). A random walk is computed by repeatedly performing the calculation xj+1 = xj + s where s is a number drawn from the standard normal distribution (randn in MATLAB). For example, a walk of N steps would be handled by the code fragment x(1) = 0; for j = 1:N x(j+1) = x(j) + randn(1,1); end There are three possible ways that the walls can "act": a. Reflecting - In this case, when the new position is outside the walls, the walker is "bounced" back by the amount that it exceeded the barrier. That is, when xj+1 > B, xj+1 = B - |B - xj+1| when xj+1 < (-B), xj+1 = (-B) + |(-B) - xj+1| If you plot the paths, you should not see any positions that are beyond |B| units from the origin. b. Absorbing - In this case, if a walker hits or exceeds the wall positions, it "dies" or is absorbed and the walk ends. For this case, it is of interest to determine the mean lifetime of a walker (i.e., the mean and distribution of the number of steps the "average" walker will take before being absorbed). c. Partially absorbing - This case is a combination of the previous two cases. When a walker encounters a wall, "a coin is flipped" to see if the walker relfects or is absorbed. Assuming a probability p (0 < p < 1) for reflection, the pseudo-code fragment that follows uses the MATLAB uniform random-number generator to make the reflect/absorb decision: if rand < p reflect else absorb end What do you do with all the walks that you generate? Compute statistics, of course. Answering questions like What is the average position of the walkers as a function of time? What is the standard deviation of the position of the walkers as a function of time? Does the absorbing or reflecting character influence these summaries? For the absorbing/partial-reflection case, a plot of the number of surviving walkers as a function of step numbers is a very interesting thing. is useful, informative and interesting, particularly if graphically displayed. 10. The properties of saturated steam are tabulated in many books but that is not a useful form for computer use. The following expressions were provided by Williamson (Chemical Engineering, May 15, 1972, p. 128): Saturation temperature, degrees F Tsat = 8576.65/(15.47538 - ln(Psat)) - 459.216 - 0.023719Psat + (0.84219e-4)Psat2 - (0.70854e-7)Psat3 Liquid specific volume, ft3/# Vliq = 0.01655 + (0.150326e-4)Psat - (0.40488e-7)Psat2 + (0.665584e-10)Psat3 - (0.4053e-13)Psat4 Vapor specific volume, ft3/# Vvap = 430.8419/(Psat + 1.66) + 0.2031 - (0.000258)Psat Liquid specific enthalpy, BTU/# Hliq = 6473.878/(14.01875 - ln(Psat)) - 391.6036 + (0.022915)Psat Vapor specific enthalpy, BTU/# Hvap = 1142.342 + (0.76833)Psat - (0.004194)Psat2 + (0.11642e-4)Psat3 - (0.157e-7)Psat4 + (0.8086e-11)Psat5 In all cases, the pressure is in psia and the equations are valid for pressures between 20 and 600 psia. Maximum error in any computation is below 1%. a. Create a function that provides these steam properties for a given vector of pressures. To save typing, you should be able to copy directly from the browser window and paste into the MATLAB editor. Some further editting will be required, though. b. Write a script that allows a user to get different steam properties as desired. The script should run until the user decides to quit. c. Extend your script to allow the user to specify the units set (e.g., SI, English, CGS) before requesting numerical values for the properties. d. Extend the function from a. to provide the heat of vaporization of water as a property. e. What if the user wants to specify the temperature rather than the pressure? 11. Write a function which computes the cumulative product of the elements in a vector. The cumulative product of the jth element of the vector x, xj, is defined by pj = (x1)(x2) ... (xj) for j = 1:length of the vector x. Create 2 different versions of this function: a. One that uses two for-loops to explicitly carry out the calculations, element by element. An "inner" loop should accumulate the product and an "outer" loop should more through the elements of the vector p. b. One that uses the built-in function prod to replace the inner loop. In each case, you can check your results with the built-in function, cumprod. 12. Follow the directions for Exercise 11 but create a function that computes the cumulative sum of the elements of a vector. The elements of the cumulative sum vector are defined by sj = x1 + x2 + ... + xj for j = 1: length of the vector x. The built-in functions sum and cumsum should be used instead of prod and cumprod, respectively. 13. Create a function that generates random arrays of integers beween a and b, inclusive. Use the definition line function A = randint(a,b,M,N) where a and b define the (inclusive) range and M and N define the size of the output array (rows and columns, respectively). a. Test your function with the following piece of code: x = randint(10,17,100000,1); hist(x,10:17) The resulting histogram should be nearly flat from 10 to 17. Pay particular attention to the endpoints of the distribution. b. Test your function with the following pieces of code: x = randint(-30,5,100000,1); hist(x,-30:5) x = randint(-45,-35,100000,1); hist(x,-45:-35) x = randint(7,-2,100000,1); hist(x,-2:7) Consider that each of these uses is valid. c. Modify your code to allow for missing/default inputs. Use the behavior of rand to help design your code. You will need to use nargin in your routine. For example, the function should be able to return - a 5-by-5 array of integers between 1 and 20: e.g., A = randint(1,20,5) - a single random integer between 10 and 50: e.g., A = randint(10,50) - an empty array if not enough inputs are provided.

[Back to main page] [Back to exercises page]
Comments? Contact Jim Maneval at maneval@bucknell.edu