User Manual

ManualsBrandsHP ManualsCalculators, organizers15c

HEWLETT-PACKARD

HP-15C

ADVANCED FUNCTIONS

HANDBOOK

Summary of content (186 pages)

PAGE 1
HEWLETT-PACKARD HP-15C ADVANCED FUNCTIONS HANDBOOK
PAGE 2
Legal Notice This manual and any examples contained herein are provided “as is” and are subject to change without notice. Hewlett-Packard Company makes no warranty of any kind with regard to this manual, including, but not limited to, the implied warranties of merchantability non-infringement and fitness for a particular purpose. In this regard, HP shall not be liable for technical or editorial errors or omissions contained in the manual.
PAGE 3
HEWLETT PACKARD HP-15C Advanced Functions Handbook August 1982 00015-90011 Printed in U.S.A.
PAGE 4
Contents Contents ...............................................................................................................4 Introduction .........................................................................................................7 Section 1: Using _ Effectively ..................................................................9 Finding Roots .........................................................................................................................9 How _ Samples ...........
PAGE 5
Contents 5 Applications ......................................................................................................................... 65 Storing and Recalling Complex Numbers Using a Matrix .............................................. 65 Calculating the nth Roots of a Complex Number ............................................................ 67 Solving an Equation for Its Complex Roots .................................................................... 69 Contour Integrals .................
PAGE 6
PAGE 7
Introduction The HP-15C provides several advanced capabilities never before combined so conveniently in a handheld calculator:  Finding the roots of equations.  Evaluating definite integrals.  Calculating with complex numbers.  Calculating with matrices. The HP-15C Owner's Handbook gives the basic information about performing these advanced operations. It also includes numerous examples that show how to use these features.
PAGE 8
PAGE 9
Section 1: Using _ Effectively The _ algorithm provides an effective method for finding a root of an equation. This section describes the numerical method used by _ and gives practical information about using _ in various situations. Finding Roots In general, no numerical technique can be guaranteed to find a root of every equation that has one.
PAGE 10
Section 1: Using _ Effectively If c isn't a root, but f(c) is closer to zero than f(b), then b is relabeled as a, c is relabeled as b, and the prediction process is repeated. Provided the graph of f(x) is smooth and provided the initial values of a and b are close to a simple root, the secant method rapidly converges to a root.
PAGE 11
Section 1: Using _ Effectively 11 If _ hasn't found a sign change and a sample value c doesn't yield a function value with diminished magnitude, then _ fits a parabola through the function values at a, b, and c. _ finds the value d at which the parabola has its maximum or minimum, relabels d as a, and then continues the search using the secant method. _ abandons the search for a root only when three successive parabolic fits yield no decrease in the function magnitude or when d = b.
PAGE 12
Section 1: Using _ Effectively In general, every equation is one of an infinite family of equivalent equations with the same real roots. And some of those equations must be easier to solve than others. While _ may fail to find a root for one of those equations, it may succeed with another. Inaccurate Equations _ can't calculate an equation's root incorrectly unless the function is incorrectly calculated. The accuracy of your function subroutine affects the accuracy of the root that you find.
PAGE 13
Section 1: Using _ Effectively 13 This function equals zero at no more than n real values of x, called zeros of the polynomial. A limit to the number of positive zeros of this function can be determined by counting the number of times the signs of the coefficients change as you scan the polynomial from left to right. Similarly, a limit to the number of negative zeros can be determined by scanning a new function obtained by substituting −x in place of x in the original polynomial.
PAGE 14
Section 1: Using _ Effectively First convert 5° 18'S to decimal degrees (press 5.18”|À), obtaining −5.3000 (using •4 display mode). (Southern latitudes are expressed as negative numbers for calculation purposes.) The solution to Coulerre's problems is the value of t satisfying −5.3000 = a4t4 + a3t3 + a2t2 + a1t + a0. Expressed in the form required by _ the equation is 0 = a4t4 + a3t3 + a2t2 + a1t − 2.8806 where the last, constant term now incorporates the value of the declination.
PAGE 15
Section 1: Using _ Effectively 15 In Run mode, key in the five coefficients: Keystrokes |¥ 4.2725 “8” O4 1.9931”“ 5” O3 1.0229“3” O2 3.7680“1” O1 2.8806” O0 Display Run mode. 4.2725 4.2725 -08 -08 -1.9931 1.0229 0.0010 3.7680 0.3768 -2.8806 -05 -03 Coefficient of t4. Coefficient of t3. Coefficient of t2. -01 Coefficient of t. Constant term. Because the desired solution should be between 1 and 32, key in these values for initial estimates. Then use _ to find the roots.
PAGE 16
Section 1: Using _ Effectively Keystrokes − ) ) |(|( ´_A − ) ) Display 278.4497 276.7942 7.8948 278.4497 Error 8 278.4398 278.4497 7.8948 Last estimate tried. A previous estimate. Nonzero value of function. Restores stack to original state. Again, no root found. Approximately same estimate. A previous estimate. Same function value. You have found a positive local minimum rather than a root. Now try to find the negative root Keystrokes 1000 ”v 1100 ” ´_A ) ) Display -1,000.0000 -1,100 -108.
PAGE 17
Section 1: Using _ Effectively 17 For some systems of equations, expressed as f1(x1, …, xn) = 0 ⋮ fn(x1, …, xn) = 0 it is possible through algebraic manipulation to eliminate all but one variable. That is, you can use the equations to derive expressions for all but one variable in terms of the remaining variable. By using these expressions, you can reduce the problem to using _ to find the root of a single equation.
PAGE 18
Section 1: Using _ Effectively A final consideration for this example is to choose the initial estimates that would be appropriate. Because the argument of the inverse cosine must be between −1 and 1, x must be more negative than about −0.1059 (found by trial and error or by using _). The initial guesses might be near but more negative than this value, −0.11 and −0.2 for example. (The complex equation used in this example is solved using an iterative procedure in the example on page 69.
PAGE 19
Section 1: Using _ Effectively Keystrokes |¥ ´CLEAR M ´ b0 \ l0 * \ l1 ® [ ÷ ® ] ÷ ” ® \ l0 * [ l0 * + l2 * |n Display Program mode. 000001-42,21, 0 00224 00345 0 00420 00524 00645 1 00730 00834 00923 01010 01134 01225 01310 01416 01534 01624 01745 0 01820 01923 02045 0 02120 02240 02345 2 02420 02543 32 In Radians mode, calculate and store the three constants. Keystrokes |¥ |R 2|$* .6*O0 \O1 ”1+ ⁄O2 Display Run mode. Specifies Radians mode. 6.2832 3.7699 -0.8090 1.8090 0.5528 Constant r0.
PAGE 20
Section 1: Using _ Effectively The relative field intensity is maximum at an angle of 90° (perpendicular to the tower). To find the minimum, use angles closer to zero as initial estimates, such as the radian equivalents of 10° and 60°. Keystrokes 10´r 60´r ´_0 )) |(|( |d Display 0.1745 1.0472 0.4899 -5.5279 0.4899 28.0680 -10 Initial estimates. Angle giving zero slope. Slope at specified angle. Restores the stack. Angle in degrees. The relative field intensity is most negative at an angle of 28.
PAGE 21
Section 1: Using _ Effectively 21 The accuracy of this approximation depends upon the increment Δ and the nature of the function. Smaller values of Δ give better approximations to the derivative, but excessively small values can cause round-off inaccuracy. A value of x at which the slope is zero is potentially a local extreme of the function. Example: Solve the previous example without using the equation for the derivative dE/dθ.
PAGE 22
Section 1: Using _ Effectively 22 Keystrokes Display ® [ ÷ l2 * |n 029030031032033034- 34 23 10 45 2 20 43 32 In the previous example, the calculator was set to Radians mode and the three constants were stored in registers R0, R1, and R2. Key in the same initial estimates as before and execute _. Keystrokes Display |¥ 10´r 60´r ´_A )) |(|( vv´B ® |d Run mode. 0.1745 1.0472 0.4899 0.0000 0.4899 -0.2043 Initial estimates. Angle given zero slope. Slope at specified angle. Restores stack.
PAGE 23
Section 1: Using _ Effectively  23 If a root of g(x) is found, either the number e is not beyond the extreme value of f(x) or else _ has found a different region where f(x) equals e. Revise e so that it is close to—but beyond—the extreme value of f(x) and try _ again. It may also be possible to modify g(x) in order to eliminate the distant root. Example: Solve the previous example without calculating the derivative of the relative field intensity E.
PAGE 24
Section 1: Using _ Effectively Based on these samples, try using an extreme estimate of −0.25 and initial _ estimates (in radians ) near 10° and 30°. Keystrokes Display .25”O9 .2v .6 ´_1 −O4 )O5 ) .9* O+9 l4 vv´B − l5 -0.2500 0.2000 0.6 Error 8 0.4849 0.4698 0.0457 0.0411 0.0411 0.4849 -0.2043 0.0000 0.4698 ´_1 − ® ® vv´B ® |d |D Error 8 0.4898 0.4893 0.4898 -0.2043 0.4898 28.0660 28.0660 Stores extreme estimate. Initial estimates. No root found. Stores θ estimate. Stores previous θ estimate.
PAGE 25
Section 1: Using _ Effectively 25 n The number of compounding periods. (For example, a 30 year loan with monthly payments has n = 12 x 30 = 360.) i The interest rate per compounding period expressed as a percent. (To calculate i, divide the annual percentage rate by the number of compounding periods in a year. That is, 12% annual interest compounded monthly equals 1% periodic interest.) PV The present value of a series of future cash flows or the initial cash flow. PMT The periodic payment amount.
PAGE 26
Section 1: Using _ Effectively 26 compounding period—months or years, for example. Vertical arrows represent exchanges of money, following the convention that an upward arrow (positive) represents money received and a downward arrow (negative) represents money paid out. (The examples that follow are illustrated using cash flow diagrams.) Pressing ´CLEARQ provides a convenient way to set up the calculator for a new problem. However, it isn't necessary to press ´CLEARQ between problems.
PAGE 27
Section 1: Using _ Effectively 27  If a problem has a specified interest rate of 0, the program generates an Error 0 display (or Error 4 when solving for PMT).  Problems with extremely large (greater than 106) or extremely small (less than 10−6) values for n and i may give invalid results:  Interest problems with balloon payments of opposite signs to the periodic payments may have more than one mathematically correct answer (or no answer at all).
PAGE 28
Section 1: Using _ Effectively Keystrokes Display v “ ” 3 |"1 02936 03026 03116 0323 033-43, 5, 1 ´_3 t4 t0 ´b4 “ 2 * O2 |n ´bC O3 ¦ G1 G2 ” O3 |n ´bÁ O4 ¦ 1 O4 G1 l3 G2 ® ÷ ” O4 |n ´bE O5 034-42,10, 3 03522 4 03622 0 037-42,21, 4 03826 0392 04020 04144 2 04243 32 043-42,21,13 04444 3 04531 04632 1 04732 2 04816 04944 3 05043 32 051-42,21,14 05244 4 05331 0541 05544 4 05632 1 05745 3 05832 2 05934 06010 06116 06244 4 06343 32 064-42,21,15 06544 5 28 Clears flag 1 for _ subroutine. Calculates i.
PAGE 29
Section 1: Using _ Effectively Keystrokes Display ¦ G1 l+3 l÷7 ” O5 |n ´b1 |F1 1 l2 |k ´b3 O8 1 O0 + |T4 t0 O6 |?0 06631 06732 1 068-45,40, 3 069-45,10, 7 07016 07144 5 07243 32 073-42,21, 1 074-43, 4, 1 0751 07645 2 07743 14 078-42,21, 3 07944 8 0801 08144 0 08240 083-43,30, 4 08422 0 08544 6 086-43, 6, 0 O0 l1 ” Y O7 1 ® |~ t0 l*0 l4 l÷8 * |?1 |n 08744 0 08845 1 08916 09014 09144 7 0921 09334 09430 09543 20 09622 0 097-45,20, 0 09845 4 099-45,10, 8 10020 101-43, 6, 1 10243 32 Calculates FV.
PAGE 30
Section 1: Using _ Effectively 30 Keystrokes l+3 ´b2 l5 l*7 + |n Display 103-45,40, 3 104-42,21, 2 10545 5 106-45,20, 7 10740 10843 32 _ subroutine continues. Calculates FV(1 + i/100)−n. _ subroutine ends. Labels used: A, B, C, D, E, 0, 1, 2, 3, and 4. Registers used: R0 (A), R1 (n), R2 (i), R3 (PV), R4 (PMT), R5 (FV), R6, R7, and R8. To use the program: 1. Press 8´m% to reserve R0 through R8. 2. Press ´U to activate User mode. 3. If necessary, press ´CLEARQ to clear all of the financial variables.
PAGE 31
Section 1: Using _ Effectively 31 Example: You place $155 in a savings account paying 5¾% compounded monthly. What sum of money can you withdraw at the end of 9 years? Keystrokes |¥ ´CLEARQ ´•2 ´U |"0 9v12*A 5.75v12÷B 155”C E¦ Display Run mode. Clears financial variables. 108.00 0.48 -155.00 259.74 Activates User mode. Ordinary annuity. Enters n = 9 × 12. Enters i = 5.75 / 12. Enters PV = −155 (money paid out). Calculates FV.
PAGE 32
Section 1: Using _ Effectively Keystrokes ´CLEARQ 30v12*A 13v12÷B 30000C Á¦ Display 360.00 1.08 30,000.00 -331.86 Clears financial variables Enters n = 30 × 12. Enters i = 13/12. Enters PV = 30,000. Calculates PMT (money paid out). Example: You offer a loan of $3,600 that is to be repaid in 36 monthly payments of $100 with an annual interest rate of 10%.
PAGE 33
Section 1: Using _ Effectively Keystrokes ´CLEARQ 360A 14v12÷B 50000”C Á¦ 33 Display 360.00 1.17 -50,000.00 592.44 Clears financial variables Enters n = 360. Enters i = 14/12. Enters PV = −50,000. Calculates PMT. Now calculate the remaining balance at month 24. Keystrokes 24A E¦ Display 24.00 49,749.56 Enters n = 24. Calculates FV at month 24. Store this remaining balance, then calculate the remaining balance at month 12 and the principal reduction between payments 12 and 24.
PAGE 34
Section 1: Using _ Effectively Keystrokes Display ´CLEARQ |F0 5v12*A 13v12÷B 63000”C 10000E Á¦ 60.00 1.08 -63,000.00 10,000.00 1,300.16 Clears financial variables. Specifies beginning of period payments. Enters n = 5 × 12. Enters i = 13/12. Enters PV = −63,000. Enters FV = 10,000. Calculates PMT. If the prices of the computer increases to $70,000, what should the payments be? Keystrokes 70000”C Á¦ Display -70,000.00 1,457.73 Enters PV = −70,000. Calculates PMT.
PAGE 35
Section 1: Using _ Effectively 35 flows must occur at equal intervals; if no cash flow occurs for several time periods, enter 0 for the cash flow amount and the number of zero cash flows in that group. After all the cash flows have been stored in matrix C, you can enter an assumed interest rate and calculate the net present value (NPV) of the investment. Alternatively, you can calculate the internal rate of return (IRR).
PAGE 36
Section 1: Using _ Effectively Keystrokes “ ” 3 ´_2 t1 t0 ´b1 “ 2 * ¦ ´b2 |"0 O2 1 O4 + |T4 t0 O3 0 O5 ´>1 ´b3 |?0 t7 G6 l2 |~ t4 1 + G6 ” Y O4 1 ® Display 01026 01116 0123 013-42,10, 2 01422 1 01522 0 016-42,21, 1 01726 0182 01920 02031 021-42,21, 2 022-43, 5, 0 02344 2 0241 02544 4 02640 027-43,30, 4 02822 0 02944 3 0300 03144 5 032-42,16, 1 033-42,21, 3 034-43, 6, 0 03522 7 03632 6 03745 2 03843 20 03922 4 0401 04140 04232 6 04316 04414 04544 4 0461 04734 36 Branch for no IRR solution.
PAGE 37
Section 1: Using _ Effectively Keystrokes l÷2 l*3 t5 ´b4 ® G6 ´b5 * O+5 l4 O*3 t3 ´b6 ´UlC ´U |n | F0 |n ´b7 l5 |n 37 Display 04830 049-45,10, 2 050-45,20, 3 05122 5 052-42,21, 4 05334 05432 6 055-42,21, 5 05620 057-44,40, 5 05845 4 059-44,20, 3 06022 3 061-42,21, 6 062u 45 13 06343 32 064-43, 4, 0 06543 32 066-42,21, 7 06745 5 06843 32 Recalls cash flow element. Sets flag 0 if last element. Recalls NPV. Labels used: A, B, and 0 through 7. Registers used: R0 through R5. Matrix used: C.
PAGE 38
Section 1: Using _ Effectively  To calculate NPV, enter periodic interest rate i in percent and press A. Repeat for as many interest rates as needed. 7. Repeat steps 3 through 6 for other sets of cash flows. Example: An investor pays $80,000 for a duplex that he intends to sell after 7 years. He must spend some money the first year for repairs. At the end of the seventh year the duplex is sold for $91,000.
PAGE 39
Section 1: Using _ Effectively 39 Since the NPV is negative, the investment does not achieve the desired 9% yield. Calculate the IRR. Keystrokes B Display 8.04 IRR (after about 5 seconds). The IRR is less than the desired 9% yield. Example: An investment of $620,000,000 is expected to have an annual income stream for the next 15 years as shown in the diagram. What is the expected rate of return? Keystrokes 3v2 ´mC ´>1 620000000” OC 1OC 100000000OC 10OC 5000000OC 5OC B ´•4 ´U Display 2 2.00 2.
PAGE 40
Section 2: Working with f The HP-15C gives you the ability to perform numerical integration using f. This section shows you how to use f effectively and describes techniques that enable you to handle difficult integrals. Numerical Integration Using f A calculator using numerical integration can almost never calculate an integral precisely. But the f function asks you in a convenient way to specify how much error is tolerable.
PAGE 41
Section 2: Working with f 41 The HP-15C doesn't prevent you from declaring that f(x) is far more accurate than it really is. You can specify the display setting after a careful error analysis, or you can just offer a guess. You may leave the display set to i or • 4 without much further thought. You will get an estimate of the integral and its uncertainty, enabling you to interpret the result more intelligently than if you got the answer with no idea of its accuracy or inaccuracy.
PAGE 42
Section 2: Working with f Functions Related to Physical Situations Functions like cos(4 - sin) are pure mathematical functions. In this context, this means that the functions do not contain any empirical constants, and neither the variables nor the limits of integration represent actual physical quantities. For such functions, you can specify as many digits as you want in the display format (up to nine) to achieve the desired degree of accuracy in the integral.
PAGE 43
Section 2: Working with f 43 mathematical operations—may not be accurate to all 10 digits that can be displayed. Note that round-off error affects the evaluation of any mathematical expression, not just the evaluation of a function to be integrated using f. (Refer to the appendix for additional information.) If f(x) is a function relating to a physical situation, its inaccuracy due to round-off typically is insignificant compared to the inaccuracy due to empirical constants, etc.
PAGE 44
Section 2: Working with f 4. To get the integral over the entire interval of integration, add together the approximations and their uncertainties from the integrals calculated over each subinterval. You can do this easily using the z key. Before subdividing the integration, check whether the calculator underflows when evaluating the function around the upper (or lower) limit of integration.
PAGE 45
Section 2: Working with f Keystrokes 45 Display 200v vv G1 2.000 2.000 2.768 00 02 -85 225v vv G1 2.250 2.250 4.324 02 02 -96 Try a smaller value of x. Calculator doesn’t underflow at x = 200; try a number between 200 and 250. Calculator is close to underflow. At this point, you can use _ to pinpoint the smallest value of x at which the calculator underflows. Keystrokes Display ) 2.250 02 ´_1 2.280 02 Roll down stack until the last value tried is in the X- and Yregisters.
PAGE 46
Section 2: Working with f Keystrokes Display ® )) 1.841 1.000 228 228 ´i0 2. 02 ´f1 5. -04 ´i3 5.328 -04 ® 7.568 -05 ® 5.328 -04 z 2.000 00 lz 1.000 00 ® 2.598 -04 -04 01 46 Uncertainty of approximation. Roll down stack until upper limit of first integral appears in X-register. Keys upper limit of second integral into X-register. Upper limit of first integral is lifted into Y-register, becoming lower limit of second integral.
PAGE 47
Section 2: Working with f 47 Transformation of Variables In many problems where the function changes very slowly over most of a very wide interval of integration, a suitable transformation of variables may decrease the time required to calculate the integral. For example, consider again the integral  . Let e x  u 3 Then x  3 ln u And dx  3  0 xe  x dx du .
PAGE 48
Section 2: Working with f 48 The approximation agrees with the value calculated In the previous problem for the same integral. Evaluating Difficult Integrals Certain conditions can prolong the time required to evaluate an integral or can cause inaccurate results. As discussed in the HP-15C Owner's Handbook, these conditions are related to the nature of the integrand over the interval of integration. One class of integrals that are difficult to calculate is improper integrals.
PAGE 49
Section 2: Working with f 49 Although the branch for u=1 adds four steps to your subroutine, integration near x = 0 becomes more accurate. As a second example, consider the integral 1  x 1    x  1  ln x dx . 0   The derivative of the integrand approaches ∞ as x approaches 0, as shown in the illustration below. By substituting x = u2, the function becomes more well behaved, as shown in the second illustration.
PAGE 50
Section 2: Working with f  t 0 f ( x)dx  tan 1 t 0 e  tan u (1  tan 2 u )du , 2 which is calculated readily even with t as large as 1010. Using the same substitution with g(x), values near a = 0 and b = 10−5 provide good results. This example involves subdividing the interval of integration. Although a function may have features that look extreme over the entire interval of integration, over portions of that interval the function may look more well-behaved.
PAGE 51
Section 2: Working with f 51 r  c / a  0.01 However, this integral is nearly improper because q and r are both so nearly zero. But by using an integral in closed form that sufficiently resembles the troublesome part of V, the difficulty can be avoided. Try 1 1 r r W  p  du / u 2  q  p ln(u  u 2  q )  p ln((1  1  q ) /( r  r 2  q ))  8.40181880708  10 6. Then 1 V  W  p  ( (1  u 2 ) /(u 2  q)  1 / u 2  q ) du r  1 W / p u2 du.
PAGE 52
Section 2: Working with f 52 The program has the following characteristics:  The display format specifies the accuracy of the integrand in the same way as it does for f itself. However, if you specify an unnecessarily large number of display digits, the calculation will be prolonged.  Small function values, such as Q(20), P(−20), and erfc(10), are accurately computed as quickly as moderate values.
PAGE 53
Section 2: Working with f Keystrokes 53 Display ´b4 |"1 O1 ® O0 |a 1 . 6 |T8 t6 0 l0 ´f0 2 * ´b3 |$ ¤ ÷ |n ´b6 030-42,21, 4 031-43, 5, 1 03244 1 03334 03444 0 03543 16 0361 03748 0386 039-43,30, 8 04022 6 0410 04245 0 043-42,20, 0 0442 04520 046-42,21, 3 04743,26 04811 04910 05043 32 051-42,21, 6 |F1 0 l0 |x ” ' ´f1 052-43, 4, 1 0530 05445 0 05543 11 05616 05712 058-42,20, 1 G3 l0 v |a ÷ * l1 |K 059060061062063064065066- 32 45 3 0 36 43 16 10 20 45 1 43 36 Subroutine for erf(x) or erfc(x).
PAGE 54
Section 2: Working with f 54 Keystrokes Display + 067068- ” |n ´b0 |x ” ' |n ´b1 06916 07043 32 071-42,21, 0 07243 11 07316 07412 07543 32 076-42,21, 1 |~ |n |N ” ¤ ⁄ |n 077078079080081082083- 30 40 Adjusts integral for sign of x and function. Subroutine to calculate Subroutine to calculate (-ln u)−1/2. 43 20 43 32 43 12 16 11 15 43 32 Labels used: A, B, C, E, 0, 1, 2, 3, 4, 5, and 6. Registers used: R0 (x), R1, R2. Flag used: 1. To use this program: 1. Enter the argument x into the display.
PAGE 55
Section 2: Working with f Keystrokes Display 1.234´A .5´E 8.914 5.205 -01 -01 55 P(1.234). erf(0.5). Example: For a Normally distributed random variable X with mean 2.151 and standard deviation 1.085, calculate the probability Pr [2 < X  3].  2  2.151 X   3  2.151 Pr[2  X  3]  Pr     1.085   1.085  3  2.151   2  2.151   P   P   1.085   1.085  Keystrokes 2v 2.1511.085÷ ´A O3 3v 2.1511.085÷ ´A l3 ´•4 Display 2.000 -1.510 -1.392 4.447 4.447 3.000 8.490 7.825 7.830 4.
PAGE 56
Section 3: Calculating in Complex Mode Physically important problems involving real data are often solved by performing relatively simple calculations using complex numbers. This section gives important insights into complex computation and shows several examples of solving problems involving complex numbers. Using Complex Mode Complex mode in the HP-15C enables you to evaluate complex-valued expressions simply.
PAGE 57
Section 3: Calculating in Complex Mode Keystrokes 6 OV ® v v v l6 ´b1 + l% ® ÷ ´eV t1 l0 + ® |K |N |K . 5 * + |n Display 0026 00344 25 00434 00536 00636 00736 00845, 6 009-42,21, 1 01040 01145 24 01234 01310 014-42, 5,25 01522 1 01645 0 01740 01834 01930 02043 36 02143 12 02243 36 02348 0245 02530 02620 02740 02843 32 Stores counter in Index register. Fills stack with z. Recalls a6. Loop for continued fraction Recalls ai. Restores z. Decrements counter. Recalls a0. Restores z. Recalls z.
PAGE 58
Section 3: Calculating in Complex Mode Keystrokes 195v371÷ O4 1.011523068O5 1.517473649 O6 Display 0.5256 0.5256 1.0115 1.5175 Stores a4. Stores a5. Stores a6. Use this program to calculate ln((4.2)), then compare it with ln(3.2!) calculated with the ! function. Also calculate ln((1 + 5i)). Keystrokes Display 4.2´A ´•9 3.2´! |N 1v 5´V 2.0486 2.048555637 7.756689536 2.048555637 1.000000000 1.000000000 ´A ´} -6.130324145 3.815898575 ´•4 Calculates ln((4.2)). Displays 10 digits.
PAGE 59
Section 3: Calculating in Complex Mode 59 Definitions of Math Functions The lists that follow define the operation of the HP-15C in Complex mode. In these definitions, a complex number is denoted by z = x + iy (rectangular form) or z = rei (polar form). Also z  x2  y 2 .
PAGE 60
Section 3: Calculating in Complex Mode The illustrations that follow show the principal branches of the inverse relations. The lefthand graph in each figure represents the cut domain of the inverse function; the right-hand graph shows the range of the principal branch. For some inverse relations, the definitions of the principal branches are not universally agreed upon. The principal branches used by the HP-15C were carefully chosen.
PAGE 61
Section 3: Calculating in Complex Mode 61
PAGE 62
Section 3: Calculating in Complex Mode The principal branches in the last four graphs above are obtained from the equations shown, but don't necessarily use the principal branches of ln(z) and z . The remaining inverse functions may be determined from the illustrations above and the following equations: LOG(z) = LN(z) / LN(10) SINH−1(z) = −i SIN−1(iz) TANH−1(z) = −i TAN−1(iz) wz = ez LN(w).
PAGE 63
Section 3: Calculating in Complex Mode 63 Using _ and f in Complex Mode The _ and f functions use algorithms that sample your function at values along the real axis. In Complex mode, the _ and f functions operate with only the real stack, even though your function subroutine may have complex computations in it. For example, _ will not search for the roots of a complex function, but rather will sample the function on the real axis and search for a zero of the function's real part.
PAGE 64
Section 3: Calculating in Complex Mode To require approximations with accurate components is to demand more than keeping relative errors small. For example, consider this problem in which the calculations are done with four significant digits. It illustrates the limitations imposed on a complex calculation by finite accuracy. z1 = Z1 = 37.1 + 37.3i z2 = Z2 = 37.5 + 37.3i and Z1 × Z2 = (37.10 × 37.50 − 37.30 × 37.30) + i(37.10 × 37.30 + 37.30 × 37.50) = (1391. − 1391.) + i(1384. + 1399.) = 0 + i(2783.
PAGE 65
Section 3: Calculating in Complex Mode 65 then think of Z as 0.0000001234567890 × 10−3 + i(2.222222222 × 10−3). then the accurate digits are 0.000000123 × 10−3 + i(2.222222222 × 10−3). Applications The capability of the HP-15C to work with complex numbers enables you to solve problems that extend beyond the realm of real-valued numbers. On the following pages are several programs that illustrate the usefulness of complex calculations—and the HP-15C.
PAGE 66
Section 3: Calculating in Complex Mode Keystrokes Display |n ´bE O0 |` 2 O1 ) 0 + 01143 32 012-42,21,15 01344 0 01443 35 0152 01644 1 01733 0180 01940 lC ´} ´s1 |` 02045 13 02142 30 022-42, 5, 1 02343 35 lC |n 024025- 45 13 43 32 “Recall” program. R0 = k. Disables stack. Sets R1 = 2. Sets stack for another try if Error 3 occurs next. Recalls b (imaginary part). Decrements R1 = 1. Disables stack and clears real X-register. Recalls a (real part).
PAGE 67
Section 3: Calculating in Complex Mode 67 Calculating the nth Roots of a Complex Number This program calculates the nth roots of a complex number. The roots are zk for k = 0, 1, 2, ... , n - 1. You can also use the program to calculate z1/r, where r isn't necessarily an integer. The program operates the same way except that there may be infinitely many roots zk for k = 0, ±1, ±2, ... .
PAGE 68
Section 3: Calculating in Complex Mode Keystrokes Display ´V * 029030- 42 25 20 lV ® 031032- 45 25 34 1 O+V 0331 034-44,40,25 ) 035- 33 ¦ t0 036037- 31 22 0 Forms complex z0. Calculates z0eik360 / n, root number k. Recalls number k. Places zk in X-registers, k in Yregister. Increments number k in Index register. Restores zk and k to X- and Yregisters. Halts execution. Branch for next root. Labels used: A and 0. Registers used: R2, R3, R4, and Index register. To use this program: 1.
PAGE 69
Section 3: Calculating in Complex Mode Keystrokes ¦ ´% (hold) 50O V ¦ ´% (hold) 69 Display 0.9980 0.0628 50.0000 Calculates z1 (real part). Imaginary part of z1. Stores root number in Index register. Calculates z50 (real part). Imaginary part of z50 -1.0000 0.0000 Solving an Equation for Its Complex Roots A common method for solving the complex equation f(z) = 0 numerically is Newton's iteration.
PAGE 70
Section 3: Calculating in Complex Mode negative values of Re(z). This would slow convergence considerably unless the first guess z0 were extremely close to a root. In addition, this f(z) vanishes infinitely often, so it's difficult to determine when all desired roots have been calculated. But by rewriting this equation as ez = −8/(z + 9) and taking logarithms, you obtain an equivalent equation z = ln(−8/(z + 9)) ± i2nπ for n = 0, 1, 2, ….
PAGE 71
Section 3: Calculating in Complex Mode Keystrokes Display O0 ® ´V |n ´bB v v l1 ´} ® 9 + ÷ |K ® |N * ® l1 ´} l+0 |K ) + ® 1 0 + ÷ + |n ´bC 01644 0 01734 01842 25 01943 32 020-42,21,12 02136 02236 02345 1 02442 30 02534 0269 02740 02810 02943 36 03034 03143 12 03220 03334 03445 1 03542 30 036-45,40, 0 03730 03843 36 03933 040 40 04134 0421 0430 04440 04510 04640 04743 32 048-42,21,13 v ' 9 |K + 049050051052053- 36 12 9 43 36 40 Forms complex A(n). Program for zk + 1. Creates i Im(A(n)).
PAGE 72
Section 3: Calculating in Complex Mode Keystrokes Display 054055056057058059- 8 ® ÷ + |a |n 8 34 10 40 43 16 43 32 Calculates |ez + 8/(z + 9)|. Labels used: A, B, and C. Registers used: R0 and R1. Now run the program. For each root, press B until the displayed real part doesn't change. (You might also check that the imaginary part doesn't change.) Keystrokes |¥ ´U 0A B B B ´% (hold) C ® Display 1.6279 -0.1487 -0.1497 -0.1497 2.8319 1.0000 -0.1497 Run mode. Activates User mode.
PAGE 73
Section 3: Calculating in Complex Mode 73 Since all roots have negative real parts, the system is stable, but the margin of stability (the smallest in magnitude among the real parts, namely -0.1497) is small enough to cause concern if the system must withstand much noise. Contour Integrals You can use f to evaluate the contour integral  f ( z)dz , where C is a curve in the complex C plane. First parameterize the curve C by z(t)= x(t) + i y(t) for t1 ≤ t ≤ t2. Let G(t)=f(z(t))z’(t).
PAGE 74
Section 3: Calculating in Complex Mode Keystrokes v 1 ´f0 O2 ) O3 ) ´f1 l2 ´V ® l3 ´V ® |n ´b0 G1 ´} |n ´b1 l4 l5 ´V * l6 l7 ´V + GB l4 l5 ´V * |n Display 01236 0131 014-42,20, 0 01544 2 01633 01744 3 01833 019-42,20, 1 02045 2 02142 25 02234 02345 3 02442 25 02534 02643 32 027-42,21, 0 02832 1 02942 30 03043 32 031-42,21, 1 03245 4 03345 5 03442 25 03520 03645 6 03745 7 03842 25 03940 04032 12 04145 4 04245 5 04342 25 04420 04543 32 Labels used: A, 0, and 1.
PAGE 75
Section 3: Calculating in Complex Mode 75 To use this program: 1. Enter your function subroutine labeled "B" into program memory. 2. Press 7´m% to reserve registers R0 through R7. (Your subroutine may require additional registers.) 3. Set the display format for f. 4. Enter the two complex points that define the ends of the straight line that your function will be integrated along. The lower limit should be in the Y-registers; the upper limit should be in the X-registers. 5.
PAGE 76
Section 3: Calculating in Complex Mode Keystrokes ' ® ÷ |n Display 054055056057- 12 34 10 43 32 Calculates eiz. Calculates f(z). Approximate the complex integral by integrating the function from 1 + 0i to 1 + 6i using a i2 display format to obtain three significant digits. (The integral beyond 1 + 6i doesn't affect the first three digits.) Keystrokes Display |¥ ´i2 1v 1.00 00 1v6 ´V 6 1.00 00 -3.24 -01 ´A ´% (hold) ® ´% (hold) ´•4 3.82 7.87 1.23 0.0008 -01 -04 -03 Run mode.
PAGE 77
Section 3: Calculating in Complex Mode 77 The program listed below is set up to compute the values of yk from evenly spaced values of xk. You must provide a subroutine labeled "B" that places Im(P(z)) in the real X-register. The program uses inputs that specify the step size h, the number of points n along the real axis, and z0 = x0 + iy0, the initial point on the streamline. You must enter n, h, and z0 into the Z-, Y-, and X-registers before running the program.
PAGE 78
Section 3: Calculating in Complex Mode Keystrokes Display O+V 027-44,40,25 l4 lV * l2 + O6 l3 v 028029030031032033034035- ´_3 t4 1 O-V 4 O÷4 O*V t2 ´b4 l6 ´© ´UOA ´U ) ´© O3 ´UOA ´U t2 036-42,10, 3 03722 4 0381 039-44,30,25 0404 041-44,10, 4 042-44,20,25 04322 2 044-42,21, 4 04545 6 04642 31 047u 44 11 t0 ´b3 l6 ® ´V GB 05322 0 054-42,21, 3 05545 6 05634 05742 25 05832 12 45 4 45 25 20 45 2 40 44 6 45 3 36 048049050051u 33 42 31 44 3 44 11 052- 22 2 78 Increments counter k in Index reg
PAGE 79
Section 3: Calculating in Complex Mode Keystrokes l5 |n 79 Display 059060061- 45 5 30 43 32 Recalls c. Calculates Im(P(zk)) – c. Labels used: A, B, 0, 1, 2, 3, and 4. Registers used: R0, R1, R2 (x0), R3 (y0), R4 (h), R5 (c), R6 (xk), and Index register (k). Matrix used: A.
PAGE 80
Section 3: Calculating in Complex Mode Keystrokes Display |n 067- 43 32 Determine the streamline using z0 = −2 + 0.1 i, step size h = 0.5, and number of points n = 9. Keystrokes |¥ 9v .5v 2”v .1´V ´A Display Run mode. Enters n. Enters h. 9.0000 0.5000 -2.0000 -2.0000 -1.5000 0.1343 Enters z0. x1. y1. ⋮ x9. y9. Descriptor for matrix A. Deactivates Complex mode. ⋮ |"8 2.0000 0.1000 A 9 A 9 2 2 Matrix A contains the following values of xk and yk: xk yk −2.0 −1.5 −1.0 −0.5 0.0 0.5 1.0 1.
PAGE 81
Section 3: Calculating in Complex Mode 81 Example: For the same potential as the previous example, P(z) = 1/z + z, compute the velocity equipotential line starting at z = 2 + i and proceeding to the left. First change subroutine "B" so that it returns Re(P(z))—that is, remove the } instruction from "B". Try n = 6 and h = −0.5.
PAGE 82
Section 4: Using Matrix Operations Matrix algebra is a powerful tool. It allows you to more easily formulate and solve many complicated problems, simplifying otherwise intricate computations. In this section you will find information about how the HP-15C performs certain matrix operations and about using matrix operations in your applications. Several results from numerical linear algebra theory are summarized in this section. This material is not meant to be self-contained.
PAGE 83
Section 4: Using Matrix Operations 83 Row interchanges can also reduce rounding errors that can occur during the calculation of the decomposition. The HP-15C uses the Doolittle method with extended-precision arithmetic to construct the LU decomposition. It generates the decomposition entirely within the result matrix. The LU decomposition is stored in the form U  L    It is not necessary to save the diagonal elements of L since they are always equal to 1.
PAGE 84
Section 4: Using Matrix Operations ILL-Conditioned Matrices and the Condition Number In order to discuss errors in matrix calculations, it's useful to define a measure of distance between two matrices. One measure of the distance between matrices A and B is the norm of their difference, denoted ||A − B||. The norm can also be used to define the condition number of a matrix, which indicates how the relative error of a calculation compares to the relative error of the matrix itself.
PAGE 85
Section 4: Using Matrix Operations 85 that X and B are nonzero vectors satisfying AX = B for some square matrix A. Suppose A is perturbed by ΔA and we compute B + ΔB = (A + ΔA)X. Then  ΔB  ΔA   K (A) , A B with equality for some perturbation ΔA. This measures how much the relative uncertainty in A can be magnified when propagated into the product.
PAGE 86
Section 4: Using Matrix Operations 1 K (A)  min  A  S and A 1 A 1  min  A  S  , where the minimum is taken over all singular matrices S. That is, if K(A) is large, then the relative difference between A and the closest singular matrix S is small. If the norm of A−1 is large, the difference between A and the closest singular matrix S is small. For example, let 1 1 A  1 .9999999999 Then  9,999,999,999 1010  A 1    1010  1010   and ||A−1|| = 2 × 1010.
PAGE 87
Section 4: Using Matrix Operations 87 In the left diagram, ||ΔA|| < 1/||A−1||. If ||ΔA|| << 1 / ||ΔA−1|| (or K(A) ||ΔA||/||A|| << 1), then relative variation in A−1 = ||change in A−1||/||A−1|| ≈ (||ΔA||/||A||) K(A) = ||ΔA||/(1/||A−1||) = (radius of sphere)/(distance to surface) In the right diagram, ||ΔA||>1/||A−1||. In this case, there exists a singular matrix that is indistinguishable from A, and it may not even be reasonable to try to compute the inverse of A.
PAGE 88
Section 4: Using Matrix Operations    number of correct   number of        log A A 1  log(10n)  decimal digits   digits carried  where n is the dimension of A. For the HP-15C, which carries 10 accurate digits, (number of correct decimal digits) ≥ 9 − log(||A|| ||A−1||) − log(10n). In many applications, this accuracy may be adequate. When additional accuracy is desired, the computed solution Z can usually be improved by iterative refinement (also known as residual correction).
PAGE 89
Section 4: Using Matrix Operations 89 For example, let matrix A be 3  10 40 1  A 1 1  2 1  2  1 .  1  The HP-15C correctly calculates A−1 to 10-digit accuracy as A 1 3  1  2    3  4 2 .   1 2  1 Now let  10 20  LR 0  0  0 10  20 0 0   0  10  20  so that 3 E  1 2 1 10 40 10 40 2  10 40  .  10 40  E is very near a singular matrix  3 1 2 S  1 0 0 2 0 0 and ||E – S|| / ||E|| = ⅓ × 10–40.
PAGE 90
Section 4: Using Matrix Operations Multiplying the calculated inverse and the original matrix verifies that the calculated inverse is poor. The trouble is that E is badly scaled. A well-scaled matrix, like A, has all its rows and columns comparable in norm and the same must hold true for its inverse. The rows and columns of E are about as comparable in norm as those of A, but the first row and column of E−1 are small in norm compared with the others.
PAGE 91
Section 4: Using Matrix Operations 91 If (LER) −1 is verifiably poor, you can repeat the scaling, using LER in place of E and using new scaling matrices suggested by LER and the calculated (LER)−1. You can also apply scaling to solving a system of equations, for example EX = B, where E is poorly scaled. When solving for X, replace the system EX = B by a system (LER)Y = LB to be solved for Y. The diagonal scaling matrices L and R are chosen as before to make the matrix LER well-scaled.
PAGE 92
Section 4: Using Matrix Operations 2014.6 1 2014.6 1    X  2014.6 and E -1  2014.61    2014.6 1 2014.6 1 1 1 1 1 1 1 1 1 1 1 1 1 .  1 1 1 1 1 1 1 1 Substituting, you find 1.00146  0.00146   EX  0.00146 .   0.00146 0.00147 Upon checking (using >7), you find that 1/||E−1|| ≈ 9.9 × 10−5 which is very small compared with ||E|| ≈ 1.6 × 104 (or that the calculated condition number is large— ||E|| ||E−1|| ≈ 1.6 × 108).
PAGE 93
Section 4: Using Matrix Operations 93 Note that rTE was scaled by 107 so that each row of E and A has roughly the same norm as every other. Using this new system, the HP-15C calculates the solution  10 7  2000.000080  1999.999980  5    10    X  1999.999980  , with AX   9  10 6  .      0  1999.999980   0  1999.999980    This solution differs from the earlier solution and is correct to 10 digits.
PAGE 94
Section 4: Using Matrix Operations Wr 2 F n   wi2 ri 2 i 1 where W is a diagonal n × n matrix with positive diagonal elements w1, w2, ... , wn. Then Wr 2 F  (y  Xb )T WT W(y  XB) and any solution b also satisfies the weighted normal equations XTWTWXb = XTWTWy. These are the normal equations with X and y replaced by WX and Wy. Consequentially, these equations are sensitive to rounding errors also.
PAGE 95
Section 4: Using Matrix Operations 95 and  10,000.03  XT y   .  9,999.97 However, when rounded to 10 digits,  1010 XT X   10  10  1010  , 1010  which is the same as what would be calculated if X were rounded to five significant digits to the largest element: 100,000  100,000  0  0  . X  0  0   0  0  The HP-15C solves XTXb = XTy (perturbing the singular matrix as described on page 99) and gets  0.060001 b  0.060000 With 0.03 X T y  X T Xb   . 0.
PAGE 96
Section 4: Using Matrix Operations Any n × p matrix X can be factored as X = QTU, where Q is an n × n orthogonal matrix characterized by QT = Q−l and U is an n × p upper-triangular matrix. The essential property of orthogonal matrices is that they preserve length in the sense that Qr 2 F  (Qr) T (Qr )  r T Q T Qr  rT r  r 2 F . Therefore, if r = y – Xb, it has the same length as Qr = Qy – QXb = Qy – Ub.
PAGE 97
Section 4: Using Matrix Operations 97 Only the first p + 1 rows (and columns) of V need to be retained. (Note that Q here is not the same as that mentioned earlier, since this Q must also transform y.) 2. Solve the following system for b: ˆ g  b   0  U       .  0 q   1  q  (If q = 0, replace it by any small nonzero number, say 10−99.) The −1 in the solution matrix automatically appears; it requires no additional calculations.
PAGE 98
Section 4: Using Matrix Operations columns in the original data and refactor the weighted constraint equations. Repeat this procedure if necessary. For example, if the factorization of wC gives 0.3 1.0 2.0 0.5  1.5  U   0 0.02 0.5 3.0 0.1 ,  0 0 2.5 1.5  1.2 then the second diagonal element is much smaller than the value 2.0 above it. This indicates that the first and second columns in the original constraints are nearly dependent.
PAGE 99
Section 4: Using Matrix Operations 99 1 0  3 3   , LU     10  .3333333333 1 0 10  which is nonsingular. The singular matrix B can't be distinguished from the nonsingular matrix 3 3  D  .9999999999 1 since they both have identical calculated LU decompositions. On the other hand, the matrix 3 3  3   1 0 3 A  1    10   LU 1 .9999999999  3 1 0  10  is nonsingular. Using 10-digit accuracy, matrix A is decomposed as 1 0 3 3  LU    . .
PAGE 100
Section 4: Using Matrix Operations matrix A + ΔA represented by the LU decomposition. It can be zero only if the product's magnitude becomes smaller than 10−99 (underflow). Applications The following programs illustrate how you can use matrix operations to solve many types of advanced problems. Constructing an Identity Matrix This program creates an identity matrix In in the matrix whose descriptor is in the Index register. The program assumes that the matrix is already dimensioned to n × n.
PAGE 101
Section 4: Using Matrix Operations 101 One-Step Residual Correction The following program solves the system of equations AX = B for X, then performs one stage iterative refinement to improve the solution. The program uses four matrices: Matrix A B Input System Matrix Right-Hand Matrix Output System Matrix Corrected Solution Keystrokes C D Uncorrected Solution LU Form of A Display |¥ ´CLEARM ´bA l>A O>Á l>B l>Á ´6 l>Á ÷ l>C + |n |¥ Program mode.
PAGE 102
Section 4: Using Matrix Operations Example: Use the residual correction program to calculate the inverse of matrix A for 16 72  33  A   24  10  57 .   8  4  17 The theoretical inverse of A is  29 / 3  8 / 3  32  A   8  5 / 2 51 / 2 .  8 / 3 2/3 9  1 Find the inverse by solving AX = B for X, where B is a 3 × 3 identity matrix. First, enter the program from above.
PAGE 103
Section 4: Using Matrix Operations 103 Let  x1   f 1 ( x)   F11 (x)  F1 p (x)  x   f ( x)    2 2  , and F(x)   F21 (x)  F2 p (x)  , x    , f ( x)                  x p   f p (x)  F p1 (x)  F pp (x) where Fij (x)   f i ( x) x j for i, j = 1, 2, …, p. The system of equations can be expressed as f(x) = 0.
PAGE 104
Section 4: Using Matrix Operations x2 f 2 (x)   ( w / 2) m exp( w / 2)dw  2(1  a)(m  1) . x1 Here x2 > xl > 0, a and n are known (n > 1), and m = (n − 1)/2 − l. An initial guess for (xl, x2) is x1( 0)  xn21, a / 2 and x2( 0)  xn21, 1a / 2 Where x d2, p is the pth percentile of the chi-square distribution with d degrees of freedom. For this example, 1  n  1 x1  F ( x)   m  x 2 2 exp  x1 2 (n  1) x2  1 x 2 2 m  .
PAGE 105
Section 4: Using Matrix Operations Keystrokes 105 Display ´UO4 023- 44 ´UlA ´U O5 024u 45 11 Skips next line for last element.
PAGE 106
Section 4: Using Matrix Operations Keystrokes l2 3 2 ÷ ´! * + OB l4 GC ” ´UOC ´U l5 GC ´UOC ´U |n |n ´bC 2 ÷ ” ' |K ” l2 3 2 ÷ Y * |n Display 059 060061062063064065066067068069 070071u 45 072073 074u 45 5 32 13 44 13 42 44 45 32 44 2 3 30 2 10 0 20 40 12 4 13 16 13 07543 32 07643 32 077-42,21,13 0782 07910 08016 08112 08243 36 08316 08445 2 0853 08630 0872 08810 08914 09020 091 43 32 Calculates m. Calculates Γ(m + 1). Calculates f2. Stores f2 in B. Calculates F21. Stores F21 in C.
PAGE 107
Section 4: Using Matrix Operations 107 Now run the program. For example, choose the values n = 11 and a = 0.05. The suggested initial guesses are x1(0) = 3.25 and x2(0) = 20.5. Remember that the display format affects the uncertainty of the integral calculation. Keystrokes Display |¥ 5´m% 11O2 .05O3 2v1 ´mA ´U ´>1 3.25OA 5.0000 11.0000 0.0500 1 1.0000 1.0000 1.0000 3.2500 20.5OA 20.5000 ´i4 A ) lA 2.0500 1.1677 1.0980 3.5519 01 00 00 00 lA 2.1556 01 Run mode. Reserves R0 through R5.
PAGE 108
Section 4: Using Matrix Operations Solving a Large System of Complex Equations Example: Find the output voltage at a radian frequency of ω = 15 × 103 rad/s for the filter network shown below.
PAGE 109
Section 4: Using Matrix Operations Keystrokes ´mA ´>1 ´U 100OA 150v 800v3÷ -OA ⋮ 150v 8v3÷ -OA l>A ´p ´>2 O< ´⁄ 109 Display 8.0000 Dimensions matrix A to 4 × 8. 8.0000 8.0000 100.0000 150.0000 266.6667 -116.6667 150.0000 2.6667 147.3333 A 4 A 8 A 8 A 8 A 8 Activates User mode. Stores Re(a11). Stores Im(a11). Stores Im(a44). 8 4 8 8 8 Transforms AC to AP. ~ Transforms AP to A. ~ Calculates inverse of A in A. Delete the second half of the rows of A to provide space to store the right-hand matrix B.
PAGE 110
Section 4: Using Matrix Operations Keystrokes Display 1v8 ´mC 8 8.0000 l< ´>4 |c C C C Redimensions matrix C to 1 × 8. 1 8 4 8 1 2 Calculates transpose. Transforms CP to CC. Matrix C contains the desired values of I1, I2, I3, and I4 in rectangular form. Their phasor forms are easy to compute: Keystrokes Display ´>1 ´i4 lC lC ®|: ® lC lC ®|: ® lC lC ®|: ® lC lC ®|: ® ®“5* ´•4 ´U C 4 C 4 1.9550 4.0964 4.1013 8.7212 -1.4489 -3.5633 3.5662 -9.2328 -1.4541 -3.5633 3.5662 -9.2337 5.3446 -2.
PAGE 111
Section 4: Using Matrix Operations 111 Least-Squares Using Normal Equations The unconstrained least-squares problem is known in statistical literature as multiple linear regression. It uses the linear model p y  bj x j  r j 1 Here, b1, …, bp are the unknown parameters, xl, ..., xp are the independent (or explanatory) variables, y is the dependent (or response) variable, and r is the random error having expected value E(r) = 0, variance σ2. After making n observations of y and x1, x2, ...
PAGE 112
Section 4: Using Matrix Operations and   E || y  Xbˆ || 2F (n  p)   2 . For b ≠ 0. When the simpler model y = r is correct, both of these expectations equal σ2. You can test the hypothesis that the simpler model is correct (against the alternative that the original model is correct) by calculating the F ratio || Xbˆ || 2F p F || y  Xbˆ || 2F (n  p) F will tend to be larger when the original model is true (b ≠ 0) than when the simpler model is true (b = 0).
PAGE 113
Section 4: Using Matrix Operations 113 The analysis of variance (ANOVA) table below partitions the total sum of squares (Tot SS) into the regression and the residual sums of squares. You can use the table to calculate the F ratio.
PAGE 114
Section 4: Using Matrix Operations You will need the original dependent variable data for each regression. If there is not enough room to store the original data in matrix E, you can compute it from the output of the regression fit. A subroutine has been included to do this. This program has the following characteristics:  If the entire program is keyed into program memory, the sizes of n and p are required to satisfy n ≥ p and (n + p)(p + 1) ≤ 56.
PAGE 115
Section 4: Using Matrix Operations Keystrokes 115 Display ´>6 ´>8 |x lmA ÷ v v l>C ´B ´>8 |x |K |n ´bB 018-42,16, 6 019-42,16, 8 02043 11 021-45,23,11 022 30 02310 02436 02536 026-45,16,13 027-42,26,13 02810 02943 33 030-45,16,12 031-42,16, 8 03243 11 03330 03443 36 03543 32 036-42,21,12 l>A l>Á ” ´6 l>Á ” |n 037-45,16,11 038-45,16,14 03916 040-42,26,12 041-42,16, 6 042-45,16,14 04316 04443 32 Calculates residuals of fit in B. Calculates Res SS. Calculates σ2 estimate.
PAGE 116
Section 4: Using Matrix Operations 4. Press ´> 1 to set registers R0 and R1. 5. Press ´U to activate User mode. 6. For each observation, store the values of the p variables in a row of matrix A. Repeat this for the n observations. 7. Store the values of the dependent variable in matrix B. 8. Press A to calculate and display the Res SS. The Y-register contains the Reg SS and the T-register contains the σ2 estimate. 9. Press lÁ to observe each of the p parameter estimates. 10.
PAGE 117
Section 4: Using Matrix Operations ´>1 ´U 1OA 3.9OA 3.5OA 1.0000 1.0000 1.0000 3.9000 3.5000 ⋮ ⋮ 1OA 19.3OA 5.8OA 5.4OB 5.9OB ⋮ 11.5OB A´•9 ) )) lÁ lÁ lÁ B l>A ´>4 2v11 ´mA l>A ´>4 A ) )) lÁ lÁ B l>A ´>4 1v11 ´mA l>A ´>4 A ) )) 1.0000 19.3000 5.8000 5.4000 5.9000 Enters independent variable data. Enters dependent variable data. ⋮ 11.5000 13.51217504 587.9878252 1.689021880 1.245864326 0.379758235 0.413552218 d 3 1 A 11 3 A 3 11 11 11.00000000 A 2 11 A 11 2 16.78680552 584.7131947 1.865200613 3.
PAGE 118
Section 4: Using Matrix Operations lÁ ´U ´•4 6.963636364 6.963636364 6.9636 b1 estimate. Deactivates User mode. The Reg SS for the PPI variable adjusted for the constant term is (Reg SS for reduced model) − (Reg SS for constant) = 51.29864900. The Reg SS for the UR variable adjusted for the PPI variable and the constant term is (Reg SS for full model) – (Reg SS for reduced model) = 3.274630500.
PAGE 119
Section 4: Using Matrix Operations 119 ˆ g  ( p rows) U   V   0 q  (1 row ),  0 0  (n  p  1 rows)   and Û is an upper-triangular matrix. If this factorization results from including n rows rm = (xm1, xm2, …, xmp, ym ) for m = 1, 2, ... , n in [X y], consider how to advance to n + 1 rows by appending row rn+1 to[X y]: X y  Q T r    n 1   0 0  V   . 1 rn 1  The zero rows of V are discarded.
PAGE 120
Section 4: Using Matrix Operations U *  A*   0 0  g *  ( p rows)  q*  (1 row ), 0  (1 row ) Where U* is also an upper-triangular matrix. You can obtain the solution b augmented system of p + 1 rows by solving U *  0 (n+1) to the g *  b ( n1)   0      * q *    1   q  By replacing the last row of A* by rn+2 and repeating the factorization, you can continue including additional rows of data in the system.
PAGE 121
Section 4: Using Matrix Operations Keystrokes 121 Display lmA 018-45,23,11 ® O2 ´>1 ´b1 |"0 l2 l0 l|A lA |T2 |F0 01934 02044 2 021-42,16, 1 022-42,21, 1 023-43, 5, 0 02445 2 02545 0 026-45,43,11 02745 11 028-43,30, 2 029-43, 4, 0 |a |: |` 1 ´; 030031032033034- |?0 ” ´V ) ´b2 |( lA l2 l1 l|A ´V * l2 l1 O|A ´} ´UOA ´U l1 035-43, 6, 0 03616 03742 25 03833 039-42,21, 2 04043 33 04145 11 04245 2 04345 1 044-45,43,11 04542 25 04620 04745 2 04845 1 049-44,43,11 05042 30 051u 44 11 052- 43 16 43 1 43 35
PAGE 122
Section 4: Using Matrix Operations Keystrokes Display l0 |£ t2 053054055- |"8 O1 l2 |£ |n t1 ´bC 056-43 5 8 05744 1 05845 2 05943 10 06043 32 06122 1 062-42,21,13 lmA v ´mA O0 O1 1 ´mC 063-45,23,11 06436 065-42,23,11 06644 0 06744 1 0681 069-42,23,13 0 O>C “ 9 9 ” lA |~ ) ” l0 1 O|C l>C l>A ´
PAGE 123
Section 4: Using Matrix Operations Keystrokes 123 Display + l0 ´mA 08940 09045 0 091-42,23,11 1 1 ´mC lA ´>1 |n 0921 09330 0941 095-42,23,13 09645 11 097-42,16, 1 09843 32 Dimensions matrix A as (p + 2) × (p + 1) Dimensions matrix C as p × 1. Recalls q. Sets k = l = 1. Labels used: A, B, C, and 1 through 5. Registers used: R0, R1, and R2 (p+2 and w) Matrices used: A (working matrix) and C (parameter estimates). Flags used: 0 and 8.
PAGE 124
Section 4: Using Matrix Operations 8. Optionally, press C|x to calculate and display the residual sum of squares q2 and to calculate the current solution b. Then press lCp times to display b1, b2, …, bp in turn. 9. Repeat steps 5 through 8 for each additional row. Example: Use this program and the CPI data from the previous example to fit the model y = bl + b2x2 + b3x3 + r, where y, x2, and x3 represent the CPI, PPI, and UR (all as percentages).
PAGE 125
Section 4: Using Matrix Operations 125 These estimates agree (to within 3 in the ninth significant digit) with the results of the preceding example, which uses the normal equations. In addition, you can include additional data and update the parameter estimates. For example, add this data from 1968: CPI = 4.2, PPI = 2.5 and UR = 3.6. Keystrokes Display 1A 1¦ 2.5¦ 3.6¦ 4.2¦ B C |x 1.000000000 2.000000000 3.000000000 4.000000000 1.000000000 1.000000000 3.700256908 13.691900119 lC 1.
PAGE 126
Section 4: Using Matrix Operations An orthogonal change of variables x = Qz, which is equivalent to rotating the coordinate axes, changes the equation of a family of quadratic surfaces (xTAx = constant) into the form k z T (Q T AQ )z   λ j z 2j  constant. j With the equation in this form, you can recognize what kind of surfaces these are (ellipsoids, hyperboloids, paraboloids, cones, cylinders, planes) because the surface's semi-axes lie along the new coordinate axes.
PAGE 127
Section 4: Using Matrix Operations Keystrokes |¥ ´CLEARM ´bA l>A O>B O>C ´>4 l>B O< + 2 ÷ O>A ´>8 O2 |` O3 O>C ´>1 ´b0 l0 l1 |T5 t3 |T7 t1 ® l|B ” ´UOB ´U t0 ´b1 l|A |a O+3 |K v + 127 Display Program mode.
PAGE 128
Section 4: Using Matrix Operations Keystrokes Display l0 v l|A l1 v l|A |T3 t2 ” 03745 0 03836 039-45,43,11 04045 1 04136 042-45,43,11 04330 044-43,30, 3 04522 2 04616 ® ” ® ´b2 |: |` 4 ÷ ] ´UOB ´U t0 ´b3 1 OC ´UOB ´U t0 l3 l÷2 ¦ 2 l>B ÷ l>C l>A ´
PAGE 129
Section 4: Using Matrix Operations Keystrokes 129 Display ´>5 l>B ´
PAGE 130
Section 4: Using Matrix Operations Keystrokes Display ⋮ 3OA 4OA A ¦ ¦ ¦ ¦ ¦ lA lA lA lA lA lA lA lA lA ´U 3.0000 4.0000 0.8660 0.2304 0.1039 0.0060 3.0463 5.8257 -0.8730 -9.0006 -2.0637 -9.0006 9.3429 1.0725 -2.0637 1.0725 6.8730 6.8730 Enters a32. Enters a33. Calculates ratio－too large. Again, too large. Again, too large. Again, too large. Again, too large. Negligible ratio. Recalls a11=λ1. Recalls a12. Recalls a13. Recalls a21. Recalls a22=λ2. Recalls a23. Recalls a31. Recalls a32.
PAGE 131
Section 4: Using Matrix Operations 131 The value used for λk need not be exact; the calculated eigenvector is determined accurately in spite of small inaccuracies in λk. Furthermore, don't be concerned about having too accurate an approximation to λk; the HP-15C can calculate the eigenvector even when A − λk I is very ill-conditioned. This technique requires that vector z(0) have a nonzero component along the unknown eigenvector qk..
PAGE 132
Section 4: Using Matrix Operations Keystrokes Display l>C O>Á O< l>B ÷ v ´>7 ÷ ´>1 ´b7 024-45,16,13 025-44,16,14 02644 26 027-45,16,12 02810 02936 030-42,16, 7 03110 032-42,16, 1 033-42,21, 7 ´UlC ´U v 034u 45 13 035- 36 |a 1 |T6 t7 l>C |K ÷ l>Á O< ´>7 ´>1 ¦ 03643 16 0371 038-43,30, 6 03922 7 040-45,16,13 04143 36 04210 043-45,16,14 04444 26 04530 046-42,16, 7 047-42,16, 1 04831 t6 049- 22 Stores z(n) in D. Calculates w(n+1) in C. Calculates ±z(n+1) in C.
PAGE 133
Section 4: Using Matrix Operations 133 3. Dimension and enter the elements into matrix A using ´mA, ´>1, and OA. 4. Key in the eigenvalue and press C. The display shows the correction parameter ||z(1) − z(0)||R. 5. Press ¦ repeatedly until the correction parameter is negligibly small. 6. Press lC repeatedly to view the components of qk, the eigenvector. Example: For matrix A of the previous example, 0 A 1  2 1 2 3 2 3  4  Calculate the eigenvectors q1, q2, and q3.
PAGE 134
Section 4: Using Matrix Operations Keystrokes Display lC lC lC 6.8730C -0.5000 1.0000 -0.5000 0.7371 ¦ ¦ ¦ lC lC lC ´U 1.9372 1.0000 0.0000 0.3923 0.6961 1.0000 1.0000 Eigenvector forλ2. Uses λ3=6.8730 (approximation). -06 -10 Eigenvector forλ3. Deactivates User mode. If matrix A is no larger than 3×3, this program can be included with the previous eigenvalue program.
PAGE 135
Section 4: Using Matrix Operations 135 Labels used: E and 8. Registers used: no additional registers. Matrices used: A (from previous program) and E (eigenvalues). To use the combined eigenvalue, eigenvalue storage, and eigenvector programs for an A matrix up to 3×3: 1. Execute the eigenvalue program as described earlier. 2. Press E to store the eigenvalues. 3. Enter again the elements of the original matrix into A. 4. Recall the desired eigenvalue from matrix E using l E. 5.
PAGE 136
Section 4: Using Matrix Operations s   f (x) for finding a minimum  f (x) for finding a maximum.  Once the direction is determined from the gradient, the program looks for the optimum distance to move from xj in the direction indicated by sj—the distance that gives the greatest improvement in f(x) toward a minimum or maximum To do this, the program finds the optimum value tj by calculating the slope of the function gj(t) = f(xj + tsj) at increasing values of t until the slope changes sign.
PAGE 137
Section 4: Using Matrix Operations N 137 Determines the maximum number of iterations that the program will attempt in each of two procedures: the bounding search and the overall optimization procedure. That is, the program halts if the bounding search finds no change of sign within N iterations. Also, the program halts if the norm of the gradient is still too large at xN. Each of these situations results in an Error 1 display. (They can be distinguished by pressing −.
PAGE 138
Section 4: Using Matrix Operations Keystrokes Display |n ´b7 l4 l÷6 O8 GE l>E O>Á l>Á |?0 ” 00843 32 009-42,21, 7 01045 4 011-45,10, 6 01244 8 01332 15 014-45,16,15 015-44,16,14 016-45,16,14 017-43, 6, 0 01816 ´>8 |~ |n ⁄ l*8 O.1 0 O.0 l5 O7 ´b6 l.1 G3 ´© |?0 ” |T4 t5 G8 l.1 O.0 l2 O*.1 ´s7 t6 l>A |a 019-42,16, 8 02043 20 02143 32 02215 023-45,20, 8 02444 .1 0250 02644 .0 02745 5 02844 7 029-42,21, 6 03045 .1 03132 3 03242 31 033-43, 6, 0 03416 035-43,30, 4 03622 5 03732 8 03845 .1 03944 .
PAGE 139
Section 4: Using Matrix Operations Keystrokes 139 Display t6 ´b5 l6 O7 ´b4 G8 l.0 l+.1 2 ÷ O8 G3 |?0 ” 1 1 OV ) |T1 ´eV l8 O% ´e7 t4 |n ´b3 l>Á ´A + G8 04622 6 047-42,21, 5 04845 6 04944 7 050-42,21, 4 05132 8 05245 .0 053-45,40,.
PAGE 140
Section 4: Using Matrix Operations Keystrokes Display ´>5 1 v l|B |n ´bA 0 O6 ´b2 1 O+6 ´i3 G7 l6 ´•0 ´© ´>1 ´i3 l9 ¦ l3 l>E ´>8 |£ tB ´© l5 l6 |T8 t2 l>C |a 083-42,16, 5 0841 08536 086-45,43,12 08743 32 088-42,21,11 0890 09044 6 091-42,21, 2 0921 093-44,40, 6 094-42, 8, 3 09532 7 09645 6 097-42, 7, 0 09842 31 099-42,16, 1 100-42, 8, 3 10145 9 10231 10345 3 104-45,16,15 105-42,16, 8 10643 10 10722 12 10842 31 10945 5 11045 6 111-43,30, 8 11222 2 113-45,16,13 11443 16 t2 ´bB |F9 ¦ tB 11522 2 116-4
PAGE 141
Section 4: Using Matrix Operations 141 Labels used: A, B, and 2 through 8. Registers used: R2 through R9, R.0, R.1, and Index register. Matrices used: A, B, C, D, and E. Your subroutine, labeled "E", may use any labels and registers not listed above, plus the Index register, matrix B, and matrix E (which should contain your calculated gradient). To use the program: 1. Enter your subroutine into program memory. 2. Press 11 ´m% to reserve registers R0 through R.
PAGE 142
Section 4: Using Matrix Operations Example: Use the optimization program to find the dimensions of the box of largest volume with the sum of the length and girth (perimeter of cross section) equaling 100 centimeters. For this problem l + (2h + 2w) = 100 v = whl v(w,h) = wh(100－2h－2w) = 100wh － 2wh2 － 2hw2  2h(50  h  2w)  v( w, h)   . 2w(50  w  2h) The solution should satisfy w + h < 50, w > 0, and h > 0. First enter a subroutine to calculate the gradient and the volume.
PAGE 143
Section 4: Using Matrix Operations Keystrokes 143 Display l.3 l+.3 - 14445 .3 145-45,40,.3 14630 l.2 l*.3 |n 14745 .2 148-45,20,.3 14943 32 Replaces ei with lei − 2wh, the gradient elements. Calculates lwh. Now enter the necessary information and run the program. Keystrokes |¥ 13´m% |"0 ´U ´>1 2v1 ´mA 15OA OA 3O2 0.1O3 0.05O4 4O5 A − Display Run mode. Reserves R0 through R.3. Finds local maximum. Activates User mode. 13.0000 13.0000 13.0000 13.0000 1 1.0000 15.0000 15.0000 3.0000 0.1000 0.
PAGE 144
Section 4: Using Matrix Operations Keystrokes ¦ ¦ ¦ − ´•4 lA lA ´U ´>0 Display 9.253 3.480 1.121 9.431 4.126 -1.139 2. 9.259 5.479 -6.127 3. 9.259 7.726 7.726 0.0773 16.6661 16.6661 16.6661 16.6661 03 01 03 02 02 03 03 -01 -01 03 -02 -02 Volume at this iteration. Gradient. Slope at u1. Slope at u2. Slope at u3. Slope at u4 (sign change). j + 1. Volume at this iteration. Gradient. Slope at u1 (sign change). j + 1. Volume at this iteration. Gradient less than e. Stops blinking. Recalls h form a1.
PAGE 145
Appendix: Accuracy of Numerical Calculations Misconceptions About Errors Error is not sin, nor is it always a mistake. Numerical error is merely the difference between what you wish to calculate and what you get. The difference matters only if it is too big. Usually it is negligible; but sometimes error is distressingly big, hard to explain, and harder to correct. This appendix focuses on errors, especially those that might be large—however rare. Here are some examples. Example 1: A Broken Calculator.
PAGE 146
Appendix: Accuracy of Numerical Calculations where payment = $0.01 = one penny per second, i = 0.1125 = 11.25 percent per annum interest rate, n = 60×60×24×365 = number of seconds in a year. Using her HP-15C, Susan reckons that the total will be $376,877.67. But at year's end the bank account is found to hold $333,783.35. Is Susan entitled to the $43,094.32 difference? In both examples the discrepancies are caused by rounding errors that could have been avoided. This appendix explains how.
PAGE 147
Appendix: Accuracy of Numerical Calculations 147 needed to line up the bridge girders before she can rivet them together. Both see the same discrepancy d, but what looks negligible to one person can seem awfully big to another. Whether large or small, errors must have sources which, if understood, usually permit us to compensate for the errors or to circumvent them altogether.
PAGE 148
Appendix: Accuracy of Numerical Calculations 0 < x < 1 above. When R(x) is squared 50 times to produce F(x) = S(R(x)), the result is clearly 1 for x ≥ 1, but why is F(x) = 0 for 0 ≤ x < 1? When x <1, sRx   s0.9999999999  1  10 10    250  6.14  10 48898. This value is so small that the calculated value F(x) = S(R(x)) underflows to 0. So the HP-15C isn't broken; it is doing the best that can be done with 10 significant digits of precision and 2 exponent digits.
PAGE 149
Appendix: Accuracy of Numerical Calculations 149 Regarding the first misconception, example 1 would behave in the same perverse way if it suffered only one rounding error, the one that produces R(x) = 1 or 0.9999999999, in error by less than one unit in its last (10th) significant digit. Regarding the second misconception, cancellation is what happens when two nearly equal numbers are subtracted.
PAGE 150
Appendix: Accuracy of Numerical Calculations A Hierarchy of Errors Some errors are easier to explain and to tolerate than others. Therefore, the functions delivered by single keystrokes on the HP-15C have been categorized, for the purposes of easier exposition, according to how difficult their errors are to estimate. The estimates should be regarded as goals set by the calculator's designers rather than as specifications that guarantee some stated level of accuracy.
PAGE 151
Appendix: Accuracy of Numerical Calculations 151 significantly smaller than one unit in the 10th significant digit of the result, include ∆, À,r,d,p,and c; N,o, @,P]for real arguments; :, ,, {, /,H[, H\and H]for real and complex arguments; a, ¤, and⁄for complex arguments; matrix norms >7, >8; and finally [, \, and ] for real arguments in Degrees and Grads modes (but not in Radians mode − refer to Level 2, page 154) A function that grows to ∞ or decays to 0 exponentially fast as its argument approaches ± ∞
PAGE 152
Appendix: Accuracy of Numerical Calculations She calculated $376,877.67 on her HP-15C, but the bank's total was $333,783.35, and this latter total agrees with the results calculated on good, modern financial calculators like the HP-12C, HP-37E, HP-38E/38C, and HP-92. Where did Susan's calculation go awry? No severe cancellation, no vast accumulation of errors; just one rounding error that grew insidiously caused the damage: i/n = 0.000000003567351598 1 + i/n = 1.
PAGE 153
Appendix: Accuracy of Numerical Calculations Keystrokes ÷ * |n |¥ 153 Display 011012013014- Calculates u − 1 when u≠1. Calculates x/(u − 1) or 1/1. Calculates λ(x). 30 10 20 43 32 The calculated value of u, correctly rounded by the HP-15C, is u=(1+ε)(1+x), where |ε| <5×10-10. If u=1, then |x| = |1/(1+ε)-1| ≤ 5×10-10 too, in which case the Taylor series λ(x) = x (1 − ½ x + ⅓ x2 − ... ) tells us that the correctly rounded value of λ(x) must be just x.
PAGE 154
Appendix: Accuracy of Numerical Calculations numbers. In other words, every complex function f in Level 1C will produce a calculated complex value F = (1 + ε) f whose small complex relative error ε must satisfy |ε| < 10-9. The complex functions in Level lC are *,÷,x,N,o,,,{,/, H[, H\, and H]. Therefore, a function like λ(z) = ln(1+z) can be calculated accurately for all z by the same program given above with the same explanation.
PAGE 155
Appendix: Accuracy of Numerical Calculations 155 with [ (1014 $) = 0.7990550814, although the true sin (1014 $) = −0.78387… The wrong sign is an error too serious to ignore; it seems to suggest a defect in the calculator. To understand the error in trigonometric functions we must pay attention to small differences among π and two approximations to π: true π key $ internal p = 3.1415926535897932384626433 ... = 3.141592654 (matches π to 10 digits) = 3.
PAGE 156
Appendix: Accuracy of Numerical Calculations [(2x) = −0.00001100815000 2[(x) \(x) = −0.00001100815000. Note the close agreement even though for this x, sin(2x) = 2sin(x)cos(x) = −0.0000110150176 ... disagrees with [(2x) in its fourth significant digit. The same identities are satisfied by Æ(x) values as by trig(x) values even though Æ(x) and trig(x) may disagree.
PAGE 157
Appendix: Accuracy of Numerical Calculations 157 It all seems like much ado about very little. After a blizzard of formulas and examples, we conclude that the error caused by p ≠ π is negligible for engineering purposes, so we need not have bothered to know about it. That is the burden that conscientious error analysts must bear; if they merely took for granted that small errors are negligible, they might be wrong.
PAGE 158
Appendix: Accuracy of Numerical Calculations Some transformations f are stable in the presence of input noise; they keep Δy relatively small as long as Δx is relatively small. Other transformations f may be unstable in the presence of noise because certain relatively small input noises Δx cause relatively huge perturbations Δy in the output.
PAGE 159
Appendix: Accuracy of Numerical Calculations 159 introduced to explain diverse noise sources actually distributed throughout F. Some of the noise appears as a tolerably small perturbation δx to the input—hence the term "backward error analysis." Such a system F, whose noise can be accounted for by two tolerably small perturbations, is therefore classified into Level 2 for purposes of exposition.
PAGE 160
Appendix: Accuracy of Numerical Calculations Example 6: The Smaller Root of a Quadratic. The two roots x and y of the quadratic equation c − 2bz + az2 = 0 are real whenever d = b2 − ac is nonnegative. Then the root y of smaller magnitude can be regarded as a function y = f(a,b,c) of the quadratic's coefficients (b  d sgn( b)) / a f (a, b, c)    (c / b) / 2 if a  0 otherwise.
PAGE 161
Appendix: Accuracy of Numerical Calculations 161 when we wish to calculate f(x), the best we could hope to get is f(x + Δx), but we actually get F(x+Δx)=(f+δf)(x+Δx+δx), where |δf| ≤ ε|f| and |δx| ≤ η|x|. What we get is scarcely worse than the best we could hope for provided the tolerances ε and η are small enough, particularly if |Δx| is likely to be at least roughly as big as η|x|.
PAGE 162
Appendix: Accuracy of Numerical Calculations all calculated accurately because N is in Level 1. What might happen if N were in Level 2 instead? If N were in Level 2, then "successful" backward error analysis would show that, for arguments u near 1, N (u) = ln(u +δu) with |δu| < 10-9.
PAGE 163
Appendix: Accuracy of Numerical Calculations 163 Example 7: The Angle in a Triangle. The cosine law for triangles says r2 = p2 + q2 – 2pq cos θ for the figure shown below. Engineering and scientific calculations often require that the angle θ be calculated from given values p, q, and r for the length of the triangle's sides.
PAGE 164
Appendix: Accuracy of Numerical Calculations Disparate Results from Three Programs FA, FB, FC Case 1 Case 2 p 1. 9.999999996 q 1. 9.999999994 r 1.00005 × 10-5 3 × 10-9 0. 0. FA -4 FB 5.73072 × 10 FC 5.72986 × 10-4 10. 5.000000001 15. 180. Error 0 1.28117 × 10-8 Case 4 p Case 3 180. 179.9985965 Case 5 Case 6 0.527864055 9.999999996 9.999999999 q 9.472135941 -9 3 × 10 9.999999999 r 9.999999996 9.999999994 FA Error 0 FB Error 0 FC 48.18968509 Error 0 180. 48.
PAGE 165
Appendix: Accuracy of Numerical Calculations Keystrokes * ® |K |x + |( ® v + ÷ |{ |n ´bB O1 v |( O+1 |( O+1 2 O÷1 ) l-1 ® l-1 * ¤ ® l-1 l*1 ” ¤ |: ) * |n ´bC O0 Display 00720 00834 00943 36 01043 11 01140 01243 33 01330 01434 01536 01640 01710 01843 24 01943 32 020-42,21,12 02144 1 02236 02343 33 024-44,40, 1 02543 33 026-44,40, 1 0272 028-44,10, 1 02933 030-45,30, 1 03134 032-45,30, 1 03320 03411 03534 036-45,30, 1 037-45,20, 1 03816 03911 04043 1 04133 04220 04343 32 044-42,21,13 04544 0 165
PAGE 166
Appendix: Accuracy of Numerical Calculations Keystrokes ) |£ ® O1 O+0 ® O+0 |( O-1 |K v l+1 ¤ ´X0 ¤ O*0 |` + ) + ´X1 |( |K |£ t.9 ) |T2 ¤ ® t.8 ´b.9 |T2 ¤ |( ´b.8 ¤ l1 Display 04633 04743 10 04834 04944 1 050-44,40, 0 05134 052-44,40, 0 05330 05443 33 055-44,30, 1 05643 36 05736 058-45,40, 1 05911 060-42, 4, 0 06111 062-44,20, 0 06343 35 06440 06533 06640 067-42, 4, 1 06843 33 06943 36 07043 10 07122 .9 07233 073-43,30, 2 07411 07534 07622 .8 077-42,21,.9 078-43,30, 2 07911 08043 33 081-42,21,.
PAGE 167
Appendix: Accuracy of Numerical Calculations Keystrokes ¤ * l0 |: |~ ÷ ® v + |n |¥ 167 Display 085086087088089090091092093094- 45 43 43 43 11 20 0 1 20 10 34 36 40 32 The results FC(p, q , r) are correct to at least nine significant digits. They are obtained from a program "C" that is utterly reliable though rather longer than the unreliable programs "A" and "B". The method underlying program "C" is: 1. If p < q, then swap them to ensure p ≥ q. 2. Calculate b=(p−q)+r, c=(p−r)+q, and s=(p+r)+q. 3.
PAGE 168
Appendix: Accuracy of Numerical Calculations constitute the edge lengths of a feasible triangle, so FC might produce an error message when it shouldn't, or vice-versa, on those machines. Backward Error Analysis of Matrix Inversion The usual measure of the magnitude of a matrix X is a norm ||X|| such as is calculated by either >7 or >8; we shall use the former norm, the row norm X  max  xij i j in what follows.
PAGE 169
Appendix: Accuracy of Numerical Calculations 169 (A + ΔA)-1 is at least about as far from A-1 in norm as the calculated ⁄ (A). The figure below illustrates the situation. As A + ΔA runs through matrices with ||ΔA|| at least about as large as roundoff in ||A||, its inverse (A + ΔA)-1 must roam at least about as far from A-1 as the distance from A-1 to the computed ⁄ (A).
PAGE 170
Appendix: Accuracy of Numerical Calculations 50,000.03  45  0.00002  50 ,000 0 50,000  50,000.03 45  X  0 0 0.00002  50,000.03   0 0 52,000  0 and p q 50,000 50 ,000    0 0.00002 50 ,000.03 48,076.98077  -1 X   0 0 50 ,000 48,076.95192    0 0 0 0.00001923076923  Ideally, p = q = 0, but the HP-15C's approximation to X-1, namely ⁄ (X), has q = 9.643.269231 instead, a relative error X 1  ⁄ ( X) X 1  0.0964, Nearly 10 percent.
PAGE 171
Appendix: Accuracy of Numerical Calculations 171 Of course, none of these horrible things could happen if X were not so nearly singular. Because ||X|| ||X-1|| > 1010, a change in X amounting to less than one unit in the 10th significant digit of ||X|| could make X singular; such a change might replace one of the diagonal elements 0.00002 of X by zero. Since X is so nearly singular, the accuracy ⁄(X) in this case rather exceeds what might be expected in general.
PAGE 172
Appendix: Accuracy of Numerical Calculations Example 6 Continued. The program listed below solves the real quadratic equation c − 2 bz + az2 = 0 for real or complex roots. To use the program, key the real constants into the stack (c v b v a) and run program "A". The roots x and y will appear in the X- and Y-registers. If the roots are complex, the C annunciator turns on, indicating that Complex mode has been activated. The program uses labels "A" and ".
PAGE 173
Appendix: Accuracy of Numerical Calculations Keystrokes ÷ |K |( ÷ |n ´b.9 ¤ lV |( ÷ ® |K ÷ ´V v ´} ” ´} |n |¥ 173 Display 02110 02243 36 02343 33 02410 02543 32 026-42,21,.9 02711 02845 25 02943 33 03010 03134 03243 36 03310 03442 25 03536 03642 30 03716 03842 30 03943 32 The method uses d = b2 − ac. If d < 0, then the roots are a complex conjugate pair b a  i  d a. If d ≥ 0, then the roots are real numbers x and y calculated by s  b  d sgn (b) x  s/a c / s if s  0 y  0 if s  0.
PAGE 174
Appendix: Accuracy of Numerical Calculations The results of several cases are summarized below. Case 1 Case 2 Case 3 Case 4 c 3 4 1 654,321 b 2 0 1 654,322 a 1 1 Roots Real Complex 10 654,323 Real 0 ± 2i 3 -13 2 × 10 1 0.5 Real 13 0.9999984717 0.9999984717 Case 5 Case 6 c 46,152,709 12,066,163 b 735,246 987,644 a 11,713 80,841 Roots Real Complex 12.21711755 ± i0.001377461 62.77179203 62.
PAGE 175
Appendix: Accuracy of Numerical Calculations 175 Subroutine "B" below, which uses such a trick,* is a very short program that guarantees nine correct significant digits on a 10-digit calculator. It uses labels "B", ".7", and ".8" and registers R0 through R9 and the Index register. To use it, key in c v b v a, run subroutine "B", and wait for results as before. Keystrokes |¥ ´CLEARM ´bB OV ) O0 O8 ® O1 O9 ´i2 ´b.
PAGE 176
Appendix: Accuracy of Numerical Calculations Keystrokes |( O8 l7 O9 |a “ 2 0 * l1 |a |£ t.8 ´bB ´•9 l8 |x O7 lV l9 |w l7 |T2 T .7 ¤ ´X0 |T2 l-0 |T3 l+0 ´X1 |T0 l÷1 l1 l÷V |n ´b.7 ” ¤ Display 03243 33 03344 8 03445 7 03544 9 03643 16 03726 0382 0390 04020 04145 1 04243 16 04343 10 04422 .8 045-42,21,12 046-42, 7, 9 04745 8 04843 11 04944 7 05045 25 05145 9 05243 49 05345 7 054-43,30, 2 05522 .
PAGE 177
Appendix: Accuracy of Numerical Calculations Keystrokes l÷V v ” l0 lV ÷ ® ´V v |( ´V |n |¥ 177 Display 071-45,10,25 07236 07316 07445 0 07545 25 07610 07734 07842 25 07936 08043 33 08142 25 08243 32 This program's accuracy is phenomenal: better than nine significant digits even for the imaginary parts of nearly indistinguishable complex roots (as when c = 4,877,163,849 and b = 4,877,262,613 and a = 4,877,361,379); if the roots are integers, real or complex, and if a = 1, then the roots are calculated ex
PAGE 178
Index Page numbers in bold type indicate primary references; page numbers in regular type indicate secondary references.
PAGE 179
Index Complementary normal distribution function · 51–55 Complex components, accurate · 64 Complex equations, solving large system · 107–110 Complex math functions · 59–62 Complex mode · 56–81 _ and f · 63 accuracy · 63–65 Complex multivalued functions · 59–62 Complex number, n th roots · 59, 67–69 Complex number, storing and recalling · 65–66 Complex potential function · 76–81 Complex relative error · 154 Complex roots of equation · 69–73, 17–18 Complex roots of quadratic equation · 171–177 Complex single
PAGE 180
Index Electrostatic field · 50 Endpoint, f sampling at · 41, 48 Equations complex, solving large system · 107–110 equivalent · 11–12 solving inaccurate · 12 solving nonlinear system · 102–107 with several roots · 12 Equipotential lines · 76–81 Equivalent equations · 11–12 Error · 145 absolute · 146, 153 hierarchy · 150 in matrix elements · 84 misconceptions · 145–149 relative · 146, 151 Error 0 · 27, 164, 167, 173 Error 1 · 137, 141 Error 4 · 26, 27, 35 Error 8 · 11, 22 Error analysis, backward · 157–
PAGE 181
Index Extended precision · 41 Extremes of function · 18–24 F F actorization, orthogonal · 95–98, 118–25 F ratio · 112–18 Field Intensity · 18–24 Financial equation · 26 Financial problems · 24–39 Format, display · 40–41, 42 Frobenius norm · 84 Functions, complex · 59–62 G Gamma function, complex · 56–58 Gradient · 135, 136, 137, 141, 142 Grandmother's expensive chinaware · 177 H Hierarchy of error · 150 Horner's method · 13, 14 Hyperbolic cylinder · 129–30 I Identity matrix · 100 Ill-conditioned matrix · 8
PAGE 182
Index backward error analysis · 168–171 IRR · 34–39 Iterative refinement · 88, 100–102 J Jordon canonical form · 131 L Large system of complex equations, solving · 107–110 Least-squares · 93–98, 110–25, 157 linearly constrained · 94, 120 weighted · 93, 94, 97, 98, 120 Least-Squares linearly constrained · 98 Level 0 · 150 Level 1 · 150–53, 159, 162 Level 1C · 153–54 Level 2 · 154–177 Level ∞ · 150 Line search · 136 Linear model · 111 Linear regression, multiple .
PAGE 183
Index Multivalued functions, complex · 59–62 N n th roots of complex number · 59, 67–69 Nearly singular matrix · 88, 93 Net present value · 34–39 equation · 34 Network, filter · 107–10 Newton's iteration method · 69–70, 103 Noise, input and output · 157–61 Nonlinear equations, solving system · 102–7 Nonsingular matrix · 85–86, 99 Norm · 84, 85, 86, 90, 91 Normal distribution · 42, 103, 112 Normal distribution function · 42, 51–55 complementary · 51–55 Normal equations · 93–94, 110–18 augmented · 93–94 weig
PAGE 184
Index Preconditioning a system · 91–93 Present value · 24–39 Principal branch · 59–62 Principal value · 59–62 Q Quadratic equation, roots · 160, 171–77 Quadratic surface · 126, 129–30 R Radians, in complex mode · 58 Rate of return · 34–39 Recalling complex numbers · 65–66 Rectangular form · 59 Refinement, iterative · 88 Regression sum of squares · 111–18 mean-adjusted · 113 Regression, multiple linear .
PAGE 185
Index Sampling, f · 41, 48, 63 Scaling a matrix · 88–91 Scaling a system · 90–91 Secant method · 9 Sign change · 10 Sign symmetry · 151, 156 Single-valued functions, complex · 59 Singular matrix · 85–86, 98–100, 168 Singularity and backward error analysis · 161–62 Skew-symmetric matrix · 126 Slope · 20–22 Smaller root of quadratic equation · 160, 171–77 Solutions to linear system, accuracy · 87–88 _ · 9–39 algorithm · 9–11, 63 in Complex mode · 63 Solving a system of equations · 16–18, 83, 99, 102–107 Solv
PAGE 186
Index V Variables, transforming · 47–48 W Weighted least-squares · 93–94, 97–98, 120 Weighted normal equations · 94 Y Yield · 35 Z Zero of polynomial · 13 186