+2005-01-24 Vangelis Rokas <vrokas AT otenet.gr>
+
+ * device/lib/pic16/libm: added Math library sources
+
2005-01-24 Raphael Neider <rneider AT web.de>
* src/pic16/pcode.h: added second memory operand to pCodeOpReg
--- /dev/null
+#
+# Makefile - Makefile to build pic16 Math Library
+#
+# This file is part of the GNU PIC Library.
+#
+# January, 2004
+# The GNU PIC Library is maintained by,
+# Vangelis Rokas <vrokas@otenet.gr>
+#
+# This program is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Library General Public License
+# as published by the Free Software Foundation; either version 2
+# 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 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, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
+#
+#
+# $Id$
+#
+#
+
+include ../Makefile.common
+
+PRJDIR = ../../../..
+
+LIB = libm18f.lib
+
+SRCS = acosf \
+ asincosf \
+ asinf \
+ atan2f \
+ atanf \
+ ceilf \
+ cosf \
+ coshf \
+ cotf \
+ errno \
+ expf \
+ fabsf \
+ floorf \
+ frexpf \
+ ldexpf \
+ log10f \
+ logf \
+ modff \
+ powf \
+ sincosf \
+ sincoshf \
+ sinf \
+ sinhf \
+ sqrtf \
+ tancotf \
+ tanf \
+ tanhf
+
+COMPILE_FLAGS += $(MODELFLAGS) $(OPT_FLAGS)
+
+#CFLAGS += -I$(LIBC_INC_DIR)
+
+
+CFILES = $(patsubst %,%.c,$(SRCS))
+OFILES = $(patsubst %.c,%.o,$(CFILES))
+
+%.o: %.c
+ $(CC) $(CFLAGS) $(COMPILE_FLAGS) -c $<
+
+all: build-library
+
+build-library: $(LIB)
+
+$(LIB): $(OFILES)
+ @echo Creating $(LIB) ...
+ @for object in $(OFILES) ; do \
+ if [ ! -e $(LIB) ]; then \
+ $(AR) -c $(LIB) $$object ; \
+ else \
+ $(AR) -r $(LIB) $$object ; \
+ fi; \
+ echo adding $$object ; \
+ done ;
+ mv -v $(LIB) ../bin
+
+
+all-clean: clean
+
+clean-intermediate:
+ @echo Removing intermediate files ...
+ $(RM) -f *.lst *.asm
+
+clean: clean-intermediate
+ $(RM) -f $(LIB) *.o
+
+dep .depend:
+ -rm .depend
+ @for source in $(CFILES); do \
+ $(CC) $(MM) $(CFLAGS) $$source > .tmpdepend ; \
+ $(SED) s/.rel/.o/g .tmpdepend >> .depend ; \
+ $(RM) -f .tmpdepend; \
+ done
+
+include .depend
+
--- /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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float asincosf(const float x, const int isacos);
+
+float acosf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+
+#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)
+
+#ifdef SDCC_mcs51
+ #define myconst code
+#else
+ #define myconst const
+#endif
+
+float asincosf(const float x, const int isacos)
+{
+ float y, g, r;
+ int i;
+
+ static myconst float a[2]={ 0.0, QUART_PI };
+ static myconst 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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float asincosf(const float x, const int isacos);
+
+float asinf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#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 */
+
+/*
+** $Id$
+*/
+
+#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) */
+
+#ifdef SDCC_mcs51
+ #define myconst code
+#else
+ #define myconst const
+#endif
+
+float atanf(const float x) _MATH_REENTRANT
+{
+ float f, r, g;
+ int n=0;
+ static myconst 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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float ceilf(float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float sincosf(float x, int iscos);
+
+float cosf(float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float sincoshf(const float x, const int iscosh);
+
+float coshf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+#include <errno.h>
+
+float tancotf(const float x, const int iscot);
+
+float cotf(const float x) _MATH_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
+/*-------------------------------------------------------------------------
+
+ errno.c :- just declares errno as a variable
+
+ This library is free software; you can redistribute it and/or modify it
+ under the terms of the GNU Library General Public License as published by the
+ Free Software Foundation; either version 2, 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 Library General Public License for more details.
+
+ You should have received a copy of the GNU Library 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.
+
+-------------------------------------------------------------------------*/
+
+/*
+** $Id$
+*/
+
+int errno;
--- /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 */
+
+/*
+** $Id$
+*/
+
+#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;
+ char 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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+#include <errno.h>
+
+float fabsf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float floorf (float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#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 */
+
+/*
+** $Id$
+*/
+
+#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);
+}
--- /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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+#include <errno.h>
+
+float log10f(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#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) _MATH_REENTRANT
+{
+#if defined(SDCC_mcs51) && defined(SDCC_MODEL_SMALL) \
+ && !defined(SDCC_NOOVERLAY)
+ volatile
+#endif
+ float Rz;
+ float 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 */
+
+/*
+** $Id$
+*/
+
+#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 */
+
+/*
+** $Id$
+*/
+
+#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 */
+
+/*
+** $Id$
+*/
+
+#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(float x, int iscos)
+{
+ float y, f, r, g, XN;
+ int N;
+ char 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 */
+
+/*
+** $Id$
+*/
+
+#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;
+ char 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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float sincosf(float x, int iscos);
+
+float sinf(float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float sincoshf(const float x, const int iscosh);
+
+float sinhf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+#include <errno.h>
+
+float sqrtf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#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 */
+
+/*
+** $Id$
+*/
+
+#include <math.h>
+
+float tancotf(const float x, const int iscot);
+
+float tanf(const float x) _MATH_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 */
+
+/*
+** $Id$
+*/
+
+#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) _MATH_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;
+}
+
projects / fw / sdcc / commitdiff