Morgan - Numerical Methods (523161), страница 36

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 36)

If everything proceeds smoothly, the routine exits when theestimates stop changing.In some circumstances, however, the routine will hang, toggling between twopossible roots. Another escape is provided for that contingency. A counter, cntr, isloaded with the maximum number of iterations. If that number is exceeded, theroutine leaves with the last best estimate, which is probably close enough. An253NUMERICAL METHODSalternative would be to use another variable to define an error band and compare itwith the difference between each new estimate and the last; exiting when thedifference is less than the error (this-estimate -last_estimate<error).fx_sqr: Algorithm1.

Establish a limiton the number of iterations possible in cntr.Check for negative or zero input.If true, exit through step 3.leave radicand in the register to be justified and make our firstestimate.2. Decrement the limit counter, cntr.If there is a carry, exit with current estimate and the carry setthrough step 3.If there is no carry, continue.Test the estimate to see that it fits within sixteen bits.If not, shift right until it does.Store the estimate.Divide the radicand by the estimate.Add the result to the estimate.Divide that by two.Compare last estimate with current estimate.If is different continue with the beginning of step 2.Otherwise, go to step 3.3.

Write the result to the output and leave.fx_sqr: Listing; *****; accepts integers; Remember that the powers follow the powers of two (the root of a double word; is a word, the root of a word is a byte, the root of a byte is a nibble, etc.).; new_estimate = (radicand/last_estimate+last_estimate)/2, last-estimate=new-estimate.fx_sqr proc uses bx cx dx di si, radicand:dword, root:word254THE ELEMENTARY FUNCTIONSlocalmovsubmovmovorjsjejmpzero_exit:orjnesigr_exit:stcsubmovjmpestimate:word, cntr:bytebyte ptr cntr, 16bx, bxax, word ptr radicanddx, word ptr radicand [2]dx, dxsign_exitzero_exitfind_rootax, axfind_root;to test radicand;not zero;no negatives or zeros;indicate error in the;operationax, axdx, axroot_exitfind_root:subjcbyte ptr cntr, 1root_exitfind_root1:orjeshrrcrjmpdx, dxfitsdx, 1ax, 1find_root1;must be zero;some kind of estimatemovsubmovdivmovword ptr estimate, axdx, dxax, word ptr radicand [2]word ptr estimatebx, ax;store first estimate of rootmovdivmovax, word ptr radicandword ptr estimatedx, bx;will exit with carry;set and an approximate;root;cannot have a root;greater than 16 bits;for a 32-bit radicand!fits:;save quotient from division;of upper word;divide lower word;concatenate quotients255NUMERICAL METHODSaddax, word ptr estimateadcshrrcrdx, 0dx, 1ax, 1orjnecmpjneclcdx, dxfind_rootax, word ptr estimatefind_rootroot_exit:movmovmovretfx_sqr endp;(radicand/estimate +;estimate)/2;to prevent any modular aliasing;is the estimate still changing?;clear the carry to indicate;successdi, word ptr rootword ptr [di], axword ptr [di][2], dxThe next approach is based on the technique taught in school for doing squareroots by hand.

This method turns out to be much simpler in binary than in decimalbecause of its modulus of 2.8It may not be faster than Newton’s Method, but it’s agood alternative for those processors without hardware division instructions.school_sqr: Algorithm1. Determine that the radicand is positive and not zero.If so, continue with step 2.If not, signal the error and exit through step 5.2.

