918 Chapter 20.Less-Numerical Algorithms Full multiplication of two digit strings,if done by the traditional hand method. is not a fast operation:In multiplying two strings of length N,the multiplicand would be short-multiplied in turn by each byte of the multiplier,requiring O(N2) operations in all.We will see,however,that all the arithmetic operations on numbers of length N can in fact be done in O(N x log N x log log N)operations. The trick is to recognize that multiplication is essentially a comolution(813.1) of the digits of the multiplicand and multiplier,followed by some kind of carry operation.Consider,for example,two ways of writing the calculation 456 x 789: 456 456 ×789 ×789 4104 364554 3648 324048 3192 359784 283542 28671189354 359784 The tableau on the left shows the conventional method of multiplication,in which three separate short multiplications of the full multiplicand(by 9,8,and 7)are added to obtain the final result.The tableau on the right shows a different method (sometimes taught for mental arithmetic),where the single-digit cross products are 9」 all computed (e.g.8 x 6 =48).then added in columns to obtain an incompletely carried result (here,the list 28,67,118,93,54).The final step is a single pass from OF SCIENTIFIC right to left,recording the single least-significant digit and carrying the higher digit or digits into the total to the left (e.g.93+5=98,record the 8,carry 9). You can see immediately that the column sums in the right-hand method are components of the convolution of the digit strings,for example 118=4 x 9+5x 8+6 x 7.In 813.1 we learned how to compute the convolution of two vectors by the fast Fourier transform(FFT):Each vector is FFT'd,the two complex transforms 9 are multiplied,and the result is inverse-FFT'd.Since the transforms are done with floating arithmetic,we need sufficient precision so that the exact integer value of each component of the result is discernible in the presence of roundoff error.We Recipes Numerica 10621 should therefore allow a(conservative)few times log2(log2 N)bits for roundoff 431 in the FFT.A number of length N bytes in radix 256 can generate convolution Recipes components as large as the order of(256)2N,thus requiring 16+log2 N bits of 腿 precision for exact storage.If it is the number of bits in the floating mantissa (cf.$20.1),we obtain the condition North 16+log2 N+few x log2 log2 N<it (20.6.3) We see that single precision,say with it =24,is inadequate for any interesting value of N,while double precision,say with it =53,allows N to be greater than 106,corresponding to some millions of decimal digits.The following routine therefore presumes double precision versions of realft (812.3)and four1(812.2). here called drealft and dfour1.(These routines are included on the Numerical Recipes diskettes.)918 Chapter 20. Less-Numerical Algorithms Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). Full multiplication of two digit strings, if done by the traditional hand method, is not a fast operation: In multiplying two strings of length N, the multiplicand would be short-multiplied in turn by each byte of the multiplier, requiring O(N 2) operations in all. We will see, however, that all the arithmetic operations on numbers of length N can in fact be done in O(N × log N × log log N) operations. The trick is to recognize that multiplication is essentially a convolution (§13.1) of the digits of the multiplicand and multiplier, followed by some kind of carry operation. Consider, for example, two ways of writing the calculation 456 × 789: 456 × 789 4104 3648 3192 359784 456 × 789 36 45 54 32 40 48 28 35 42 28 67 118 93 54 359 784 The tableau on the left shows the conventional method of multiplication, in which three separate short multiplications of the full multiplicand (by 9, 8, and 7) are added to obtain the final result. The tableau on the right shows a different method (sometimes taught for mental arithmetic), where the single-digit cross products are all computed (e.g. 8 × 6 = 48), then added in columns to obtain an incompletely carried result (here, the list 28, 67, 118, 93, 54). The final step is a single pass from right to left, recording the single least-significant digit and carrying the higher digit or digits into the total to the left (e.g. 93 + 5 = 98, record the 8, carry 9). You can see immediately that the column sums in the right-hand method are components of the convolution of the digit strings, for example 118 = 4 × 9+5 × 8+6 × 7. In §13.1 we learned how to compute the convolution of two vectors by the fast Fourier transform (FFT): Each vector is FFT’d, the two complex transforms are multiplied, and the result is inverse-FFT’d. Since the transforms are done with floating arithmetic, we need sufficient precision so that the exact integer value of each component of the result is discernible in the presence of roundoff error. We should therefore allow a (conservative) few times log 2(log2 N) bits for roundoff in the FFT. A number of length N bytes in radix 256 can generate convolution components as large as the order of (256)2N, thus requiring 16 + log2 N bits of precision for exact storage. If it is the number of bits in the floating mantissa (cf. §20.1), we obtain the condition 16 + log2 N + few × log2 log2 N < it (20.6.3) We see that single precision, say with it = 24, is inadequate for any interesting value of N, while double precision, say with it = 53, allows N to be greater than 106, corresponding to some millions of decimal digits. The following routine therefore presumes double precision versions of realft (§12.3) and four1 (§12.2), here called drealft and dfour1. (These routines are included on the Numerical Recipes diskettes.)