home *** CD-ROM | disk | FTP | other *** search
/ Frozen Fish 1: Amiga / FrozenFish-Apr94.iso / bbs / gnu / libm-src.lha / src / amiga / libm / vax / support.s < prev    next >
Text File  |  1993-09-23  |  6KB  |  229 lines
  1. /*
  2.  * Copyright (c) 1985 Regents of the University of California.
  3.  * All rights reserved.
  4.  *
  5.  * Redistribution and use in source and binary forms, with or without
  6.  * modification, are permitted provided that the following conditions
  7.  * are met:
  8.  * 1. Redistributions of source code must retain the above copyright
  9.  *    notice, this list of conditions and the following disclaimer.
  10.  * 2. Redistributions in binary form must reproduce the above copyright
  11.  *    notice, this list of conditions and the following disclaimer in the
  12.  *    documentation and/or other materials provided with the distribution.
  13.  * 3. All advertising materials mentioning features or use of this software
  14.  *    must display the following acknowledgement:
  15.  *    This product includes software developed by the University of
  16.  *    California, Berkeley and its contributors.
  17.  * 4. Neither the name of the University nor the names of its contributors
  18.  *    may be used to endorse or promote products derived from this software
  19.  *    without specific prior written permission.
  20.  *
  21.  * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
  22.  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  23.  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  24.  * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
  25.  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
  26.  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
  27.  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
  28.  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
  29.  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
  30.  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
  31.  * SUCH DAMAGE.
  32.  *
  33.  *    @(#)support.s    5.5 (Berkeley) 10/9/90
  34.  */
  35.     .data
  36.     .align    2
  37. _sccsid:
  38. .asciz    "@(#)support.s    1.3 (Berkeley) 8/21/85; 5.5 (ucb.elefunt) 10/9/90"
  39.  
  40. /*
  41.  * copysign(x,y),
  42.  * logb(x),
  43.  * scalb(x,N),
  44.  * finite(x),
  45.  * drem(x,y),
  46.  * Coded in vax assembly language by K.C. Ng,  3/14/85.
  47.  * Revised by K.C. Ng on 4/9/85.
  48.  */
  49.  
  50. /*
  51.  * double copysign(x,y)
  52.  * double x,y;
  53.  */
  54.     .globl    _copysign
  55.     .text
  56.     .align    1
  57. _copysign:
  58.     .word    0x4
  59.     movq    4(ap),r0        # load x into r0
  60.     bicw3    $0x807f,r0,r2        # mask off the exponent of x
  61.     beql    Lz            # if zero or reserved op then return x
  62.     bicw3    $0x7fff,12(ap),r2    # copy the sign bit of y into r2
  63.     bicw2    $0x8000,r0        # replace x by |x|
  64.     bisw2    r2,r0            # copy the sign bit of y to x
  65. Lz:    ret
  66.  
  67. /*
  68.  * double logb(x)
  69.  * double x;
  70.  */
  71.     .globl    _logb
  72.     .text
  73.     .align    1
  74. _logb:
  75.     .word    0x0
  76.     bicl3    $0xffff807f,4(ap),r0    # mask off the exponent of x
  77.     beql    Ln
  78.     ashl    $-7,r0,r0        # get the bias exponent
  79.     subl2    $129,r0            # get the unbias exponent
  80.     cvtld    r0,r0            # return the answer in double
  81.     ret
  82. Ln:    movq    4(ap),r0        # r0:1 = x (zero or reserved op)
  83.     bneq    1f            # simply return if reserved op
  84.     movq     $0x0000fe00ffffcfff,r0  # -2147483647.0
  85. 1:    ret
  86.  
  87. /*
  88.  * long finite(x)
  89.  * double x;
  90.  */
  91.     .globl    _finite
  92.     .text
  93.     .align    1
  94. _finite:
  95.     .word    0x0000
  96.     bicw3    $0x7f,4(ap),r0        # mask off the mantissa
  97.     cmpw    r0,$0x8000        # to see if x is the reserved op
  98.     beql    1f            # if so, return FALSE (0)
  99.     movl    $1,r0            # else return TRUE (1)
  100.     ret
  101. 1:    clrl    r0
  102.     ret
  103.  
  104. /*
  105.  * double scalb(x,N)
  106.  * double x; int N;
  107.  */
  108.     .globl    _scalb
  109.     .set    ERANGE,34
  110.     .text
  111.     .align    1
  112. _scalb:
  113.     .word    0xc
  114.     movq    4(ap),r0
  115.     bicl3    $0xffff807f,r0,r3
  116.     beql    ret1            # 0 or reserved operand
  117.     movl    12(ap),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,_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
  136.  
  137. /*
  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.     .globl    _drem
  144.     .set    EDOM,33
  145.     .text
  146.     .align    1
  147. _drem:
  148.     .word    0xffc
  149.     subl2    $12,sp    
  150.     movq    4(ap),r0        #r0=x
  151.     movq    12(ap),r2        #r2=y
  152.     jeql    Rop            #if y=0 then generate reserved op fault
  153.     bicw3    $0x007f,r0,r4        #check if x is Rop
  154.     cmpw    r4,$0x8000
  155.     jeql    Ret            #if x is Rop then return Rop
  156.     bicl3    $0x007f,r2,r4        #check if y is Rop
  157.     cmpw    r4,$0x8000
  158.     jeql    Ret            #if y is Rop then return Rop
  159.     bicw2    $0x8000,r2        #y  := |y|
  160.     movw    $0,-4(fp)        #-4(fp) = nx := 0
  161.     cmpw    r2,$0x1c80        #yexp ? 57 
  162.     bgtr    C1            #if yexp > 57 goto C1
  163.     addw2    $0x1c80,r2        #scale up y by 2**57
  164.     movw    $0x1c80,-4(fp)        #nx := 57 (exponent field)
  165. C1:
  166.     movw    -4(fp),-8(fp)        #-8(fp) = nf := nx
  167.     bicw3    $0x7fff,r0,-12(fp)    #-12(fp) = sign of x
  168.     bicw2    $0x8000,r0        #x  := |x|
  169.     movq    r2,r10            #y1 := y
  170.     bicl2    $0xffff07ff,r11        #clear the last 27 bits of y1
  171. loop:
  172.     cmpd    r0,r2            #x ? y
  173.     bleq    E1            #if x <= y goto E1
  174.  /* begin argument reduction */
  175.     movq    r2,r4            #t =y
  176.     movq    r10,r6            #t1=y1
  177.     bicw3    $0x807f,r0,r8        #xexp= exponent of x
  178.     bicw3    $0x807f,r2,r9        #yexp= exponent fo y
  179.     subw2    r9,r8            #xexp-yexp
  180.     subw2    $0x0c80,r8        #k=xexp-yexp-25(exponent bit field)
  181.     blss    C2            #if k<0 goto C2
  182.     addw2    r8,r4            #t +=k    
  183.     addw2    r8,r6            #t1+=k, scale up t and t1
  184. C2:
  185.     divd3    r4,r0,r8        #x/t
  186.     cvtdl    r8,r8            #n=[x/t] truncated
  187.     cvtld    r8,r8            #float(n)
  188.     subd2    r6,r4            #t:=t-t1
  189.     muld2    r8,r4            #n*(t-t1)
  190.     muld2    r8,r6            #n*t1
  191.     subd2    r6,r0            #x-n*t1
  192.     subd2    r4,r0            #(x-n*t1)-n*(t-t1)
  193.     brb    loop
  194. E1:
  195.     movw    -4(fp),r6        #r6=nx
  196.     beql    C3            #if nx=0 goto C3
  197.     addw2    r6,r0            #x:=x*2**57 scale up x by nx
  198.     movw    $0,-4(fp)        #clear nx
  199.     brb    loop
  200. C3:
  201.     movq    r2,r4            #r4 = y
  202.     subw2    $0x80,r4        #r4 = y/2
  203.     cmpd    r0,r4            #x:y/2
  204.     blss    E2            #if x < y/2 goto E2
  205.     bgtr    C4            #if x > y/2 goto C4
  206.     cvtdl    r8,r8            #ifix(float(n))
  207.     blbc    r8,E2            #if the last bit is zero, goto E2
  208. C4:
  209.     subd2    r2,r0            #x-y
  210. E2:
  211.     xorw2    -12(fp),r0        #x^sign (exclusive or)
  212.     movw    -8(fp),r6        #r6=nf
  213.     bicw3    $0x807f,r0,r8        #r8=exponent of x
  214.     bicw2    $0x7f80,r0        #clear the exponent of x
  215.     subw2    r6,r8            #r8=xexp-nf
  216.     bgtr    C5            #if xexp-nf is positive goto C5
  217.     movw    $0,r8            #clear r8
  218.     movq    $0,r0            #x underflow to zero
  219. C5:
  220.     bisw2    r8,r0            #put r8 into x's exponent field
  221.     ret
  222. Rop:                    #Reserved operand
  223.     pushl    $EDOM
  224.     calls    $1,_infnan        #generate reserved op fault
  225.     ret
  226. Ret:
  227.     movq    $0x8000,r0        #propagate reserved op
  228.     ret
  229.