Patrick Welche <prlw1@cam.ac.uk>
[netbsd-mini2440.git] / lib / libm / arch / vax / n_support.S
blob02f0e8a407e4ead75cc38dada32c92f59c94fe9d
1 /*      $NetBSD: n_support.S,v 1.5 2002/06/23 21:55:12 matt Exp $       */
2 /*
3  * Copyright (c) 1985, 1993
4  *      The Regents of the University of California.  All rights reserved.
5  *
6  * Redistribution and use in source and binary forms, with or without
7  * modification, are permitted provided that the following conditions
8  * are met:
9  * 1. Redistributions of source code must retain the above copyright
10  *    notice, this list of conditions and the following disclaimer.
11  * 2. Redistributions in binary form must reproduce the above copyright
12  *    notice, this list of conditions and the following disclaimer in the
13  *    documentation and/or other materials provided with the distribution.
14  * 3. Neither the name of the University nor the names of its contributors
15  *    may be used to endorse or promote products derived from this software
16  *    without specific prior written permission.
17  *
18  * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
19  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
21  * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
22  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
23  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
24  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
25  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
26  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
27  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
28  * SUCH DAMAGE.
29  *
30  *      @(#)support.s   8.1 (Berkeley) 6/4/93
31  */
32 #include <machine/asm.h>
34         .text
35 _sccsid:
36         .asciz "@(#)support.s\t1.3 (Berkeley) 8/21/85; 8.1 (ucb.elefunt) 6/4/93"
39  * copysign(x,y),
40  * logb(x),
41  * scalb(x,N),
42  * finite(x),
43  * drem(x,y),
44  * Coded in vax assembly language by K.C. Ng,  3/14/85.
45  * Revised by K.C. Ng on 4/9/85.
46  */
49  * double copysign(double x,double y)
50  */
52 ENTRY(copysign, 0)
53         movq    4(%ap),%r0              # load x into %r0
54         bicw3   $0x807f,%r0,%r2         # mask off the exponent of x
55         beql    Lz                      # if zero or reserved op then return x
56         bicw3   $0x7fff,12(%ap),%r2     # copy the sign bit of y into %r2
57         bicw2   $0x8000,%r0             # replace x by |x|
58         bisw2   %r2,%r0                 # copy the sign bit of y to x
59 Lz:     ret
62  * double logb(double x);
63  */
64 ENTRY(logb, 0)
65         bicl3   $0xffff807f,4(%ap),%r0  # mask off the exponent of x
66         beql    Ln
67         ashl    $-7,%r0,%r0             # get the bias exponent
68         subl2   $129,%r0                        # get the unbias exponent
69         cvtld   %r0,%r0                 # return the answer in double
70         ret
71 Ln:     movq    4(%ap),%r0              # %r0:1 = x (zero or reserved op)
72         bneq    1f                      # simply return if reserved op
73         movq    $0x0000fe00ffffcfff,%r0  # -2147483647.0
74 1:      ret
77  * long finite(double x);
78  */
79 #ifndef __GFLOAT__
80         .globl finitef
81 finitef = finite
82 #endif
83 ENTRY(finite, 0)
84         bicw3   $0x7f,4(%ap),%r0        # mask off the mantissa
85         cmpw    %r0,$0x8000             # to see if x is the reserved op
86         beql    1f                      # if so, return FALSE (0)
87         movl    $1,%r0                  # else return TRUE (1)
88         ret
89 1:      clrl    %r0
90         ret
92 /* int isnan(double x);
93  */
94 #if 0
95 ENTRY(isnan, 0)
96         clrl    %r0
97         ret
98 #endif
100 /* int isnanf(float x);
101  */
102 ENTRY(isnanf, 0)
103         clrl    %r0
104         ret
107  * double scalb(x,N)
108  * double x; double N;
109  */
110         .set    ERANGE,34
112 ENTRY(scalb, 0)
113         movq    4(%ap),%r0
114         bicl3   $0xffff807f,%r0,%r3
115         beql    ret1                    # 0 or reserved operand
116         movq    12(%ap),%r4
117         cvtdl   %r4, %r2
118         cmpl    %r2,$0x12c
119         bgeq    ovfl
120         cmpl    %r2,$-0x12c
121         bleq    unfl
122         ashl    $7,%r2,%r2
123         addl2   %r2,%r3
124         bleq    unfl
125         cmpl    %r3,$0x8000
126         bgeq    ovfl
127         addl2   %r2,%r0
128         ret
129 ovfl:   pushl   $ERANGE
130         calls   $1,_C_LABEL(infnan)     # if it returns
131         bicw3   $0x7fff,4(%ap),%r2      # get the sign of input arg
132         bisw2   %r2,%r0                 # re-attach the sign to %r0/1
133         ret
134 unfl:   movq    $0,%r0
135 ret1:   ret
138  * DREM(X,Y)
139  * RETURN X REM Y =X-N*Y, N=[X/Y] ROUNDED (ROUNDED TO EVEN IN THE HALF WAY CASE)
140  * DOUBLE PRECISION (VAX D format 56 bits)
141  * CODED IN VAX ASSEMBLY LANGUAGE BY K.C. NG, 4/8/85.
142  */
143         .set    EDOM,33
145 ENTRY(drem, 0x0fc0)
146         subl2   $12,%sp
147         movq    4(%ap),%r0              #%r0=x
148         movq    12(%ap),%r2             #%r2=y
149         jeql    Rop                     #if y=0 then generate reserved op fault
150         bicw3   $0x007f,%r0,%r4         #check if x is Rop
151         cmpw    %r4,$0x8000
152         jeql    Ret                     #if x is Rop then return Rop
153         bicl3   $0x007f,%r2,%r4         #check if y is Rop
154         cmpw    %r4,$0x8000
155         jeql    Ret                     #if y is Rop then return Rop
156         bicw2   $0x8000,%r2             #y  := |y|
157         movw    $0,-4(%fp)              #-4(%fp) = nx := 0
158         cmpw    %r2,$0x1c80             #yexp ? 57
159         bgtr    C1                      #if yexp > 57 goto C1
160         addw2   $0x1c80,%r2             #scale up y by 2**57
161         movw    $0x1c80,-4(%fp)         #nx := 57 (exponent field)
163         movw    -4(%fp),-8(%fp)         #-8(%fp) = nf := nx
164         bicw3   $0x7fff,%r0,-12(%fp)    #-12(%fp) = sign of x
165         bicw2   $0x8000,%r0             #x  := |x|
166         movq    %r2,%r10                        #y1 := y
167         bicl2   $0xffff07ff,%r11                #clear the last 27 bits of y1
168 loop:
169         cmpd    %r0,%r2                 #x ? y
170         bleq    E1                      #if x <= y goto E1
171  /* begin argument reduction */
172         movq    %r2,%r4                 #t =y
173         movq    %r10,%r6                        #t1=y1
174         bicw3   $0x807f,%r0,%r8         #xexp= exponent of x
175         bicw3   $0x807f,%r2,%r9         #yexp= exponent fo y
176         subw2   %r9,%r8                 #xexp-yexp
177         subw2   $0x0c80,%r8             #k=xexp-yexp-25(exponent bit field)
178         blss    C2                      #if k<0 goto C2
179         addw2   %r8,%r4                 #t +=k
180         addw2   %r8,%r6                 #t1+=k, scale up t and t1
182         divd3   %r4,%r0,%r8             #x/t
183         cvtdl   %r8,%r8                 #n=[x/t] truncated
184         cvtld   %r8,%r8                 #float(n)
185         subd2   %r6,%r4                 #t:=t-t1
186         muld2   %r8,%r4                 #n*(t-t1)
187         muld2   %r8,%r6                 #n*t1
188         subd2   %r6,%r0                 #x-n*t1
189         subd2   %r4,%r0                 #(x-n*t1)-n*(t-t1)
190         jbr     loop
192         movw    -4(%fp),%r6             #%r6=nx
193         beql    C3                      #if nx=0 goto C3
194         addw2   %r6,%r0                 #x:=x*2**57 scale up x by nx
195         movw    $0,-4(%fp)              #clear nx
196         jbr     loop
198         movq    %r2,%r4                 #%r4 = y
199         subw2   $0x80,%r4               #%r4 = y/2
200         cmpd    %r0,%r4                 #x:y/2
201         blss    E2                      #if x < y/2 goto E2
202         bgtr    C4                      #if x > y/2 goto C4
203         cvtdl   %r8,%r8                 #ifix(float(n))
204         blbc    %r8,E2                  #if the last bit is zero, goto E2
206         subd2   %r2,%r0                 #x-y
208         xorw2   -12(%fp),%r0            #x^sign (exclusive or)
209         movw    -8(%fp),%r6             #%r6=nf
210         bicw3   $0x807f,%r0,%r8         #%r8=exponent of x
211         bicw2   $0x7f80,%r0             #clear the exponent of x
212         subw2   %r6,%r8                 #%r8=xexp-nf
213         bgtr    C5                      #if xexp-nf is positive goto C5
214         movw    $0,%r8                  #clear %r8
215         movq    $0,%r0                  #x underflow to zero
217         bisw2   %r8,%r0                 /* put %r8 into x's exponent field */
218         ret
219 Rop:                                    #Reserved operand
220         pushl   $EDOM
221         calls   $1,_C_LABEL(infnan)     #generate reserved op fault
222         ret
223 Ret:
224         movq    $0x8000,%r0             #propagate reserved op
225         ret