1 /*
   2  * CDDL HEADER START
   3  *
   4  * The contents of this file are subject to the terms of the
   5  * Common Development and Distribution License (the "License").
   6  * You may not use this file except in compliance with the License.
   7  *
   8  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
   9  * or http://www.opensolaris.org/os/licensing.
  10  * See the License for the specific language governing permissions
  11  * and limitations under the License.
  12  *
  13  * When distributing Covered Code, include this CDDL HEADER in each
  14  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
  15  * If applicable, add the following below this CDDL HEADER, with the
  16  * fields enclosed by brackets "[]" replaced with your own identifying
  17  * information: Portions Copyright [yyyy] [name of copyright owner]
  18  *
  19  * CDDL HEADER END
  20  */
  21 /*
  22  * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
  23  */
  24 /*
  25  * Copyright 2006 Sun Microsystems, Inc.  All rights reserved.
  26  * Use is subject to license terms.
  27  */
  28 
  29         .file "hypot.s"
  30 
  31 #include "libm.h"
  32 LIBM_ANSI_PRAGMA_WEAK(hypot,function)
  33 #include "libm_protos.h"
  34 
  35 #undef fabs
  36 
  37         .data
  38         .align  4
  39 inf:
  40         .long   0x7f800000
  41 
  42         ENTRY(hypot)
  43         movl    8(%esp),%eax            / eax <-- hi_32(x)
  44         andl    $0x7fffffff,%eax        / eax <-- hi_32(|x|)
  45         jz      .x_maybe_0              / if x = +/-0, return |y|
  46         subl    $0x7ff00000,%eax        / eax <-- hi_32(|x|) - hi_32(INF)
  47         jz      .x_maybe_inf
  48 .check_y:
  49         movl    16(%esp),%eax           / eax <-- hi_32(y)
  50         andl    $0x7fffffff,%eax        / eax <-- hi_32(|y|)
  51         jz      .y_maybe_0              / if y = +/-0, return |x|
  52         subl    $0x7ff00000,%eax        / eax <-- hi_32(|y|) - hi_32(INF)
  53         jz      .y_maybe_inf
  54 .do_hypot:
  55         fldl    12(%esp)                / ,y
  56         fmul    %st(0),%st              / ,y*y
  57         fldl    4(%esp)                 / x,y*y
  58         fmul    %st(0),%st              / x*x,y*y
  59         faddp   %st,%st(1)              / x*x+y*y
  60         fsqrt                           / sqrt(x*x+y*y)
  61         subl    $8,%esp
  62         fstpl   (%esp)                  / round to double
  63         fldl    (%esp)                  / sqrt(x*x+y*y) rounded to double
  64         PIC_SETUP(1)
  65         flds    PIC_L(inf)              / inf , sqrt(x*x+y*y)
  66         PIC_WRAPUP
  67         addl    $8,%esp
  68         fucomp
  69         fstsw   %ax                     / store status in %ax
  70         sahf                            / 80387 flags in %ah to 80386 flags
  71         jz      .maybe_ovflw
  72         ret
  73 
  74 .maybe_ovflw:
  75         jnp     .ovflw
  76         ret
  77 
  78 .ovflw:
  79         /       overflow occurred
  80         fstp    %st(0)                  / stack empty
  81         pushl   %ebp
  82         movl    %esp,%ebp
  83         PIC_SETUP(2)
  84         pushl   $4
  85         pushl   20(%ebp)                / high y
  86         pushl   16(%ebp)                / low y
  87         pushl   12(%ebp)                / high x
  88         pushl   8(%ebp)                 / low x
  89         call    PIC_F(_SVID_libm_err)
  90         addl    $20,%esp
  91         PIC_WRAPUP
  92         leave
  93         ret
  94 
  95 .x_maybe_0:
  96         movl    4(%esp),%ecx            / ecx <-- lo_32(x)
  97         orl     %ecx,%eax               / is x = +/-0?
  98         jnz     .check_y                / branch if x is denormal
  99         /  x = +/-0, so return |y|
 100         fldl    12(%esp)
 101         fabs
 102         ret
 103 
 104 .x_maybe_inf:
 105         movl    4(%esp),%ecx            / ecx <-- lo_32(x)
 106         orl     %ecx,%eax               / is x = +/-INF?
 107         jnz     .check_y                / branch if x is NaN
 108         / push&pop y in case y is a SNaN
 109         fldl    12(%esp)
 110         fstp    %st(0)
 111         / x = +/-INF, so return |x|
 112         fldl    4(%esp)
 113         fabs
 114         ret
 115 
 116 .y_maybe_0:
 117         movl    12(%esp),%ecx           / ecx <-- lo_32(y)
 118         orl     %ecx,%eax               / is y = +/-0?
 119         jnz     .do_hypot               / branch if y is denormal
 120         /  y = +/-0, so return |x|
 121         fldl    4(%esp)
 122         fabs
 123         ret
 124 
 125 .y_maybe_inf:
 126         movl    12(%esp),%ecx           / ecx <-- lo_32(y)
 127         orl     %ecx,%eax               / is y = +/-INF?
 128         jnz     .do_hypot               / branch if y is NaN
 129         / push&pop x in case x is a SNaN
 130         fldl    4(%esp)
 131         fstp    %st(0)
 132         /  y = +/-INF, so return |y|
 133         fldl    12(%esp)
 134         fabs
 135         ret
 136         .align  4
 137         SET_SIZE(hypot)