Subversion Repositories HelenOS-historic

Rev

Rev 660 | Go to most recent revision | Blame | 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. #include<sftypes.h>
  30. #include<arithmetic.h>
  31. #include<comparison.h>
  32.  
  33. /** Add two Float32 numbers with same signs
  34.  */
  35. float32 addFloat32(float32 a, float32 b)
  36. {
  37.     int expdiff;
  38.     __u32 exp1,exp2,mant1,mant2;
  39.    
  40.     expdiff=a.parts.exp - b.parts.exp;
  41.     if (expdiff<0) {
  42.         if (isFloat32NaN(b)) {
  43.             //TODO: fix SigNaN
  44.             if (isFloat32SigNaN(b)) {
  45.             };
  46.             return b;
  47.         };
  48.        
  49.         if (b.parts.exp==0xFF) {
  50.             return b;
  51.         }
  52.        
  53.         mant1=b.parts.mantisa;
  54.         exp1=b.parts.exp;
  55.         mant2=a.parts.mantisa;
  56.         exp2=a.parts.exp;
  57.         expdiff*=-1;
  58.     } else {
  59.         if (isFloat32NaN(a)) {
  60.             //TODO: fix SigNaN
  61.             if ((isFloat32SigNaN(a))||(isFloat32SigNaN(b))) {
  62.             };
  63.             return a;
  64.         };
  65.        
  66.         if (a.parts.exp==0xFF) {
  67.             return a;
  68.         }
  69.        
  70.         mant1=a.parts.mantisa;
  71.         exp1=a.parts.exp;
  72.         mant2=b.parts.mantisa;
  73.         exp2=b.parts.exp;
  74.     };
  75.    
  76.     if (exp1==0) {
  77.         //both are denormalized
  78.         mant1+=mant2;
  79.         if (mant1&0xF00000) {
  80.             a.parts.exp=1;
  81.         };
  82.         a.parts.mantisa=mant1;
  83.         return a;
  84.     };
  85.    
  86.     // create some space for rounding
  87.     mant1<<=6;
  88.     mant2<<=6;
  89.    
  90.     mant1|=0x20000000; //add hidden bit
  91.    
  92.    
  93.     if (exp2==0) {
  94.         --expdiff; 
  95.     } else {
  96.         mant2|=0x20000000; //hidden bit
  97.     };
  98.    
  99.     if (expdiff>24) {
  100.          goto done;
  101.          };
  102.    
  103.     mant2>>=expdiff;
  104.     mant1+=mant2;
  105. done:
  106.     if (mant1&0x40000000) {
  107.         ++exp1;
  108.         mant1>>=1;
  109.     };
  110.    
  111.     //rounding - if first bit after mantisa is set then round up
  112.     mant1+=0x20;
  113.    
  114.     if (mant1&0x40000000) {
  115.         ++exp1;
  116.         mant1>>=1;
  117.     };
  118.    
  119.     a.parts.exp=exp1;
  120.     a.parts.mantisa = ((mant1&(~0x20000000))>>6); /*Clear hidden bit and shift */
  121.     return a;
  122. };
  123.  
  124. /** Subtract two float32 numbers with same signs
  125.  */
  126. float32 subFloat32(float32 a, float32 b)
  127. {
  128.     int expdiff;
  129.     __u32 exp1,exp2,mant1,mant2;
  130.     float32 result;
  131.  
  132.     result.f = 0;
  133.    
  134.     expdiff=a.parts.exp - b.parts.exp;
  135.     if ((expdiff<0)||((expdiff==0)&&(a.parts.mantisa<b.parts.mantisa))) {
  136.         if (isFloat32NaN(b)) {
  137.             //TODO: fix SigNaN
  138.             if (isFloat32SigNaN(b)) {
  139.             };
  140.             return b;
  141.         };
  142.        
  143.         if (b.parts.exp==0xFF) {
  144.             b.parts.sign=!b.parts.sign; /* num -(+-inf) = -+inf */
  145.             return b;
  146.         }
  147.        
  148.         result.parts.sign = !a.parts.sign;
  149.        
  150.         mant1=b.parts.mantisa;
  151.         exp1=b.parts.exp;
  152.         mant2=a.parts.mantisa;
  153.         exp2=a.parts.exp;
  154.         expdiff*=-1;
  155.     } else {
  156.         if (isFloat32NaN(a)) {
  157.             //TODO: fix SigNaN
  158.             if ((isFloat32SigNaN(a))||(isFloat32SigNaN(b))) {
  159.             };
  160.             return a;
  161.         };
  162.        
  163.         if (a.parts.exp==0xFF) {
  164.             if (b.parts.exp==0xFF) {
  165.                 /* inf - inf => nan */
  166.                 //TODO: fix exception
  167.                 result.binary = FLOAT32_NAN;
  168.                 return result;
  169.             };
  170.             return a;
  171.         }
  172.        
  173.         result.parts.sign = a.parts.sign;
  174.        
  175.         mant1=a.parts.mantisa;
  176.         exp1=a.parts.exp;
  177.         mant2=b.parts.mantisa;
  178.         exp2=b.parts.exp;
  179.        
  180.  
  181.        
  182.     };
  183.    
  184.     if (exp1==0) {
  185.         //both are denormalized
  186.         result.parts.mantisa=mant1-mant2;
  187.         if (result.parts.mantisa>mant1) {
  188.             //TODO: underflow exception
  189.             return result;
  190.         };
  191.         result.parts.exp=0;
  192.         return result;
  193.     };
  194.    
  195.     // create some space for rounding
  196.     mant1<<=6;
  197.     mant2<<=6;
  198.    
  199.     mant1|=0x20000000; //add hidden bit
  200.    
  201.    
  202.     if (exp2==0) {
  203.         --expdiff; 
  204.     } else {
  205.         mant2|=0x20000000; //hidden bit
  206.     };
  207.    
  208.     if (expdiff>24) {
  209.          goto done;
  210.          };
  211.    
  212.     mant1 = mant1-(mant2>>expdiff);
  213. done:
  214.    
  215.     //TODO: find first nonzero digit and shift result and detect possibly underflow
  216.     while ((exp1>0)&&(!(mant1&0x20000000))) {
  217.         exp1--;
  218.         mant1 <<= 1;
  219.         if(mant1 == 0) {
  220.             /* Realy is it an underflow? ... */
  221.             //TODO: fix underflow
  222.         };
  223.     };
  224.    
  225.     //rounding - if first bit after mantisa is set then round up   
  226.     mant1 += 0x20;
  227.  
  228.     if (mant1&0x40000000) {
  229.         ++exp1;
  230.         mant1>>=1;
  231.     };
  232.    
  233.     result.parts.mantisa = ((mant1&(~0x20000000))>>6); /*Clear hidden bit and shift */
  234.     result.parts.exp = exp1;
  235.    
  236.     return result;
  237. };
  238.  
  239.