1# 2# Copyright (c) 1985 Regents of the University of California. 3# 4# Use and reproduction of this software are granted in accordance with 5# the terms and conditions specified in the Berkeley Software License 6# Agreement (in particular, this entails acknowledgement of the programs' 7# source, and inclusion of this notice) with the additional understanding 8# that all recipients should regard themselves as participants in an 9# ongoing research project and hence should feel obligated to report 10# their experiences (good or bad) with these elementary function codes, 11# using "sendbug 4bsd-bugs@BERKELEY", to the authors. 12# 13 14# @(#)argred.s 1.1 (Berkeley) 08/21/85 15 16# libm$argred implements Bob Corbett's argument reduction and 17# libm$sincos implements Peter Tang's double precision sin/cos. 18# 19# Note: The two entry points libm$argred and libm$sincos are meant 20# to be used only by _sin, _cos and _tan. 21# 22# method: true range reduction to [-pi/4,pi/4], P. Tang & B. Corbett 23# S. McDonald, April 4, 1985 24# 25 .globl libm$argred 26 .globl libm$sincos 27 .text 28 .align 1 29 30libm$argred: 31# 32# Compare the argument with the largest possible that can 33# be reduced by table lookup. r3 := |x| will be used in table_lookup . 34# 35 movd r0,r3 36 bgeq abs1 37 mnegd r3,r3 38abs1: 39 cmpd r3,$0d+4.55530934770520019583e+01 40 blss small_arg 41 jsb trigred 42 rsb 43small_arg: 44 jsb table_lookup 45 rsb 46# 47# At this point, 48# r0 contains the quadrant number, 0, 1, 2, or 3; 49# r2/r1 contains the reduced argument as a D-format number; 50# r3 contains a F-format extension to the reduced argument; 51# r4 contains a 0 or 1 corresponding to a sin or cos entry. 52# 53libm$sincos: 54# 55# Compensate for a cosine entry by adding one to the quadrant number. 56# 57 addl2 r4,r0 58# 59# Polyd clobbers r5-r0 ; save X in r7/r6 . 60# This can be avoided by rewriting trigred . 61# 62 movd r1,r6 63# 64# Likewise, save alpha in r8 . 65# This can be avoided by rewriting trigred . 66# 67 movf r3,r8 68# 69# Odd or even quadrant? cosine if odd, sine otherwise. 70# Save floor(quadrant/2) in r9 ; it determines the final sign. 71# 72 rotl $-1,r0,r9 73 blss cosine 74sine: 75 muld2 r1,r1 # Xsq = X * X 76 polyd r1,$7,sin_coef # Q = P(Xsq) , of deg 7 77 mulf3 $0f3.0,r8,r4 # beta = 3 * alpha 78 mulf2 r0,r4 # beta = Q * beta 79 addf2 r8,r4 # beta = alpha + beta 80 muld2 r6,r0 # S(X) = X * Q 81# cvtfd r4,r4 ... r5 = 0 after a polyd. 82 addd2 r4,r0 # S(X) = beta + S(X) 83 addd2 r6,r0 # S(X) = X + S(X) 84 brb done 85cosine: 86 muld2 r6,r6 # Xsq = X * X 87 beql zero_arg 88 mulf2 r1,r8 # beta = X * alpha 89 polyd r6,$7,cos_coef # Q = P'(Xsq) , of deg 7 90 subd3 r0,r8,r0 # beta = beta - Q 91 subw2 $0x80,r6 # Xsq = Xsq / 2 92 addd2 r0,r6 # Xsq = Xsq + beta 93zero_arg: 94 subd3 r6,$0d1.0,r0 # C(X) = 1 - Xsq 95done: 96 blbc r9,even 97 mnegd r0,r0 98even: 99 rsb 100 101.data 102.align 2 103 104sin_coef: 105 .double 0d-7.53080332264191085773e-13 # s7 = 2^-29 -1.a7f2504ffc49f8.. 106 .double 0d+1.60573519267703489121e-10 # s6 = 2^-21 1.611adaede473c8.. 107 .double 0d-2.50520965150706067211e-08 # s5 = 2^-1a -1.ae644921ed8382.. 108 .double 0d+2.75573191800593885716e-06 # s4 = 2^-13 1.71de3a4b884278.. 109 .double 0d-1.98412698411850507950e-04 # s3 = 2^-0d -1.a01a01a0125e7d.. 110 .double 0d+8.33333333333325688985e-03 # s2 = 2^-07 1.11111111110e50 111 .double 0d-1.66666666666666664354e-01 # s1 = 2^-03 -1.55555555555554 112 .double 0d+0.00000000000000000000e+00 # s0 = 0 113 114cos_coef: 115 .double 0d-1.13006966202629430300e-11 # s7 = 2^-25 -1.8D9BA04D1374BE.. 116 .double 0d+2.08746646574796004700e-09 # s6 = 2^-1D 1.1EE632650350BA.. 117 .double 0d-2.75573073031284417300e-07 # s5 = 2^-16 -1.27E4F31411719E.. 118 .double 0d+2.48015872682668025200e-05 # s4 = 2^-10 1.A01A0196B902E8.. 119 .double 0d-1.38888888888464709200e-03 # s3 = 2^-0A -1.6C16C16C11FACE.. 120 .double 0d+4.16666666666664761400e-02 # s2 = 2^-05 1.5555555555539E 121 .double 0d+0.00000000000000000000e+00 # s1 = 0 122 .double 0d+0.00000000000000000000e+00 # s0 = 0 123 124# 125# Multiples of pi/2 expressed as the sum of three doubles, 126# 127# trailing: n * pi/2 , n = 0, 1, 2, ..., 29 128# trailing[n] , 129# 130# middle: n * pi/2 , n = 0, 1, 2, ..., 29 131# middle[n] , 132# 133# leading: n * pi/2 , n = 0, 1, 2, ..., 29 134# leading[n] , 135# 136# where 137# leading[n] := (n * pi/2) rounded, 138# middle[n] := (n * pi/2 - leading[n]) rounded, 139# trailing[n] := (( n * pi/2 - leading[n]) - middle[n]) rounded . 140 141trailing: 142 .double 0d+0.00000000000000000000e+00 # 0 * pi/2 trailing 143 .double 0d+4.33590506506189049611e-35 # 1 * pi/2 trailing 144 .double 0d+8.67181013012378099223e-35 # 2 * pi/2 trailing 145 .double 0d+1.30077151951856714215e-34 # 3 * pi/2 trailing 146 .double 0d+1.73436202602475619845e-34 # 4 * pi/2 trailing 147 .double 0d-1.68390735624352669192e-34 # 5 * pi/2 trailing 148 .double 0d+2.60154303903713428430e-34 # 6 * pi/2 trailing 149 .double 0d-8.16726343231148352150e-35 # 7 * pi/2 trailing 150 .double 0d+3.46872405204951239689e-34 # 8 * pi/2 trailing 151 .double 0d+3.90231455855570147991e-34 # 9 * pi/2 trailing 152 .double 0d-3.36781471248705338384e-34 # 10 * pi/2 trailing 153 .double 0d-1.06379439835298071785e-33 # 11 * pi/2 trailing 154 .double 0d+5.20308607807426856861e-34 # 12 * pi/2 trailing 155 .double 0d+5.63667658458045770509e-34 # 13 * pi/2 trailing 156 .double 0d-1.63345268646229670430e-34 # 14 * pi/2 trailing 157 .double 0d-1.19986217995610764801e-34 # 15 * pi/2 trailing 158 .double 0d+6.93744810409902479378e-34 # 16 * pi/2 trailing 159 .double 0d-8.03640094449267300110e-34 # 17 * pi/2 trailing 160 .double 0d+7.80462911711140295982e-34 # 18 * pi/2 trailing 161 .double 0d-7.16921993148029483506e-34 # 19 * pi/2 trailing 162 .double 0d-6.73562942497410676769e-34 # 20 * pi/2 trailing 163 .double 0d-6.30203891846791677593e-34 # 21 * pi/2 trailing 164 .double 0d-2.12758879670596143570e-33 # 22 * pi/2 trailing 165 .double 0d+2.53800212047402350390e-33 # 23 * pi/2 trailing 166 .double 0d+1.04061721561485371372e-33 # 24 * pi/2 trailing 167 .double 0d+6.11729905311472319056e-32 # 25 * pi/2 trailing 168 .double 0d+1.12733531691609154102e-33 # 26 * pi/2 trailing 169 .double 0d-3.70049587943078297272e-34 # 27 * pi/2 trailing 170 .double 0d-3.26690537292459340860e-34 # 28 * pi/2 trailing 171 .double 0d-1.14812616507957271361e-34 # 29 * pi/2 trailing 172 173middle: 174 .double 0d+0.00000000000000000000e+00 # 0 * pi/2 middle 175 .double 0d+5.72118872610983179676e-18 # 1 * pi/2 middle 176 .double 0d+1.14423774522196635935e-17 # 2 * pi/2 middle 177 .double 0d-3.83475850529283316309e-17 # 3 * pi/2 middle 178 .double 0d+2.28847549044393271871e-17 # 4 * pi/2 middle 179 .double 0d-2.69052076007086676522e-17 # 5 * pi/2 middle 180 .double 0d-7.66951701058566632618e-17 # 6 * pi/2 middle 181 .double 0d-1.54628301484890040587e-17 # 7 * pi/2 middle 182 .double 0d+4.57695098088786543741e-17 # 8 * pi/2 middle 183 .double 0d+1.07001849766246313192e-16 # 9 * pi/2 middle 184 .double 0d-5.38104152014173353044e-17 # 10 * pi/2 middle 185 .double 0d-2.14622680169080983801e-16 # 11 * pi/2 middle 186 .double 0d-1.53390340211713326524e-16 # 12 * pi/2 middle 187 .double 0d-9.21580002543456677056e-17 # 13 * pi/2 middle 188 .double 0d-3.09256602969780081173e-17 # 14 * pi/2 middle 189 .double 0d+3.03066796603896507006e-17 # 15 * pi/2 middle 190 .double 0d+9.15390196177573087482e-17 # 16 * pi/2 middle 191 .double 0d+1.52771359575124969107e-16 # 17 * pi/2 middle 192 .double 0d+2.14003699532492626384e-16 # 18 * pi/2 middle 193 .double 0d-1.68853170360202329427e-16 # 19 * pi/2 middle 194 .double 0d-1.07620830402834670609e-16 # 20 * pi/2 middle 195 .double 0d+3.97700719404595604379e-16 # 21 * pi/2 middle 196 .double 0d-4.29245360338161967602e-16 # 22 * pi/2 middle 197 .double 0d-3.68013020380794313406e-16 # 23 * pi/2 middle 198 .double 0d-3.06780680423426653047e-16 # 24 * pi/2 middle 199 .double 0d-2.45548340466059054318e-16 # 25 * pi/2 middle 200 .double 0d-1.84316000508691335411e-16 # 26 * pi/2 middle 201 .double 0d-1.23083660551323675053e-16 # 27 * pi/2 middle 202 .double 0d-6.18513205939560162346e-17 # 28 * pi/2 middle 203 .double 0d-6.18980636588357585202e-19 # 29 * pi/2 middle 204 205leading: 206 .double 0d+0.00000000000000000000e+00 # 0 * pi/2 leading 207 .double 0d+1.57079632679489661351e+00 # 1 * pi/2 leading 208 .double 0d+3.14159265358979322702e+00 # 2 * pi/2 leading 209 .double 0d+4.71238898038468989604e+00 # 3 * pi/2 leading 210 .double 0d+6.28318530717958645404e+00 # 4 * pi/2 leading 211 .double 0d+7.85398163397448312306e+00 # 5 * pi/2 leading 212 .double 0d+9.42477796076937979208e+00 # 6 * pi/2 leading 213 .double 0d+1.09955742875642763501e+01 # 7 * pi/2 leading 214 .double 0d+1.25663706143591729081e+01 # 8 * pi/2 leading 215 .double 0d+1.41371669411540694661e+01 # 9 * pi/2 leading 216 .double 0d+1.57079632679489662461e+01 # 10 * pi/2 leading 217 .double 0d+1.72787595947438630262e+01 # 11 * pi/2 leading 218 .double 0d+1.88495559215387595842e+01 # 12 * pi/2 leading 219 .double 0d+2.04203522483336561422e+01 # 13 * pi/2 leading 220 .double 0d+2.19911485751285527002e+01 # 14 * pi/2 leading 221 .double 0d+2.35619449019234492582e+01 # 15 * pi/2 leading 222 .double 0d+2.51327412287183458162e+01 # 16 * pi/2 leading 223 .double 0d+2.67035375555132423742e+01 # 17 * pi/2 leading 224 .double 0d+2.82743338823081389322e+01 # 18 * pi/2 leading 225 .double 0d+2.98451302091030359342e+01 # 19 * pi/2 leading 226 .double 0d+3.14159265358979324922e+01 # 20 * pi/2 leading 227 .double 0d+3.29867228626928286062e+01 # 21 * pi/2 leading 228 .double 0d+3.45575191894877260523e+01 # 22 * pi/2 leading 229 .double 0d+3.61283155162826226103e+01 # 23 * pi/2 leading 230 .double 0d+3.76991118430775191683e+01 # 24 * pi/2 leading 231 .double 0d+3.92699081698724157263e+01 # 25 * pi/2 leading 232 .double 0d+4.08407044966673122843e+01 # 26 * pi/2 leading 233 .double 0d+4.24115008234622088423e+01 # 27 * pi/2 leading 234 .double 0d+4.39822971502571054003e+01 # 28 * pi/2 leading 235 .double 0d+4.55530934770520019583e+01 # 29 * pi/2 leading 236 237twoOverPi: 238 .double 0d+6.36619772367581343076e-01 239 .text 240 .align 1 241 242table_lookup: 243 muld3 r3,twoOverPi,r0 244 cvtrdl r0,r0 # n = nearest int to ((2/pi)*|x|) rnded 245 mull3 $8,r0,r5 246 subd2 leading(r5),r3 # p = (|x| - leading n*pi/2) exactly 247 subd3 middle(r5),r3,r1 # q = (p - middle n*pi/2) rounded 248 subd2 r1,r3 # r = (p - q) 249 subd2 middle(r5),r3 # r = r - middle n*pi/2 250 subd2 trailing(r5),r3 # r = r - trailing n*pi/2 rounded 251# 252# If the original argument was negative, 253# negate the reduce argument and 254# adjust the octant/quadrant number. 255# 256 tstw 4(ap) 257 bgeq abs2 258 mnegf r1,r1 259 mnegf r3,r3 260# subb3 r0,$8,r0 ...used for pi/4 reduction -S.McD 261 subb3 r0,$4,r0 262abs2: 263# 264# Clear all unneeded octant/quadrant bits. 265# 266# bicb2 $0xf8,r0 ...used for pi/4 reduction -S.McD 267 bicb2 $0xfc,r0 268 rsb 269# 270# p.0 271 .text 272 .align 2 273# 274# Only 256 (actually 225) bits of 2/pi are needed for VAX double 275# precision; this was determined by enumerating all the nearest 276# machine integer multiples of pi/2 using continued fractions. 277# (8a8d3673775b7ff7 required the most bits.) -S.McD 278# 279 .long 0 280 .long 0 281 .long 0xaef1586d 282 .long 0x9458eaf7 283 .long 0x10e4107f 284 .long 0xd8a5664f 285 .long 0x4d377036 286 .long 0x09d5f47d 287 .long 0x91054a7f 288 .long 0xbe60db93 289bits2opi: 290 .long 0x00000028 291 .long 0 292# 293# Note: wherever you see the word `octant', read `quadrant'. 294# Currently this code is set up for pi/2 argument reduction. 295# By uncommenting/commenting the appropriate lines, it will 296# also serve as a pi/4 argument reduction code. 297# 298 299# p.1 300# Trigred preforms argument reduction 301# for the trigonometric functions. It 302# takes one input argument, a D-format 303# number in r1/r0 . The magnitude of 304# the input argument must be greater 305# than or equal to 1/2 . Trigred produces 306# three results: the number of the octant 307# occupied by the argument, the reduced 308# argument, and an extension of the 309# reduced argument. The octant number is 310# returned in r0 . The reduced argument 311# is returned as a D-format number in 312# r2/r1 . An 8 bit extension of the 313# reduced argument is returned as an 314# F-format number in r3. 315# p.2 316trigred: 317# 318# Save the sign of the input argument. 319# 320 movw r0,-(sp) 321# 322# Extract the exponent field. 323# 324 extzv $7,$7,r0,r2 325# 326# Convert the fraction part of the input 327# argument into a quadword integer. 328# 329 bicw2 $0xff80,r0 330 bisb2 $0x80,r0 # -S.McD 331 rotl $16,r0,r0 332 rotl $16,r1,r1 333# 334# If r1 is negative, add 1 to r0 . This 335# adjustment is made so that the two's 336# complement multiplications done later 337# will produce unsigned results. 338# 339 bgeq posmid 340 incl r0 341posmid: 342# p.3 343# 344# Set r3 to the address of the first quadword 345# used to obtain the needed portion of 2/pi . 346# The address is longword aligned to ensure 347# efficient access. 348# 349 ashl $-3,r2,r3 350 bicb2 $3,r3 351 subl3 r3,$bits2opi,r3 352# 353# Set r2 to the size of the shift needed to 354# obtain the correct portion of 2/pi . 355# 356 bicb2 $0xe0,r2 357# p.4 358# 359# Move the needed 128 bits of 2/pi into 360# r11 - r8 . Adjust the numbers to allow 361# for unsigned multiplication. 362# 363 ashq r2,(r3),r10 364 365 subl2 $4,r3 366 ashq r2,(r3),r9 367 bgeq signoff1 368 incl r11 369signoff1: 370 subl2 $4,r3 371 ashq r2,(r3),r8 372 bgeq signoff2 373 incl r10 374signoff2: 375 subl2 $4,r3 376 ashq r2,(r3),r7 377 bgeq signoff3 378 incl r9 379signoff3: 380# p.5 381# 382# Multiply the contents of r0/r1 by the 383# slice of 2/pi in r11 - r8 . 384# 385 emul r0,r8,$0,r4 386 emul r0,r9,r5,r5 387 emul r0,r10,r6,r6 388 389 emul r1,r8,$0,r7 390 emul r1,r9,r8,r8 391 emul r1,r10,r9,r9 392 emul r1,r11,r10,r10 393 394 addl2 r4,r8 395 adwc r5,r9 396 adwc r6,r10 397# p.6 398# 399# If there are more than five leading zeros 400# after the first two quotient bits or if there 401# are more than five leading ones after the first 402# two quotient bits, generate more fraction bits. 403# Otherwise, branch to code to produce the result. 404# 405 bicl3 $0xc1ffffff,r10,r4 406 beql more1 407 cmpl $0x3e000000,r4 408 bneq result 409more1: 410# p.7 411# 412# generate another 32 result bits. 413# 414 subl2 $4,r3 415 ashq r2,(r3),r5 416 bgeq signoff4 417 418 emul r1,r6,$0,r4 419 addl2 r1,r5 420 emul r0,r6,r5,r5 421 addl2 r0,r6 422 brb addbits1 423 424signoff4: 425 emul r1,r6,$0,r4 426 emul r0,r6,r5,r5 427 428addbits1: 429 addl2 r5,r7 430 adwc r6,r8 431 adwc $0,r9 432 adwc $0,r10 433# p.8 434# 435# Check for massive cancellation. 436# 437 bicl3 $0xc0000000,r10,r6 438# bneq more2 -S.McD Test was backwards 439 beql more2 440 cmpl $0x3fffffff,r6 441 bneq result 442more2: 443# p.9 444# 445# If massive cancellation has occurred, 446# generate another 24 result bits. 447# Testing has shown there will always be 448# enough bits after this point. 449# 450 subl2 $4,r3 451 ashq r2,(r3),r5 452 bgeq signoff5 453 454 emul r0,r6,r4,r5 455 addl2 r0,r6 456 brb addbits2 457 458signoff5: 459 emul r0,r6,r4,r5 460 461addbits2: 462 addl2 r6,r7 463 adwc $0,r8 464 adwc $0,r9 465 adwc $0,r10 466# p.10 467# 468# The following code produces the reduced 469# argument from the product bits contained 470# in r10 - r7 . 471# 472result: 473# 474# Extract the octant number from r10 . 475# 476# extzv $29,$3,r10,r0 ...used for pi/4 reduction -S.McD 477 extzv $30,$2,r10,r0 478# 479# Clear the octant bits in r10 . 480# 481# bicl2 $0xe0000000,r10 ...used for pi/4 reduction -S.McD 482 bicl2 $0xc0000000,r10 483# 484# Zero the sign flag. 485# 486 clrl r5 487# p.11 488# 489# Check to see if the fraction is greater than 490# or equal to one-half. If it is, add one 491# to the octant number, set the sign flag 492# on, and replace the fraction with 1 minus 493# the fraction. 494# 495# bitl $0x10000000,r10 ...used for pi/4 reduction -S.McD 496 bitl $0x20000000,r10 497 beql small 498 incl r0 499 incl r5 500# subl3 r10,$0x1fffffff,r10 ...used for pi/4 reduction -S.McD 501 subl3 r10,$0x3fffffff,r10 502 mcoml r9,r9 503 mcoml r8,r8 504 mcoml r7,r7 505small: 506# p.12 507# 508## Test whether the first 29 bits of the ...used for pi/4 reduction -S.McD 509# Test whether the first 30 bits of the 510# fraction are zero. 511# 512 tstl r10 513 beql tiny 514# 515# Find the position of the first one bit in r10 . 516# 517 cvtld r10,r1 518 extzv $7,$7,r1,r1 519# 520# Compute the size of the shift needed. 521# 522 subl3 r1,$32,r6 523# 524# Shift up the high order 64 bits of the 525# product. 526# 527 ashq r6,r9,r10 528 ashq r6,r8,r9 529 brb mult 530# p.13 531# 532# Test to see if the sign bit of r9 is on. 533# 534tiny: 535 tstl r9 536 bgeq tinier 537# 538# If it is, shift the product bits up 32 bits. 539# 540 movl $32,r6 541 movq r8,r10 542 tstl r10 543 brb mult 544# p.14 545# 546# Test whether r9 is zero. It is probably 547# impossible for both r10 and r9 to be 548# zero, but until proven to be so, the test 549# must be made. 550# 551tinier: 552 beql zero 553# 554# Find the position of the first one bit in r9 . 555# 556 cvtld r9,r1 557 extzv $7,$7,r1,r1 558# 559# Compute the size of the shift needed. 560# 561 subl3 r1,$32,r1 562 addl3 $32,r1,r6 563# 564# Shift up the high order 64 bits of the 565# product. 566# 567 ashq r1,r8,r10 568 ashq r1,r7,r9 569 brb mult 570# p.15 571# 572# The following code sets the reduced 573# argument to zero. 574# 575zero: 576 clrl r1 577 clrl r2 578 clrl r3 579 brw return 580# p.16 581# 582# At this point, r0 contains the octant number, 583# r6 indicates the number of bits the fraction 584# has been shifted, r5 indicates the sign of 585# the fraction, r11/r10 contain the high order 586# 64 bits of the fraction, and the condition 587# codes indicate where the sign bit of r10 588# is on. The following code multiplies the 589# fraction by pi/2 . 590# 591mult: 592# 593# Save r11/r10 in r4/r1 . -S.McD 594 movl r11,r4 595 movl r10,r1 596# 597# If the sign bit of r10 is on, add 1 to r11 . 598# 599 bgeq signoff6 600 incl r11 601signoff6: 602# p.17 603# 604# Move pi/2 into r3/r2 . 605# 606 movq $0xc90fdaa22168c235,r2 607# 608# Multiply the fraction by the portion of pi/2 609# in r2 . 610# 611 emul r2,r10,$0,r7 612 emul r2,r11,r8,r7 613# 614# Multiply the fraction by the portion of pi/2 615# in r3 . 616 emul r3,r10,$0,r9 617 emul r3,r11,r10,r10 618# 619# Add the product bits together. 620# 621 addl2 r7,r9 622 adwc r8,r10 623 adwc $0,r11 624# 625# Compensate for not sign extending r8 above.-S.McD 626# 627 tstl r8 628 bgeq signoff6a 629 decl r11 630signoff6a: 631# 632# Compensate for r11/r10 being unsigned. -S.McD 633# 634 addl2 r2,r10 635 adwc r3,r11 636# 637# Compensate for r3/r2 being unsigned. -S.McD 638# 639 addl2 r1,r10 640 adwc r4,r11 641# p.18 642# 643# If the sign bit of r11 is zero, shift the 644# product bits up one bit and increment r6 . 645# 646 blss signon 647 incl r6 648 ashq $1,r10,r10 649 tstl r9 650 bgeq signoff7 651 incl r10 652signoff7: 653signon: 654# p.19 655# 656# Shift the 56 most significant product 657# bits into r9/r8 . The sign extension 658# will be handled later. 659# 660 ashq $-8,r10,r8 661# 662# Convert the low order 8 bits of r10 663# into an F-format number. 664# 665 cvtbf r10,r3 666# 667# If the result of the conversion was 668# negative, add 1 to r9/r8 . 669# 670 bgeq chop 671 incl r8 672 adwc $0,r9 673# 674# If r9 is now zero, branch to special 675# code to handle that possibility. 676# 677 beql carryout 678chop: 679# p.20 680# 681# Convert the number in r9/r8 into 682# D-format number in r2/r1 . 683# 684 rotl $16,r8,r2 685 rotl $16,r9,r1 686# 687# Set the exponent field to the appropriate 688# value. Note that the extra bits created by 689# sign extension are now eliminated. 690# 691 subw3 r6,$131,r6 692 insv r6,$7,$9,r1 693# 694# Set the exponent field of the F-format 695# number in r3 to the appropriate value. 696# 697 tstf r3 698 beql return 699# extzv $7,$8,r3,r4 -S.McD 700 extzv $7,$7,r3,r4 701 addw2 r4,r6 702# subw2 $217,r6 -S.McD 703 subw2 $64,r6 704 insv r6,$7,$8,r3 705 brb return 706# p.21 707# 708# The following code generates the appropriate 709# result for the unlikely possibility that 710# rounding the number in r9/r8 resulted in 711# a carry out. 712# 713carryout: 714 clrl r1 715 clrl r2 716 subw3 r6,$132,r6 717 insv r6,$7,$9,r1 718 tstf r3 719 beql return 720 extzv $7,$8,r3,r4 721 addw2 r4,r6 722 subw2 $218,r6 723 insv r6,$7,$8,r3 724# p.22 725# 726# The following code makes an needed 727# adjustments to the signs of the 728# results or to the octant number, and 729# then returns. 730# 731return: 732# 733# Test if the fraction was greater than or 734# equal to 1/2 . If so, negate the reduced 735# argument. 736# 737 blbc r5,signoff8 738 mnegf r1,r1 739 mnegf r3,r3 740signoff8: 741# p.23 742# 743# If the original argument was negative, 744# negate the reduce argument and 745# adjust the octant number. 746# 747 tstw (sp)+ 748 bgeq signoff9 749 mnegf r1,r1 750 mnegf r3,r3 751# subb3 r0,$8,r0 ...used for pi/4 reduction -S.McD 752 subb3 r0,$4,r0 753signoff9: 754# 755# Clear all unneeded octant bits. 756# 757# bicb2 $0xf8,r0 ...used for pi/4 reduction -S.McD 758 bicb2 $0xfc,r0 759# 760# Return. 761# 762 rsb 763