Subversion Repositories HelenOS

Rev

Rev 1787 | Go to most recent revision | Blame | Compare with Previous | Last modification | View Log | Download | RSS feed

  1. /*
  2.  * Copyright (c) 2005 Josef Cejka
  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.  *
  9.  * - Redistributions of source code must retain the above copyright
  10.  *   notice, this list of conditions and the following disclaimer.
  11.  * - 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.  * - The name of the author may not be used to endorse or promote products
  15.  *   derived from this software without specific prior written permission.
  16.  *
  17.  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
  18.  * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
  19.  * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
  20.  * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
  21.  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
  22.  * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
  23.  * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
  24.  * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
  25.  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
  26.  * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  27.  */
  28.  
  29. /** @addtogroup softfloat  
  30.  * @{
  31.  */
  32. /** @file
  33.  */
  34.  
  35. #include<sftypes.h>
  36. #include<add.h>
  37. #include<comparison.h>
  38.  
  39. /** Add two Float32 numbers with same signs
  40.  */
  41. float32 addFloat32(float32 a, float32 b)
  42. {
  43.     int expdiff;
  44.     uint32_t exp1, exp2,frac1, frac2;
  45.    
  46.     expdiff = a.parts.exp - b.parts.exp;
  47.     if (expdiff < 0) {
  48.         if (isFloat32NaN(b)) {
  49.             /* TODO: fix SigNaN */
  50.             if (isFloat32SigNaN(b)) {
  51.             };
  52.  
  53.             return b;
  54.         };
  55.        
  56.         if (b.parts.exp == FLOAT32_MAX_EXPONENT) {
  57.             return b;
  58.         }
  59.        
  60.         frac1 = b.parts.fraction;
  61.         exp1 = b.parts.exp;
  62.         frac2 = a.parts.fraction;
  63.         exp2 = a.parts.exp;
  64.         expdiff *= -1;
  65.     } else {
  66.         if ((isFloat32NaN(a)) || (isFloat32NaN(b))) {
  67.             /* TODO: fix SigNaN */
  68.             if (isFloat32SigNaN(a) || isFloat32SigNaN(b)) {
  69.             };
  70.             return (isFloat32NaN(a)?a:b);
  71.         };
  72.        
  73.         if (a.parts.exp == FLOAT32_MAX_EXPONENT) {
  74.             return a;
  75.         }
  76.        
  77.         frac1 = a.parts.fraction;
  78.         exp1 = a.parts.exp;
  79.         frac2 = b.parts.fraction;
  80.         exp2 = b.parts.exp;
  81.     };
  82.    
  83.     if (exp1 == 0) {
  84.         /* both are denormalized */
  85.         frac1 += frac2;
  86.         if (frac1 & FLOAT32_HIDDEN_BIT_MASK ) {
  87.             /* result is not denormalized */
  88.             a.parts.exp = 1;
  89.         };
  90.         a.parts.fraction = frac1;
  91.         return a;
  92.     };
  93.    
  94.     frac1 |= FLOAT32_HIDDEN_BIT_MASK; /* add hidden bit */
  95.  
  96.     if (exp2 == 0) {
  97.         /* second operand is denormalized */
  98.         --expdiff;
  99.     } else {
  100.         /* add hidden bit to second operand */
  101.         frac2 |= FLOAT32_HIDDEN_BIT_MASK;
  102.     };
  103.    
  104.     /* create some space for rounding */
  105.     frac1 <<= 6;
  106.     frac2 <<= 6;
  107.    
  108.     if (expdiff < (FLOAT32_FRACTION_SIZE + 2) ) {
  109.         frac2 >>= expdiff;
  110.         frac1 += frac2;
  111.     } else {
  112.         a.parts.exp = exp1;
  113.         a.parts.fraction = (frac1 >> 6) & (~(FLOAT32_HIDDEN_BIT_MASK));
  114.         return a;
  115.     }
  116.    
  117.     if (frac1 & (FLOAT32_HIDDEN_BIT_MASK << 7) ) {
  118.         ++exp1;
  119.         frac1 >>= 1;
  120.     };
  121.    
  122.     /* rounding - if first bit after fraction is set then round up */
  123.     frac1 += (0x1 << 5);
  124.    
  125.     if (frac1 & (FLOAT32_HIDDEN_BIT_MASK << 7)) {
  126.         /* rounding overflow */
  127.         ++exp1;
  128.         frac1 >>= 1;
  129.     };
  130.    
  131.    
  132.     if ((exp1 == FLOAT32_MAX_EXPONENT ) || (exp2 > exp1)) {
  133.             /* overflow - set infinity as result */
  134.             a.parts.exp = FLOAT32_MAX_EXPONENT;
  135.             a.parts.fraction = 0;
  136.             return a;
  137.             }
  138.    
  139.     a.parts.exp = exp1;
  140.    
  141.     /*Clear hidden bit and shift */
  142.     a.parts.fraction = ((frac1 >> 6) & (~FLOAT32_HIDDEN_BIT_MASK)) ;
  143.     return a;
  144. }
  145.  
  146.  
  147. /** Add two Float64 numbers with same signs
  148.  */
  149. float64 addFloat64(float64 a, float64 b)
  150. {
  151.     int expdiff;
  152.     uint32_t exp1, exp2;
  153.     uint64_t frac1, frac2;
  154.    
  155.     expdiff = ((int )a.parts.exp) - b.parts.exp;
  156.     if (expdiff < 0) {
  157.         if (isFloat64NaN(b)) {
  158.             /* TODO: fix SigNaN */
  159.             if (isFloat64SigNaN(b)) {
  160.             };
  161.  
  162.             return b;
  163.         };
  164.        
  165.         /* b is infinity and a not */  
  166.         if (b.parts.exp == FLOAT64_MAX_EXPONENT ) {
  167.             return b;
  168.         }
  169.        
  170.         frac1 = b.parts.fraction;
  171.         exp1 = b.parts.exp;
  172.         frac2 = a.parts.fraction;
  173.         exp2 = a.parts.exp;
  174.         expdiff *= -1;
  175.     } else {
  176.         if (isFloat64NaN(a)) {
  177.             /* TODO: fix SigNaN */
  178.             if (isFloat64SigNaN(a) || isFloat64SigNaN(b)) {
  179.             };
  180.             return a;
  181.         };
  182.        
  183.         /* a is infinity and b not */
  184.         if (a.parts.exp == FLOAT64_MAX_EXPONENT ) {
  185.             return a;
  186.         }
  187.        
  188.         frac1 = a.parts.fraction;
  189.         exp1 = a.parts.exp;
  190.         frac2 = b.parts.fraction;
  191.         exp2 = b.parts.exp;
  192.     };
  193.    
  194.     if (exp1 == 0) {
  195.         /* both are denormalized */
  196.         frac1 += frac2;
  197.         if (frac1 & FLOAT64_HIDDEN_BIT_MASK) {
  198.             /* result is not denormalized */
  199.             a.parts.exp = 1;
  200.         };
  201.         a.parts.fraction = frac1;
  202.         return a;
  203.     };
  204.    
  205.     /* add hidden bit - frac1 is sure not denormalized */
  206.     frac1 |= FLOAT64_HIDDEN_BIT_MASK;
  207.  
  208.     /* second operand ... */
  209.     if (exp2 == 0) {
  210.         /* ... is denormalized */
  211.         --expdiff; 
  212.     } else {
  213.         /* is not denormalized */
  214.         frac2 |= FLOAT64_HIDDEN_BIT_MASK;
  215.     };
  216.    
  217.     /* create some space for rounding */
  218.     frac1 <<= 6;
  219.     frac2 <<= 6;
  220.    
  221.     if (expdiff < (FLOAT64_FRACTION_SIZE + 2) ) {
  222.         frac2 >>= expdiff;
  223.         frac1 += frac2;
  224.     } else {
  225.         a.parts.exp = exp1;
  226.         a.parts.fraction = (frac1 >> 6) & (~(FLOAT64_HIDDEN_BIT_MASK));
  227.         return a;
  228.     }
  229.    
  230.     if (frac1 & (FLOAT64_HIDDEN_BIT_MASK << 7) ) {
  231.         ++exp1;
  232.         frac1 >>= 1;
  233.     };
  234.    
  235.     /* rounding - if first bit after fraction is set then round up */
  236.     frac1 += (0x1 << 5);
  237.    
  238.     if (frac1 & (FLOAT64_HIDDEN_BIT_MASK << 7)) {
  239.         /* rounding overflow */
  240.         ++exp1;
  241.         frac1 >>= 1;
  242.     };
  243.    
  244.     if ((exp1 == FLOAT64_MAX_EXPONENT ) || (exp2 > exp1)) {
  245.             /* overflow - set infinity as result */
  246.             a.parts.exp = FLOAT64_MAX_EXPONENT;
  247.             a.parts.fraction = 0;
  248.             return a;
  249.             }
  250.    
  251.     a.parts.exp = exp1;
  252.     /*Clear hidden bit and shift */
  253.     a.parts.fraction = ( (frac1 >> 6 ) & (~FLOAT64_HIDDEN_BIT_MASK));
  254.    
  255.     return a;
  256. }
  257.  
  258. /** @}
  259.  */
  260.