18.330 Introduction to Numerical Analysis Spring 2015 Problem Set 10 Due: May 12, 2015, in class Problem 1: Arbitrary-Precision Arithmetic. Use FFT techniques to implement computer subroutines for arbitrary-precision integer addition and multiplication. (Actually, addition is easy; only multiplication requires FFT techniques.1 ) Use your subroutine to compute the 100th Catalan number C100 , where C0 = 1, C1 = 1, C2 = 2, C3 = 5, Cn+1 = n X Ck Cn−k . k=0 Note that C100 has 57 digits. For a list of 66 different ways to interpret the Catalan numbers, see Richard Stanley’s book Enumerative Combinatorics, Volume 2, Exercise 6.19. (And you thought my problem sets were long.) Problem 2: Numerical differentiation by FFT. In previous problem sets you investigated methods of numerical differentiation—that is, methods of estimating the derivative of a function f (x) given only discrete samples of the function at discrete values of x. In this problem you will explore the use of the FFT as a technique for numerical differentiation. In this problem, let fn ≡ f (nh), n = 0, · · · , N − 1 be samples of f (x) at N evenly-spaced data points. (a) Write a computer subroutine that accepts a vector of function samples {fn } and a number h and uses the (second-order) centered-difference stencil to 1 The easiest way to do this is to represent the decimal digits of an arbitrary-size integer as an array of ordinary integers, with each array element lying between 0 and 9. Adding two numbers then amounts to a simple elementwise addition of the digit arrays, with the slight complication that the result of this procedure may involve “digits” that are greater than 10; remedying this simply involves a “rectification” step (equivalent to “carrying” in elementaryschool arithmetic) in which each digit is reduced to a number between 0 and 9, possibly by augmenting a neighboring digit. Multiplication corresponds to a discrete convolution of the digit arrays (which you will here implement using FFT techniques) followed again by a rectification step. Don’t forget zero padding! 1 compute and return a vector of derivative values {fn0 }. For differentiating at the endpoints, your routine may assume that the underlying function is periodic, so that f0 = fN and fN +1 = f1 . (Note: this subroutine should be approximately one line long.) (b) Now write a subroutine that differentiates using FFT techniques. This subroutine will again accept a vector of function samples and a spacing h. It will compute the DFT of the function samples and use the results to construct a trigonometric interpolant passing through the data points. [In testing your routine, you may find it useful to plot the trigonometric interpolant you construct and confirm that it does run through all the data points.] Then your program will differentiate this interpolant to compute the values fn0 that it will return in the output vector. (c) Compare the accuracy of your subroutines in parts (a) and (b) by testing them on the function f (x) = esin(x) . That is, construct a vector of 25 evenly-spaced samples of this function at points x between 0 and 2π. Also evaluate the derivative f 0 (x) analytically and use this expression to construct a vector of exact derivative values at the same evenly-spaced x points. Then plot, as a function of x, the error in the finite-difference derivative and the error in the FFT derivative as computed using your solutions to parts (b) and (c). (d) For the function f (x) considered above, how would you expect the (i) computational cost (CPU time), (ii) accuracy of the FFT-based differentiation to vary with the number of sample points N (where N = 25 in the above example)? How do these predictions compare with the cost and accuracy scaling of the finite-difference approach of part (a)? For the accuracy in particular, explain why the FFT algorithm exhibits the behavior it does. (e) Now consider using the FFT-based algorithm to differentiate the function g(x) = ex over the same interval. How would you expect the (i) cost, (ii) accuracy of the FFT-based differentiation algorithm to vary with N in this case? [Just state your predictions; you don’t need to implement the calculation in this case, although you may find it instructive to do so.] Problem 3: Chebyshev differentiation and boundary-value problems. In this problem you will revisit the boundary-value problem you investigated on PSet 4. Recall that this problem was d2 x + ω02 x = cos(t2 ), dt2 x(0) = 0.9, x(10) = −2.1, ω0 = 3. (1) In PSet 4 you solved this problem using finite-difference techniques, which exhibit algebraic convergence. Here you will study the same problem using a spectral method and compare the rates of convergence. 2 (a) Let f be a vector of length N + 1 whose entries are samples of a function f at the Chebyshev points, i.e. nπ fn = f (tn ), tn = cos N for n = 0, 1, · · · , N . Similarly, let f 0 be the vector whose entries are samples of the derivative of f : nπ fn0 = f 0 (tn ), tn = cos . N Write a computer program that accepts an integer argument N and computes the (N + 1) × (N + 1) matrix DN that operates on f to yield f 0 . (b) Write a computer program that uses your answer to Part (a) to solve the BVP (1). [To ensure that your code is correct, you may wish to plot the solution x(t) that it predicts and compare against the solution obtained in Problem 2 of PSet 4.] Now compare the values computed for the particular value x(5) by your Chebyshev code and by your finite-difference code. Plot, as a function of the dimension of the linear system solved in each case, the error in this quantity as computed by the two solution methods. 3 Extra credit (15%). In linguistics, a palindrome is a word or phrase that remains invariant under left-right reversal. Examples include • A Santa dog lived as a devil god at NASA. • No, Mel Gibson is a casino’s big lemon. • Yo, banana boy! On the other hand, in mathematics a palindrome is an integer that remains invariant under a left-to-right reversal of its decimal digits. Examples include the integers 2552 and 42677624. Now consider the following algorithm for producing a palindromic number given an arbitrary initial integer n: 1. Set n = n + REVERSE(n). Here REVERSE(n) is the integer whose decimal digits are the left-to-right reversal of the decimal digits of n: for example, REVERSE(259) = 952. 2. Repeat step 1 until n is a palindrome, i.e. we have n=REVERSE(n). For example, given the initial value n=194, the algorithm terminates in 3 iterations: 194 + 491 = 685 685 + 586 = 1271 1271 + 1721 = 2992. A longstanding conjecture holds that this algorithm terminates for any starting value, and evidence supporting the conjecture comes from the fact that, indeed, for most starting values the iteration converges in a small number (∼ 20) of iterations. However, confidence in the conjecture is shaken by the existence of rare starting values for which the iteration seems to take an exceedingly long time to converge. For example, there are 249 integers less than 10,000 for which the iteration fails to converge within 100 iterations. Using your implementation of arbitrary-precision arithmetic, plot the number of iterations required for the above iteration to converge—and the number of digits in the final palindromic number—for starting values of n in the range 10 < n < 300. Do you observe any regularity in the pattern here? Identify at least two integers for which the iteration fails to converge within 100 iterations. References. A lovely discussion on mathematical palindromes may be found in the article “Backward run numbers, letters, words, and sentences until boggles the mind,” by Martin Gardner in the August 1970 issue of Scientific American, http://www.nature.com/scientificamerican/journal/v223/n2/ pdf/scientificamerican0870-110.pdf. Meanwhile, a classic meditation on the crucial psychosomatic benefits of linguistic palindromes (which includes, among other gems, the instant classic Lana C. Ladaug was I ere I saw Guadalcanal ) is “Ainmosni,” by Roger Angell in the May 31, 1969 issue of the New Yorker : http://archives.newyorker. com/?i=1969-05-31#folio=032. 4

© Copyright 2018