| Linus Torvalds | 1da177e | 2005-04-16 15:20:36 -0700 | [diff] [blame] | 1 | | | 
|  | 2 | |	sasin.sa 3.3 12/19/90 | 
|  | 3 | | | 
|  | 4 | |	Description: The entry point sAsin computes the inverse sine of | 
|  | 5 | |		an input argument; sAsind does the same except for denormalized | 
|  | 6 | |		input. | 
|  | 7 | | | 
|  | 8 | |	Input: Double-extended number X in location pointed to | 
|  | 9 | |		by address register a0. | 
|  | 10 | | | 
|  | 11 | |	Output: The value arcsin(X) returned in floating-point register Fp0. | 
|  | 12 | | | 
|  | 13 | |	Accuracy and Monotonicity: The returned result is within 3 ulps in | 
|  | 14 | |		64 significant bit, i.e. within 0.5001 ulp to 53 bits if the | 
|  | 15 | |		result is subsequently rounded to double precision. The | 
|  | 16 | |		result is provably monotonic in double precision. | 
|  | 17 | | | 
|  | 18 | |	Speed: The program sASIN takes approximately 310 cycles. | 
|  | 19 | | | 
|  | 20 | |	Algorithm: | 
|  | 21 | | | 
|  | 22 | |	ASIN | 
|  | 23 | |	1. If |X| >= 1, go to 3. | 
|  | 24 | | | 
|  | 25 | |	2. (|X| < 1) Calculate asin(X) by | 
|  | 26 | |		z := sqrt( [1-X][1+X] ) | 
|  | 27 | |		asin(X) = atan( x / z ). | 
|  | 28 | |		Exit. | 
|  | 29 | | | 
|  | 30 | |	3. If |X| > 1, go to 5. | 
|  | 31 | | | 
|  | 32 | |	4. (|X| = 1) sgn := sign(X), return asin(X) := sgn * Pi/2. Exit. | 
|  | 33 | | | 
|  | 34 | |	5. (|X| > 1) Generate an invalid operation by 0 * infinity. | 
|  | 35 | |		Exit. | 
|  | 36 | | | 
|  | 37 |  | 
|  | 38 | |		Copyright (C) Motorola, Inc. 1990 | 
|  | 39 | |			All Rights Reserved | 
|  | 40 | | | 
| Matt Waddel | e00d82d | 2006-02-11 17:55:48 -0800 | [diff] [blame] | 41 | |       For details on the license for this file, please see the | 
|  | 42 | |       file, README, in this same directory. | 
| Linus Torvalds | 1da177e | 2005-04-16 15:20:36 -0700 | [diff] [blame] | 43 |  | 
|  | 44 | |SASIN	idnt	2,1 | Motorola 040 Floating Point Software Package | 
|  | 45 |  | 
|  | 46 | |section	8 | 
|  | 47 |  | 
|  | 48 | PIBY2:	.long 0x3FFF0000,0xC90FDAA2,0x2168C235,0x00000000 | 
|  | 49 |  | 
|  | 50 | |xref	t_operr | 
|  | 51 | |xref	t_frcinx | 
|  | 52 | |xref	t_extdnrm | 
|  | 53 | |xref	satan | 
|  | 54 |  | 
|  | 55 | .global	sasind | 
|  | 56 | sasind: | 
|  | 57 | |--ASIN(X) = X FOR DENORMALIZED X | 
|  | 58 |  | 
|  | 59 | bra		t_extdnrm | 
|  | 60 |  | 
|  | 61 | .global	sasin | 
|  | 62 | sasin: | 
|  | 63 | fmovex		(%a0),%fp0	| ...LOAD INPUT | 
|  | 64 |  | 
|  | 65 | movel		(%a0),%d0 | 
|  | 66 | movew		4(%a0),%d0 | 
|  | 67 | andil		#0x7FFFFFFF,%d0 | 
|  | 68 | cmpil		#0x3FFF8000,%d0 | 
|  | 69 | bges		asinbig | 
|  | 70 |  | 
|  | 71 | |--THIS IS THE USUAL CASE, |X| < 1 | 
|  | 72 | |--ASIN(X) = ATAN( X / SQRT( (1-X)(1+X) ) ) | 
|  | 73 |  | 
|  | 74 | fmoves		#0x3F800000,%fp1 | 
|  | 75 | fsubx		%fp0,%fp1		| ...1-X | 
|  | 76 | fmovemx	%fp2-%fp2,-(%a7) | 
|  | 77 | fmoves		#0x3F800000,%fp2 | 
|  | 78 | faddx		%fp0,%fp2		| ...1+X | 
|  | 79 | fmulx		%fp2,%fp1		| ...(1+X)(1-X) | 
|  | 80 | fmovemx	(%a7)+,%fp2-%fp2 | 
|  | 81 | fsqrtx		%fp1		| ...SQRT([1-X][1+X]) | 
|  | 82 | fdivx		%fp1,%fp0		| ...X/SQRT([1-X][1+X]) | 
|  | 83 | fmovemx	%fp0-%fp0,(%a0) | 
|  | 84 | bsr		satan | 
|  | 85 | bra		t_frcinx | 
|  | 86 |  | 
|  | 87 | asinbig: | 
|  | 88 | fabsx		%fp0	 | ...|X| | 
|  | 89 | fcmps		#0x3F800000,%fp0 | 
|  | 90 | fbgt		t_operr		|cause an operr exception | 
|  | 91 |  | 
|  | 92 | |--|X| = 1, ASIN(X) = +- PI/2. | 
|  | 93 |  | 
|  | 94 | fmovex		PIBY2,%fp0 | 
|  | 95 | movel		(%a0),%d0 | 
|  | 96 | andil		#0x80000000,%d0	| ...SIGN BIT OF X | 
|  | 97 | oril		#0x3F800000,%d0	| ...+-1 IN SGL FORMAT | 
|  | 98 | movel		%d0,-(%sp)	| ...push SIGN(X) IN SGL-FMT | 
|  | 99 | fmovel		%d1,%FPCR | 
|  | 100 | fmuls		(%sp)+,%fp0 | 
|  | 101 | bra		t_frcinx | 
|  | 102 |  | 
|  | 103 | |end |