-/*-------------------------------------------------------------------------
- math.h - ANSI functions forward declarations
-
- Written By - Sandeep Dutta . sandeep.dutta@usa.net (1998)
-
- This program is free software; you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by the
- Free Software Foundation; either version 2, or (at your option) any
- later version.
-
- This program is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- GNU General Public License for more details.
-
- You should have received a copy of the GNU General Public License
- along with this program; if not, write to the Free Software
- Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
-
- In other words, you are welcome to use, share and improve this program.
- You are forbidden to forbid anyone else to use, share and improve
- what you give them. Help stamp out software-hoarding!
--------------------------------------------------------------------------*/
-
-#ifndef __SDC51_MATH_H
-#define __SDC51_MATH_H 1
-#error Floating point not yet completely implemented
-#endif
+/* math.h: Floating point math function declarations
+ Copyright (C) 2001 Jesus Calvino-Fraga, jesusc@ieee.org
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+/* Version 1.0 - Initial release */
+
+#ifndef _INC_MATH
+#define _INC_MATH
+
+#define PI 3.1415926536
+#define TWO_PI 6.2831853071
+#define HALF_PI 1.5707963268
+#define QUART_PI 0.7853981634
+#define iPI 0.3183098862
+#define iTWO_PI 0.1591549431
+#define TWO_O_PI 0.6366197724
+
+// EPS=B**(-t/2), where B is the radix of the floating-point representation
+// and there are t base-B digits in the significand. Therefore, for floats
+// EPS=2**(-12). Also define EPS2=EPS*EPS.
+#define EPS 244.14062E-6
+#define EPS2 59.6046E-9
+#define XMAX 3.402823466E+38
+
+union float_long
+{
+ float f;
+ long l;
+};
+
+/**********************************************
+ * Prototypes for float ANSI C math functions *
+ **********************************************/
+
+/* Trigonometric functions */
+float sinf(const float x);
+float cosf(const float x);
+float tanf(const float x);
+float cotf(const float x);
+float asinf(const float x);
+float acosf(const float x);
+float atanf(const float x);
+float atan2f(const float x, const float y);
+
+/* Hyperbolic functions */
+float sinhf(const float x);
+float coshf(const float x);
+float tanhf(const float x);
+
+/* Exponential, logarithmic and power functions */
+float expf(const float x);
+float logf(const float x);
+float log10f(const float x);
+float powf(const float x, const float y);
+float sqrtf(const float a);
+
+/* Nearest integer, absolute value, and remainder functions */
+float fabsf(const float x);
+float frexpf(const float x, int *pw2);
+float ldexpf(const float x, const int pw2);
+float ceilf(float x);
+float floorf(float x);
+float modff(float x, float * y);
+
+#endif /* _INC_MATH */
_strstr.c _strtok.c _uchar2fs.c _uint2fs.c \
_ulong2fs.c malloc.c serial.c ser_ir.c printfl.c \
printf_large.c vprintf.c puts.c gets.c \
- assert.c _strcat.c time.c printf_fast.c bpx.c
+ assert.c _strcat.c time.c printf_fast.c bpx.c \
+ fabsf.c frexpf.c ldexpf.c expf.c powf.c sincosf.c sinf.c \
+ cosf.c logf.c log10f.c sqrtf.c tancotf.c tanf.c cotf.c \
+ asincosf.c asinf.c acosf.c atanf.c atan2f.c sincoshf.c \
+ sinhf.c coshf.c tanhf.c floorf.c ceilf.c modff.c
OBJECTS = $(patsubst %.c,$(PORTDIR)/%.rel,$(SOURCES))
--- /dev/null
+/* acosf.c: Computes arc cosine of a 32-bit float
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float asincosf(const float x, const int isacos);
+
+float acosf(const float x) reentrant
+{
+ if(x== 1.0) return 0.0;
+ else if(x==-1.0) return PI;
+ else if(x== 0.0) return HALF_PI;
+ return asincosf(x,1);
+}
--- /dev/null
+/* asincosf.c: Computes asin or acos of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P1 0.933935835E+0
+#define P2 -0.504400557E+0
+#define Q0 0.560363004E+1
+#define Q1 -0.554846723E+1
+#define Q2 0.100000000E+1
+
+#define P(g) (P2*g+P1)
+#define Q(g) ((Q2*g+Q1)*g+Q0)
+
+float asincosf(const float x, const int isacos)
+{
+ float y, g, r;
+ int i;
+ float a[2]={ 0.0, QUART_PI };
+ float b[2]={ HALF_PI, QUART_PI };
+
+ y=fabsf(x);
+ i=isacos;
+ if (y < EPS) r=y;
+ else
+ {
+ if (y > 0.5)
+ {
+ i=1-i;
+ if (y > 1.0)
+ {
+ errno=EDOM;
+ return 0.0;
+ }
+ g=(0.5-y)+0.5;
+ g=ldexpf(g,-1);
+ y=sqrtf(g);
+ y=-(y+y);
+ }
+ else
+ {
+ g=y*y;
+ }
+ r=y+y*((P(g)*g)/Q(g));
+ }
+ if (isacos)
+ {
+ if (x < 0.0)
+ r=(b[i]+r)+b[i];
+ else
+ r=(a[i]-r)+a[i];
+ }
+ else
+ {
+ r=(a[i]+r)+a[i];
+ if (x<0.0) r=-r;
+ }
+ return r;
+}
--- /dev/null
+/* asinf.c: Computes asin(x)
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float asincosf(const float x, const int isacos);
+
+float asinf(const float x) reentrant
+{
+ if(x== 1.0) return HALF_PI;
+ else if(x==-1.0) return -HALF_PI;
+ else if(x== 0.0) return 0.0;
+ else return asincosf(x,0);
+}
+
--- /dev/null
+/* atan2f.c: Computes atan2(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float atan2f(const float x, const float y)
+{
+ float r;
+
+ if ((x==0.0) && (y==0.0))
+ {
+ errno=EDOM;
+ return 0.0;
+ }
+
+ if(fabsf(y)>=fabsf(x))
+ {
+ r=atanf(x/y);
+ if(y<0.0) r+=(x>=0?PI:-PI);
+ }
+ else
+ {
+ r=-atanf(y/x);
+ r+=(x<0.0?-HALF_PI:HALF_PI);
+ }
+ return r;
+}
--- /dev/null
+/* atanf.c: Computes arctan of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P0 -0.4708325141E+0
+#define P1 -0.5090958253E-1
+#define Q0 0.1412500740E+1
+#define Q1 0.1000000000E+1
+
+#define P(g,f) ((P1*g+P0)*g*f)
+#define Q(g) (Q1*g+Q0)
+
+#define K1 0.2679491924 /* 2-sqrt(3) */
+#define K2 0.7320508076 /* sqrt(3)-1 */
+#define K3 1.7320508076 /* sqrt(3) */
+
+float atanf(const float x) reentrant
+{
+ float f, r, g;
+ int n=0;
+ static float a[]={ 0.0, 0.5235987756, 1.5707963268, 1.0471975512 };
+
+ f=fabsf(x);
+ if(f>1.0)
+ {
+ f=1.0/f;
+ n=2;
+ }
+ if(f>K1)
+ {
+ f=((K2*f-1.0)+f)/(K3+f);
+ // What it is actually wanted is this more accurate formula,
+ // but SDCC optimizes it and then it does not work:
+ // f=(((K2*f-0.5)-0.5)+f)/(K3+f);
+ n++;
+ }
+ if(fabsf(f)<EPS) r=f;
+ else
+ {
+ g=f*f;
+ r=f+P(g,f)/Q(g);
+ }
+ if(n>1) r=-r;
+ r+=a[n];
+ if(x<0.0) r=-r;
+ return r;
+}
+
--- /dev/null
+/* ceilf.c: Returns the integer larger or equal than x
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float ceilf(float x) reentrant
+{
+ long r;
+ r=x;
+ if (r<0)
+ return r;
+ else
+ return (r+((r<x)?1:0));
+}
--- /dev/null
+/* cosf.c: Computes cos(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float sincosf(const float x, const int iscos);
+
+float cosf(const float x) reentrant
+{
+ if (x==0.0) return 1.0;
+ return sincosf(x, 1);
+}
--- /dev/null
+/* coshf.c: Computes cosh(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float sincoshf(const float x, const int iscosh);
+
+float coshf(const float x) reentrant
+{
+ return sincoshf(x, 1);
+}
--- /dev/null
+/* cotf.c: Computes cot(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float tancotf(const float x, const int iscot);
+
+float cotf(const float x) reentrant
+{
+ float y;
+
+ y=fabsf(x);
+ if (y<1.0E-30) //This one requires more thinking...
+ {
+ errno = ERANGE;
+ if (x<0.0)
+ return -XMAX;
+ else
+ return XMAX;
+ }
+ return tancotf(x, 1);
+}
+
--- /dev/null
+/* expf.c: Computes e**x of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P0 0.2499999995E+0
+#define P1 0.4160288626E-2
+#define Q0 0.5000000000E+0
+#define Q1 0.4998717877E-1
+
+#define P(z) ((P1*z)+P0)
+#define Q(z) ((Q1*z)+Q0)
+
+#define C1 0.693359375
+#define C2 -2.1219444005469058277e-4
+
+#define BIGX 88.72283911 /* ln(XMAX) */
+#define EXPEPS 1.0E-7 /* exp(1.0E-7)=0.0000001 */
+#define K1 1.4426950409 /* 1/ln(2) */
+
+float expf(const float x)
+{
+ int n;
+ float xn, g, r, z, y;
+ int sign;
+
+ if(x>=0.0)
+ { y=x; sign=0; }
+ else
+ { y=-x; sign=1; }
+
+ if(y<EXPEPS) return 1.0;
+
+ if(y>BIGX)
+ {
+ if(sign)
+ {
+ errno=ERANGE;
+ return XMAX;
+ }
+ else
+ {
+ return 0.0;
+ }
+ }
+
+ z=y*K1;
+ n=z;
+
+ if(n<0) --n;
+ if(z-n>=0.5) ++n;
+ xn=n;
+ g=((y-xn*C1))-xn*C2;
+ z=g*g;
+ r=P(z)*g;
+ r=0.5+(r/(Q(z)-r));
+
+ n++;
+ z=ldexpf(r, n);
+ if(sign)
+ return 1.0/z;
+ else
+ return z;
+}
+
--- /dev/null
+/* fabsf.c: Returns the absolute value of a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float fabsf(const float x) reentrant
+{
+ union float_long fl;
+
+ fl.f = x;
+ fl.l &= 0x7fffffff;
+ return fl.f;
+}
--- /dev/null
+/* floorf.c: Returns the integer smaller or equal than x
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float floorf (float x) reentrant
+{
+ long r;
+ r=x;
+ if (r<=0)
+ return (r+((r>x)?-1:0));
+ else
+ return r;
+}
--- /dev/null
+/* frexpf.c: Returns the exponent and mantisa of a 32 bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float frexpf(const float x, int *pw2)
+{
+ union float_long fl;
+ long int i;
+
+ fl.f=x;
+ /* Find the exponent (power of 2) */
+ i = ( fl.l >> 23) & 0x000000ff;
+ i -= 0x7e;
+ *pw2 = i;
+ fl.l &= 0x807fffff; /* strip all exponent bits */
+ fl.l |= 0x3f000000; /* mantissa between 0.5 and 1 */
+ return(fl.f);
+}
--- /dev/null
+/* ldexpf.c: Build a float from a mantisa and exponent.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float ldexpf(const float x, const int pw2)
+{
+ union float_long fl;
+ long e;
+
+ fl.f = x;
+
+ e=(fl.l >> 23) & 0x000000ff;
+ e+=pw2;
+ fl.l= ((e & 0xff) << 23) | (fl.l & 0x807fffff);
+
+ return(fl.f);
+}
_fs2schar
_fs2sint
_fs2slong
+fabsf
+frexpf
+ldexpf
+expf
+powf
+sincosf
+sinf
+cosf
+logf
+log10f
+sqrtf
+tancotf
+tanf
+cotf
+asincosf
+asinf
+acosf
+atanf
+atan2f
+sincoshf
+sinhf
+coshf
+tanhf
+floorf
+ceilf
+modff
--- /dev/null
+/* log10f.c: Computes the base 10 log of a 32 bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float log10f(const float x) reentrant
+{
+ return logf(x)*0.4342944819;
+}
--- /dev/null
+/* logf.c: Computes the natural log of a 32 bit float as outlined in [1].
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+/*Constans for 24 bits or less (8 decimal digits)*/
+#define A0 -0.5527074855E+0
+#define B0 -0.6632718214E+1
+#define A(w) (A0)
+#define B(w) (w+B0)
+
+#define C0 0.70710678118654752440
+#define C1 0.693359375 /*355.0/512.0*/
+#define C2 -2.121944400546905827679E-4
+
+float logf(const float x) reentrant
+{
+ float Rz, f, z, w, znum, zden, xn;
+ int n;
+
+ if (x<=0.0)
+ {
+ errno=EDOM;
+ return 0.0;
+ }
+ f=frexpf(x, &n);
+ znum=f-0.5;
+ if (f>C0)
+ {
+ znum-=0.5;
+ zden=(f*0.5)+0.5;
+ }
+ else
+ {
+ n--;
+ zden=znum*0.5+0.5;
+ }
+ z=znum/zden;
+ w=z*z;
+
+ Rz=z+z*(w*A(w)/B(w));
+ xn=n;
+ return ((xn*C2+Rz)+xn*C1);
+}
--- /dev/null
+/* modff.c: Returns both the integer and fraction of a float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float modff(float x, float * y)
+{
+ *y=((int)x);
+ return (x-*y);
+}
--- /dev/null
+/* powf.c: Computes x**y where x and y are 32-bit floats.
+ WARNING: less that 6 digits accuracy.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float powf(const float x, const float y)
+{
+ if(y == 0.0) return 1.0;
+ if(y==1.0) return x;
+ if(x <= 0.0) return 0.0;
+ return expf(logf(x) * y);
+}
--- /dev/null
+/* sincosf.c: Computes sin or cos of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define r1 -0.1666665668E+0
+#define r2 0.8333025139E-2
+#define r3 -0.1980741872E-3
+#define r4 0.2601903036E-5
+
+/* PI=C1+C2 */
+#define C1 3.140625
+#define C2 9.676535897E-4
+
+/*A reasonable value for YMAX is the int part of PI*B**(t/2)=3.1416*2**(12)*/
+#define YMAX 12867.0
+
+float sincosf(const float x, const int iscos)
+{
+ float y, f, r, g, XN;
+ int N, sign;
+
+ if(iscos)
+ {
+ y=fabsf(x)+HALF_PI;
+ sign=0;
+ }
+ else
+ {
+ if(x<0.0)
+ { y=-x; sign=1; }
+ else
+ { y=x; sign=0; }
+ }
+
+ if(y>YMAX)
+ {
+ errno=ERANGE;
+ return 0.0;
+ }
+
+ /*Round y/PI to the nearest integer*/
+ N=((y*iPI)+0.5); /*y is positive*/
+
+ /*If N is odd change sign*/
+ if(N&1) sign=~sign;
+
+ XN=N;
+ /*Cosine required? (is done here to keep accuracy)*/
+ if(iscos) XN-=0.5;
+
+ y=fabsf(x);
+ r=(int)y;
+ g=y-r;
+ f=((r-XN*C1)+g)-XN*C2;
+
+ g=f*f;
+ if(g>EPS2) //Used to be if(fabsf(f)>EPS)
+ {
+ r=(((r4*g+r3)*g+r2)*g+r1)*g;
+ f+=f*r;
+ }
+ return (sign?-f:f);
+}
--- /dev/null
+/* sincoshf.c: Computes sinh or cosh of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P0 -0.713793159E+1
+#define P1 -0.190333999E+0
+#define Q0 -0.428277109E+2
+#define Q1 0.100000000E+1
+
+#define P(z) (P1*z+P0)
+#define Q(z) (Q1*z+Q0)
+
+#define K1 0.69316101074218750000E+0 /* ln(v) */
+#define K2 0.24999308500451499336E+0 /* v**(-2) */
+#define K3 0.13830277879601902638E-4 /* v/2-1 */
+
+//WMAX is defined as ln(XMAX)-ln(v)+0.69
+#define WMAX 44.93535952E+0
+//WBAR 0.35*(b+1)
+#define WBAR 1.05
+#define YBAR 9.0 /*Works for me*/
+
+float sincoshf(const float x, const int iscosh)
+{
+ float y, w, z;
+ int sign;
+
+ if (x<0.0) { y=-x; sign=1; }
+ else { y=x; sign=0; }
+
+ if ((y>1.0) || iscosh)
+ {
+ if(y>YBAR)
+ {
+ w=y-K1;
+ if (w>WMAX)
+ {
+ errno=ERANGE;
+ z=XMAX;
+ }
+ else
+ {
+ z=expf(w);
+ z+=K3*z;
+ }
+ }
+ else
+ {
+ z=expf(y);
+ w=1.0/z;
+ if(!iscosh) w=-w;
+ z=(z+w)*0.5;
+ }
+ if(sign) z=-z;
+ }
+ else
+ {
+ if (y<EPS)
+ z=x;
+ else
+ {
+ z=x*x;
+ z=x+x*z*P(z)/Q(z);
+ }
+ }
+ return z;
+}
--- /dev/null
+/* sinf.c: Computes sin(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float sincosf(const float x, const int iscos);
+
+float sinf(const float x) reentrant
+{
+ if (x==0.0) return 0.0;
+ return sincosf(x, 0);
+}
--- /dev/null
+/* sinhf.c: Computes sinh(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float sincoshf(const float x, const int iscosh);
+
+float sinhf(const float x) reentrant
+{
+ return sincoshf(x, 0);
+}
--- /dev/null
+/* sqrtf.c: Computes square root of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+float sqrtf(const float x) reentrant
+{
+ float f, y;
+ int n;
+
+ if (x==0.0) return x;
+ else if (x==1.0) return 1.0;
+ else if (x<0.0)
+ {
+ errno=EDOM;
+ return 0.0;
+ }
+ f=frexpf(x, &n);
+ y=0.41731+0.59016*f; /*Educated guess*/
+ /*For a 24 bit mantisa (float), two iterations are sufficient*/
+ y+=f/y;
+ y=ldexpf(y, -2) + f/y; /*Faster version of 0.25 * y + f/y*/
+
+ if (n&1)
+ {
+ y*=0.7071067812;
+ ++n;
+ }
+ return ldexpf(y, n/2);
+}
--- /dev/null
+/* tancotf.c: Computes tan or cot of a 32-bit float as outlined in [1]
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P0 0.100000000E+1
+#define P1 -0.958017723E-1
+#define Q0 0.100000000E+1
+#define Q1 -0.429135777E+0
+#define Q2 0.971685835E-2
+
+#define C1 1.5703125
+#define C2 4.83826794897E-4
+
+#define P(f,g) (P1*g*f+f)
+#define Q(g) ((Q2*g+Q1)*g+Q0)
+
+//A reasonable choice for YMAX is the integer part of B**(t/2)*PI/2:
+#define YMAX 6433.0
+
+float tancotf(const float x, const int iscotan)
+{
+ float f, g, xn, xnum, xden;
+ int n;
+
+ if (fabsf(x) > YMAX)
+ {
+ errno = ERANGE;
+ return 0.0;
+ }
+
+ /*Round x*2*PI to the nearest integer*/
+ n=(x*TWO_O_PI+(x>0.0?0.5:-0.5)); /*works for +-x*/
+ xn=n;
+
+ xnum=(int)x;
+ xden=x-xnum;
+ f=((xnum-xn*C1)+xden)-xn*C2;
+
+ if (fabsf(f) < EPS)
+ {
+ xnum = f;
+ xden = 1.0;
+ }
+ else
+ {
+ g = f*f;
+ xnum = P(f,g);
+ xden = Q(g);
+ }
+
+ if(n&1)
+ //xn is odd
+ {
+ if(iscotan) return (-xnum/xden);
+ else return (-xden/xnum);
+ }
+ else
+ {
+ if(iscotan) return (xden/xnum);
+ else return (xnum/xden);
+ }
+}
+
--- /dev/null
+/* tanf.c: Computes tan(x) where x is a 32-bit float.
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+
+float tancotf(const float x, const int iscot);
+
+float tanf(const float x) reentrant
+{
+ return tancotf(x, 0);
+}
+
--- /dev/null
+/* tanhf.c: Computes tanh(x) where x is a 32-bit float as outlined in [1].
+
+ Copyright (C) 2001, 2002 Jesus Calvino-Fraga, jesusc@ieee.org
+
+ This library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ This library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with this library; if not, write to the Free Software
+ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */
+
+/* [1] William James Cody and W. M. Waite. _Software manual for the
+ elementary functions_, Englewood Cliffs, N.J.:Prentice-Hall, 1980. */
+
+/* Version 1.0 - Initial release */
+
+#include <math.h>
+#include <errno.h>
+
+#define P0 -0.8237728127E+0
+#define P1 -0.3831010665E-2
+#define Q0 0.2471319654E+1
+#define Q1 0.1000000000E+1
+
+/* ln(3)/2 */
+#define K1 0.5493061443E+0
+/* SBIG=[ln(2)+(t+1)*ln(B)]/2 */
+#define SBIG 9.01091
+
+#define P(g) ((P1*g+P0)*g)
+#define Q(g) (Q1*g+Q0)
+
+float tanhf(const float x) reentrant
+{
+ float f, g, r;
+
+ f=fabsf(x);
+ if(f>SBIG) r=1.0;
+ else if(f>K1)
+ {
+ r=0.5-1.0/(expf(f+f)+1.0);
+ r+=r;
+ }
+ else if(f<EPS) r=f;
+ else
+ {
+ g=f*f;
+ r=f+f*(P(g)/Q(g));
+ }
+ if(x<0.0) r=-r;
+ return r;
+}
+