?login_element?

Subversion Repositories NedoOS

Rev

Blame | Last modification | View Log | Download

  1. /*-------------------------------------------------------------------------
  2.    _fsdiv.c - Floating point library in optimized assembly for 8051
  3.  
  4.    Copyright (c) 2004, Paul Stoffregen, paul@pjrc.com
  5.  
  6.    This library is free software; you can redistribute it and/or modify it
  7.    under the terms of the GNU General Public License as published by the
  8.    Free Software Foundation; either version 2, or (at your option) any
  9.    later version.
  10.  
  11.    This library is distributed in the hope that it will be useful,
  12.    but WITHOUT ANY WARRANTY; without even the implied warranty of
  13.    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14.    GNU General Public License for more details.
  15.  
  16.    You should have received a copy of the GNU General Public License
  17.    along with this library; see the file COPYING. If not, write to the
  18.    Free Software Foundation, 51 Franklin Street, Fifth Floor, Boston,
  19.    MA 02110-1301, USA.
  20.  
  21.    As a special exception, if you link this library with other files,
  22.    some of which are compiled with SDCC, to produce an executable,
  23.    this library does not by itself cause the resulting executable to
  24.    be covered by the GNU General Public License. This exception does
  25.    not however invalidate any other reasons why the executable file
  26.    might be covered by the GNU General Public License.
  27. -------------------------------------------------------------------------*/
  28.  
  29.  
  30. #define __SDCC_FLOAT_LIB
  31. #include <float.h>
  32.  
  33.  
  34. #ifdef FLOAT_ASM_MCS51
  35.  
  36. // float __fsdiv (float a, float b) __reentrant
  37. static void dummy(void) __naked
  38. {
  39.         __asm
  40.         .globl  ___fsdiv
  41. ___fsdiv:
  42.         // extract the two inputs, placing them into:
  43.         //      sign     exponent   mantiassa
  44.         //      ----     --------   ---------
  45.         //  a:  sign_a   exp_a      r4/r3/r2
  46.         //  b:  sign_b   exp_b      r7/r6/r5
  47.  
  48.         lcall   fsgetargs
  49.  
  50.         // compute final sign bit
  51.         jnb     sign_b, 00001$
  52.         cpl     sign_a
  53. 00001$:
  54.  
  55.         // if divisor is zero, ...
  56.         cjne    r7, #0, 00003$
  57.         // if dividend is also zero, return NaN
  58.         cjne    r4, #0, 00002$
  59.         ljmp    fs_return_nan
  60. 00002$:
  61.         // but dividend is non-zero, return infinity
  62.         ljmp    fs_return_inf
  63. 00003$:
  64.         // if dividend is zero, return zero
  65.         cjne    r4, #0, 00004$
  66.         ljmp    fs_return_zero
  67. 00004$:
  68.         // if divisor is infinity, ...
  69.         mov     a, exp_b
  70.         cjne    a, #0xFF, 00006$
  71.         // and dividend is also infinity, return NaN
  72.         mov     a, exp_a
  73.         cjne    a, #0xFF, 00005$
  74.         ljmp    fs_return_nan
  75. 00005$:
  76.         // but dividend is not infinity, return zero
  77.         ljmp    fs_return_zero
  78. 00006$:
  79.         // if dividend is infinity, return infinity
  80.         mov     a, exp_a
  81.         cjne    a, #0xFF, 00007$
  82.         ljmp    fs_return_inf
  83. 00007$:
  84.  
  85.         // subtract exponents
  86.         clr     c
  87.         subb    a, exp_b
  88.         // if no carry then no underflow
  89.         jnc     00008$
  90.         add     a, #127
  91.         jc      00009$
  92.         ljmp    fs_return_zero
  93.  
  94. 00008$:
  95.         add     a, #128
  96.         dec     a
  97.         jnc     00009$
  98.         ljmp    fs_return_inf
  99.  
  100. 00009$:
  101.         mov     exp_a, a
  102.  
  103.         // need extra bits on a's mantissa
  104. #ifdef FLOAT_FULL_ACCURACY
  105.         clr     c
  106.         mov     a, r5
  107.         subb    a, r2
  108.         mov     a, r6
  109.         subb    a, r3
  110.         mov     a, r7
  111.         subb    a, r4
  112.         jc      00010$
  113.         dec     exp_a
  114.         clr     c
  115.         mov     a, r2
  116.         rlc     a
  117.         mov     r1, a
  118.         mov     a, r3
  119.         rlc     a
  120.         mov     r2, a
  121.         mov     a, r4
  122.         rlc     a
  123.         mov     r3, a
  124.         clr     a
  125.         rlc     a
  126.         mov     r4, a
  127.         sjmp    00011$
  128. 00010$:
  129. #endif
  130.         clr     a
  131.         xch     a, r4
  132.         xch     a, r3
  133.         xch     a, r2
  134.         mov     r1, a
  135. 00011$:
  136.  
  137.         // begin long division
  138.         push    exp_a
  139. #ifdef FLOAT_FULL_ACCURACY
  140.         mov     b, #25
  141. #else
  142.         mov     b, #24
  143. #endif
  144. 00012$:
  145.         // compare
  146.         clr     c
  147.         mov     a, r1
  148.         subb    a, r5
  149.         mov     a, r2
  150.         subb    a, r6
  151.         mov     a, r3
  152.         subb    a, r7
  153.         mov     a, r4
  154.         subb    a, #0           // carry==0 if mant1 >= mant2
  155.  
  156. #ifdef FLOAT_FULL_ACCURACY
  157.         djnz    b, 00013$
  158.         sjmp    00015$
  159. 00013$:
  160. #endif
  161.         jc      00014$
  162.         // subtract
  163.         mov     a, r1
  164.         subb    a, r5
  165.         mov     r1, a
  166.         mov     a, r2
  167.         subb    a, r6
  168.         mov     r2, a
  169.         mov     a, r3
  170.         subb    a, r7
  171.         mov     r3, a
  172.         mov     a, r4
  173.         subb    a, #0
  174.         mov     r4, a
  175.         clr     c
  176.  
  177. 00014$:
  178.         // shift result
  179.         cpl     c
  180.         mov     a, r0
  181.         rlc     a
  182.         mov     r0, a
  183.         mov     a, dpl
  184.         rlc     a
  185.         mov     dpl, a
  186.         mov     a, dph
  187.         rlc     a
  188.         mov     dph, a
  189.  
  190.         // shift partial remainder
  191.         clr     c
  192.         mov     a, r1
  193.         rlc     a
  194.         mov     r1, a
  195.         mov     a, r2
  196.         rlc     a
  197.         mov     r2, a
  198.         mov     a, r3
  199.         rlc     a
  200.         mov     r3, a
  201.         mov     a, r4
  202.         rlc     a
  203.         mov     r4, a
  204.  
  205. #ifdef FLOAT_FULL_ACCURACY
  206.         sjmp    00012$
  207. 00015$:
  208. #else
  209.         djnz    b, 00012$
  210. #endif
  211.  
  212.         // now we've got a division result, so all we need to do
  213.         // is round off properly, normalize and output a float
  214.  
  215. #ifdef FLOAT_FULL_ACCURACY
  216.         cpl     c
  217.         clr     a
  218.         mov     r1, a
  219.         addc    a, r0
  220.         mov     r2, a
  221.         clr     a
  222.         addc    a, dpl
  223.         mov     r3, a
  224.         clr     a
  225.         addc    a, dph
  226.         mov     r4, a
  227.         pop     exp_a
  228.         jnc     00016$
  229.         inc     exp_a
  230.         // incrementing exp_a without checking carry is dangerous
  231.         mov     r4, #0x80
  232. 00016$:
  233. #else
  234.         mov     r1, #0
  235.         mov     a, r0
  236.         mov     r2, a
  237.         mov     r3, dpl
  238.         mov     r4, dph
  239.         pop     exp_a
  240. #endif
  241.  
  242.         lcall   fs_normalize_a
  243.         ljmp    fs_zerocheck_return
  244.         __endasm;
  245. }
  246.  
  247. #else
  248.  
  249. /*
  250. ** libgcc support for software floating point.
  251. ** Copyright (C) 1991 by Pipeline Associates, Inc.  All rights reserved.
  252. ** Permission is granted to do *anything* you want with this file,
  253. ** commercial or otherwise, provided this message remains intact.  So there!
  254. ** I would appreciate receiving any updates/patches/changes that anyone
  255. ** makes, and am willing to be the repository for said changes (am I
  256. ** making a big mistake?).
  257. **
  258. ** Pat Wood
  259. ** Pipeline Associates, Inc.
  260. ** pipeline!phw@motown.com or
  261. ** sun!pipeline!phw or
  262. ** uunet!motown!pipeline!phw
  263. */
  264.  
  265. /* (c)2000/2001: hacked a little by johan.knol@iduna.nl for sdcc */
  266.  
  267. union float_long
  268.   {
  269.     float f;
  270.     long l;
  271.   };
  272.  
  273. /* divide two floats */
  274. static float __fsdiv_org (float a1, float a2)
  275. {
  276.   volatile union float_long fl1, fl2;
  277.   volatile long result;
  278.   volatile unsigned long mask;
  279.   volatile long mant1, mant2;
  280.   volatile int exp;
  281.   char sign;
  282.  
  283.   fl1.f = a1;
  284.   fl2.f = a2;
  285.  
  286.   /* subtract exponents */
  287.   exp = EXP (fl1.l) ;
  288.   exp -= EXP (fl2.l);
  289.   exp += EXCESS;
  290.  
  291.   /* compute sign */
  292.   sign = SIGN (fl1.l) ^ SIGN (fl2.l);
  293.  
  294.   /* divide by zero??? */
  295.   if (!fl2.l)
  296.     {/* return NaN or -NaN */
  297.       fl2.l = 0x7FC00000;
  298.       return (fl2.f);
  299.     }
  300.  
  301.   /* numerator zero??? */
  302.   if (!fl1.l)
  303.     return (0);
  304.  
  305.   /* now get mantissas */
  306.   mant1 = MANT (fl1.l);
  307.   mant2 = MANT (fl2.l);
  308.  
  309.   /* this assures we have 25 bits of precision in the end */
  310.   if (mant1 < mant2)
  311.     {
  312.       mant1 <<= 1;
  313.       exp--;
  314.     }
  315.  
  316.   /* now we perform repeated subtraction of fl2.l from fl1.l */
  317.   mask = 0x1000000;
  318.   result = 0;
  319.   while (mask)
  320.     {
  321.       if (mant1 >= mant2)
  322.         {
  323.           result |= mask;
  324.           mant1 -= mant2;
  325.         }
  326.       mant1 <<= 1;
  327.       mask >>= 1;
  328.     }
  329.  
  330.   /* round */
  331.   result += 1;
  332.  
  333.   /* normalize down */
  334.   exp++;
  335.   result >>= 1;
  336.  
  337.   result &= ~HIDDEN;
  338.  
  339.   /* pack up and go home */
  340.   if (exp >= 0x100)
  341.     fl1.l = (sign ? SIGNBIT : 0) | __INFINITY;
  342.   else if (exp < 0)
  343.     fl1.l = 0;
  344.   else
  345.     fl1.l = PACK (sign ? SIGNBIT : 0 , exp, result);
  346.   return (fl1.f);
  347. }
  348.  
  349. float __fsdiv (float a1, float a2)
  350. {
  351.   float f;
  352.   unsigned long *p = (unsigned long *) &f;
  353.  
  354.   if (a2 == 0.0f && a1 > 0.0f)
  355.     *p = 0x7f800000; // inf
  356.   else if (a2 == 0.0f && a1 < 0.0f)
  357.     *p = 0xff800000; // -inf
  358.   else if (a2 == 0.0f && a1 == 0.0f)
  359.     *p = 0xffc00000; // nan
  360.   else
  361.     f = __fsdiv_org (a1, a2);
  362.  
  363.   return f;
  364. }
  365.  
  366. #endif
  367.