Set bit counter for 16.Set buffer to hold radicand and allow for shifts.Clear space for root.3. Shift buffer left twice.Shift root left once.Subtract 2*root+l from root.If there is an underflow, restore the subtraction by means ofaddition and continue with step 4.Otherwise, increment the root and continue with step 4.256THE ELEMENTARY FUNCTIONS4. Decrement bit counter.If zero, exit through step 5.Otherwise, continue with step 3.5. Write root to output and leave.school_sqr: Listing; ******;school_sqr;accepts integersschool_sqrlocalsubmovmovorjsjejmpzero_exit:orjnesign_exit:proc uses bx cx dx di si, radicand:dword, root:wordestimate:qword, bits:bytebx, bxax, word ptr radicanddx, word ptr radicand [2]dx, dxsign-exitzero-exitsetup;notzero;no negatives or zerosax, axsetup;indicate error in the;operation; can't do;negatives;zero for failsubax, axstcjmproot_exitmovmovmovsubmovmovmovmovmovbyte ptrword ptrword ptrax, axword ptrword ptrbx, axcx, axdx, axsetup:bits, 16estimate, axestimate [2], dxestimate [4], axestimate [6], ax;root;intermediatefindroot:257NUMERICAL METHODSshlrclrclrclwordwordwordwordptrptrptrptrestimate, 1estimate[2], 1estimate[4], 1estimate[6], 1shlrclrclrclwordwordwordwordptrptrptrptrestimate, 1estimate[2], 1estimate[4], 1estimate[6], 1shlrclax, 1bx, 1movmovshlrcladdadccx,dx,cx,dx,cx,dx,subtract_root:subsbbjncaddadcjmpr_plus_one:addadccontinue_loop:decjneclcroot-exit:movmovmovretschool_sqr endp258ax,bx,;shift rootaxbx1110word ptr estimate[4],word ptr estimate[6],r_plus_oneword ptr estimate[4],word ptr estimate[6],continue_loop10byte ptr bitsfindrootdi, word ptr rootword ptr [di], axword ptr [di][2], bx;double-shift radicand;root*2;+lcxdx;accumulator-2*root+lcxdx;r+=lTHE ELEMENTARY FUNCTIONSFloating-Point ApproximationsAll the techniques explored so far in fixed point apply to floating-pointarithmetic as well.

Floating-point arithmetic is similar to fixed point except that itdeals with real numbers with far greater range. And because of its extensive use inscientific and engineering applications, greater emphasis is placed on its ability toapproximate the real world.This section presents some concepts that can also be used in fixed-point routines,but they’re most valuable in floating point because of its attention to accuracy.

Twoapproximations will also be described- a sine function and square root-based onmaterials from Software Manual for the Elementary Functions by William J. Cody,Jr. and William Waite, published by Prentice-Hall, Inc. This small book is full ofvaluable information for those writing numerical software. The sine/cosine approximation uses a minimax polynomial approximation, and the square root usesNewton’s Method with a much improved initial estimate.Floating-Point UtilitiesThe functions in this section use similar techniques to the fixed-point routines;that is, they use tables or arrays of coefficients and Homer’s rule for evaluatingpolynomial approximations to the functions.

The floating-point format also hassome new tools and requires some new handling.Many of the manipulations require argument reduction, which takes the floatingpoint word apart and puts it back together again in a different fashion. Some newfunctions will be presented here for doing that. One is frxp, which, when passed afloat (x) returns its exponent (n) and the float (f) constrained to a value between .5f < 1, where f* 2n = x.

Because it is the power to which the fixed-point mantissa mustbe raised to represent that number, the exponent is useful in finding the square rootof a number, as you’ll see in flsqr.frxp: Algorithm1. Point to the variable for the exponent.2. Test the number to see if it's zero.If so, return zero as both the exponent and the mantissa.259NUMERICAL METHODSIf not, continue with step 3.3. Discard the sign bit and subtract 126D to get the exponent, write itto exptr, and replace the exponent in the number with 126D.4.

Realign the float and write it to fraction.5. Return.frxp: Listing; *****;Frxp performs an operation similar to the C function frexp. Used;for floating-point math routines.;Returns the exponent-bias of a floating-point number.;does not convert to floating point first, but expects a single;precision number on the stack.;;uses di, float:qword, fraction:word, exptr:wordfrxp procpushfcldmovmovmovsubororjeshlsubmovmovshrmovmovleamovmovswrepfrxp_exit:popfretmake_it_zerosub260di, word ptr exptrax, word ptr float[4]dx, word ptr float[2]cx, cxcx, axcx, dxmake_it_zeroax, 1ah, 7ehbyte ptr [di],ahah, 7ehax, 1word ptr float[4], axdi, word ptr fractionsi, word ptr floatcx, 4ax, ax;assign pointer to exponent;get upper word of float;the sign means zero;subtract bias to place;float .5<=x<l;replace sign;write out new floatTHE ELEMENTARY FUNCTIONSrepfrxpmovmovstoswjmpendpbyte ptr [di], aldi, word ptr fractionfrxp_exitLdxp performs essentially the inverse of frxp.

Характеристики

Список файлов книги