From 227c08c9cb3f6c096b068d52219ff5e5a2c9c533 Mon Sep 17 00:00:00 2001 From: vrokas Date: Mon, 24 Jan 2005 10:38:53 +0000 Subject: [PATCH] * device/lib/pic16/libm: added Math library sources git-svn-id: https://sdcc.svn.sourceforge.net/svnroot/sdcc/trunk/sdcc@3654 4a8a32a2-be11-0410-ad9d-d568d2c75423 --- ChangeLog | 4 ++ device/lib/pic16/libm/Makefile | 109 +++++++++++++++++++++++++++++++ device/lib/pic16/libm/acosf.c | 35 ++++++++++ device/lib/pic16/libm/asincosf.c | 92 ++++++++++++++++++++++++++ device/lib/pic16/libm/asinf.c | 36 ++++++++++ device/lib/pic16/libm/atan2f.c | 49 ++++++++++++++ device/lib/pic16/libm/atanf.c | 80 +++++++++++++++++++++++ device/lib/pic16/libm/ceilf.c | 35 ++++++++++ device/lib/pic16/libm/cosf.c | 33 ++++++++++ device/lib/pic16/libm/coshf.c | 32 +++++++++ device/lib/pic16/libm/cotf.c | 45 +++++++++++++ device/lib/pic16/libm/errno.c | 25 +++++++ device/lib/pic16/libm/expf.c | 90 +++++++++++++++++++++++++ device/lib/pic16/libm/fabsf.c | 35 ++++++++++ device/lib/pic16/libm/floorf.c | 35 ++++++++++ device/lib/pic16/libm/frexpf.c | 41 ++++++++++++ device/lib/pic16/libm/ldexpf.c | 40 ++++++++++++ device/lib/pic16/libm/log10f.c | 31 +++++++++ device/lib/pic16/libm/logf.c | 74 +++++++++++++++++++++ device/lib/pic16/libm/modff.c | 31 +++++++++ device/lib/pic16/libm/powf.c | 35 ++++++++++ device/lib/pic16/libm/sincosf.c | 90 +++++++++++++++++++++++++ device/lib/pic16/libm/sincoshf.c | 93 ++++++++++++++++++++++++++ device/lib/pic16/libm/sinf.c | 33 ++++++++++ device/lib/pic16/libm/sinhf.c | 32 +++++++++ device/lib/pic16/libm/sqrtf.c | 55 ++++++++++++++++ device/lib/pic16/libm/tancotf.c | 89 +++++++++++++++++++++++++ device/lib/pic16/libm/tanf.c | 33 ++++++++++ device/lib/pic16/libm/tanhf.c | 64 ++++++++++++++++++ 29 files changed, 1476 insertions(+) create mode 100644 device/lib/pic16/libm/Makefile create mode 100644 device/lib/pic16/libm/acosf.c create mode 100644 device/lib/pic16/libm/asincosf.c create mode 100644 device/lib/pic16/libm/asinf.c create mode 100644 device/lib/pic16/libm/atan2f.c create mode 100644 device/lib/pic16/libm/atanf.c create mode 100644 device/lib/pic16/libm/ceilf.c create mode 100644 device/lib/pic16/libm/cosf.c create mode 100644 device/lib/pic16/libm/coshf.c create mode 100644 device/lib/pic16/libm/cotf.c create mode 100644 device/lib/pic16/libm/errno.c create mode 100644 device/lib/pic16/libm/expf.c create mode 100644 device/lib/pic16/libm/fabsf.c create mode 100644 device/lib/pic16/libm/floorf.c create mode 100644 device/lib/pic16/libm/frexpf.c create mode 100644 device/lib/pic16/libm/ldexpf.c create mode 100644 device/lib/pic16/libm/log10f.c create mode 100644 device/lib/pic16/libm/logf.c create mode 100644 device/lib/pic16/libm/modff.c create mode 100644 device/lib/pic16/libm/powf.c create mode 100644 device/lib/pic16/libm/sincosf.c create mode 100644 device/lib/pic16/libm/sincoshf.c create mode 100644 device/lib/pic16/libm/sinf.c create mode 100644 device/lib/pic16/libm/sinhf.c create mode 100644 device/lib/pic16/libm/sqrtf.c create mode 100644 device/lib/pic16/libm/tancotf.c create mode 100644 device/lib/pic16/libm/tanf.c create mode 100644 device/lib/pic16/libm/tanhf.c diff --git a/ChangeLog b/ChangeLog index 8c365448..97560ae2 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,7 @@ +2005-01-24 Vangelis Rokas + + * device/lib/pic16/libm: added Math library sources + 2005-01-24 Raphael Neider * src/pic16/pcode.h: added second memory operand to pCodeOpReg diff --git a/device/lib/pic16/libm/Makefile b/device/lib/pic16/libm/Makefile new file mode 100644 index 00000000..629407cf --- /dev/null +++ b/device/lib/pic16/libm/Makefile @@ -0,0 +1,109 @@ +# +# 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 +# +# 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 + diff --git a/device/lib/pic16/libm/acosf.c b/device/lib/pic16/libm/acosf.c new file mode 100644 index 00000000..134c5fa3 --- /dev/null +++ b/device/lib/pic16/libm/acosf.c @@ -0,0 +1,35 @@ +/* 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 + +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); +} diff --git a/device/lib/pic16/libm/asincosf.c b/device/lib/pic16/libm/asincosf.c new file mode 100644 index 00000000..f4bed6d5 --- /dev/null +++ b/device/lib/pic16/libm/asincosf.c @@ -0,0 +1,92 @@ +/* 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 +#include + +#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; +} diff --git a/device/lib/pic16/libm/asinf.c b/device/lib/pic16/libm/asinf.c new file mode 100644 index 00000000..c2928221 --- /dev/null +++ b/device/lib/pic16/libm/asinf.c @@ -0,0 +1,36 @@ +/* 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 + +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); +} + diff --git a/device/lib/pic16/libm/atan2f.c b/device/lib/pic16/libm/atan2f.c new file mode 100644 index 00000000..f5a87483 --- /dev/null +++ b/device/lib/pic16/libm/atan2f.c @@ -0,0 +1,49 @@ +/* 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 +#include + +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; +} diff --git a/device/lib/pic16/libm/atanf.c b/device/lib/pic16/libm/atanf.c new file mode 100644 index 00000000..d143fda6 --- /dev/null +++ b/device/lib/pic16/libm/atanf.c @@ -0,0 +1,80 @@ +/* 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 +#include + +#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)1) r=-r; + r+=a[n]; + if(x<0.0) r=-r; + return r; +} + diff --git a/device/lib/pic16/libm/ceilf.c b/device/lib/pic16/libm/ceilf.c new file mode 100644 index 00000000..a6343385 --- /dev/null +++ b/device/lib/pic16/libm/ceilf.c @@ -0,0 +1,35 @@ +/* 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 + +float ceilf(float x) _MATH_REENTRANT +{ + long r; + r=x; + if (r<0) + return r; + else + return (r+((r + +float sincosf(float x, int iscos); + +float cosf(float x) _MATH_REENTRANT +{ + if (x==0.0) return 1.0; + return sincosf(x, 1); +} diff --git a/device/lib/pic16/libm/coshf.c b/device/lib/pic16/libm/coshf.c new file mode 100644 index 00000000..4d1d9738 --- /dev/null +++ b/device/lib/pic16/libm/coshf.c @@ -0,0 +1,32 @@ +/* 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 + +float sincoshf(const float x, const int iscosh); + +float coshf(const float x) _MATH_REENTRANT +{ + return sincoshf(x, 1); +} diff --git a/device/lib/pic16/libm/cotf.c b/device/lib/pic16/libm/cotf.c new file mode 100644 index 00000000..cd52c68e --- /dev/null +++ b/device/lib/pic16/libm/cotf.c @@ -0,0 +1,45 @@ +/* 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 +#include + +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); +} + diff --git a/device/lib/pic16/libm/errno.c b/device/lib/pic16/libm/errno.c new file mode 100644 index 00000000..031d2bef --- /dev/null +++ b/device/lib/pic16/libm/errno.c @@ -0,0 +1,25 @@ +/*------------------------------------------------------------------------- + + 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; diff --git a/device/lib/pic16/libm/expf.c b/device/lib/pic16/libm/expf.c new file mode 100644 index 00000000..0a59823b --- /dev/null +++ b/device/lib/pic16/libm/expf.c @@ -0,0 +1,90 @@ +/* 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 +#include + +#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(yBIGX) + { + 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; +} + diff --git a/device/lib/pic16/libm/fabsf.c b/device/lib/pic16/libm/fabsf.c new file mode 100644 index 00000000..2b9fe84f --- /dev/null +++ b/device/lib/pic16/libm/fabsf.c @@ -0,0 +1,35 @@ +/* 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 +#include + +float fabsf(const float x) _MATH_REENTRANT +{ + union float_long fl; + + fl.f = x; + fl.l &= 0x7fffffff; + return fl.f; +} diff --git a/device/lib/pic16/libm/floorf.c b/device/lib/pic16/libm/floorf.c new file mode 100644 index 00000000..a9eea702 --- /dev/null +++ b/device/lib/pic16/libm/floorf.c @@ -0,0 +1,35 @@ +/* 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 + +float floorf (float x) _MATH_REENTRANT +{ + long r; + r=x; + if (r<=0) + return (r+((r>x)?-1:0)); + else + return r; +} diff --git a/device/lib/pic16/libm/frexpf.c b/device/lib/pic16/libm/frexpf.c new file mode 100644 index 00000000..40f95c34 --- /dev/null +++ b/device/lib/pic16/libm/frexpf.c @@ -0,0 +1,41 @@ +/* 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 +#include + +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); +} diff --git a/device/lib/pic16/libm/ldexpf.c b/device/lib/pic16/libm/ldexpf.c new file mode 100644 index 00000000..eeb5a059 --- /dev/null +++ b/device/lib/pic16/libm/ldexpf.c @@ -0,0 +1,40 @@ +/* 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 +#include + +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); +} diff --git a/device/lib/pic16/libm/log10f.c b/device/lib/pic16/libm/log10f.c new file mode 100644 index 00000000..7c18d7b6 --- /dev/null +++ b/device/lib/pic16/libm/log10f.c @@ -0,0 +1,31 @@ +/* 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 +#include + +float log10f(const float x) _MATH_REENTRANT +{ + return logf(x)*0.4342944819; +} diff --git a/device/lib/pic16/libm/logf.c b/device/lib/pic16/libm/logf.c new file mode 100644 index 00000000..aee57811 --- /dev/null +++ b/device/lib/pic16/libm/logf.c @@ -0,0 +1,74 @@ +/* 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 +#include + +/*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); +} diff --git a/device/lib/pic16/libm/modff.c b/device/lib/pic16/libm/modff.c new file mode 100644 index 00000000..544062d1 --- /dev/null +++ b/device/lib/pic16/libm/modff.c @@ -0,0 +1,31 @@ +/* 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 + +float modff(float x, float * y) +{ + *y=((int)x); + return (x-*y); +} diff --git a/device/lib/pic16/libm/powf.c b/device/lib/pic16/libm/powf.c new file mode 100644 index 00000000..77f768f5 --- /dev/null +++ b/device/lib/pic16/libm/powf.c @@ -0,0 +1,35 @@ +/* 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 +#include + +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); +} diff --git a/device/lib/pic16/libm/sincosf.c b/device/lib/pic16/libm/sincosf.c new file mode 100644 index 00000000..965de6bd --- /dev/null +++ b/device/lib/pic16/libm/sincosf.c @@ -0,0 +1,90 @@ +/* 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 +#include + +#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); +} diff --git a/device/lib/pic16/libm/sincoshf.c b/device/lib/pic16/libm/sincoshf.c new file mode 100644 index 00000000..87de5d13 --- /dev/null +++ b/device/lib/pic16/libm/sincoshf.c @@ -0,0 +1,93 @@ +/* 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 +#include + +#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 + +float sincosf(float x, int iscos); + +float sinf(float x) _MATH_REENTRANT +{ + if (x==0.0) return 0.0; + return sincosf(x, 0); +} diff --git a/device/lib/pic16/libm/sinhf.c b/device/lib/pic16/libm/sinhf.c new file mode 100644 index 00000000..d3386a5c --- /dev/null +++ b/device/lib/pic16/libm/sinhf.c @@ -0,0 +1,32 @@ +/* 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 + +float sincoshf(const float x, const int iscosh); + +float sinhf(const float x) _MATH_REENTRANT +{ + return sincoshf(x, 0); +} diff --git a/device/lib/pic16/libm/sqrtf.c b/device/lib/pic16/libm/sqrtf.c new file mode 100644 index 00000000..ff39339f --- /dev/null +++ b/device/lib/pic16/libm/sqrtf.c @@ -0,0 +1,55 @@ +/* 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 +#include + +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); +} diff --git a/device/lib/pic16/libm/tancotf.c b/device/lib/pic16/libm/tancotf.c new file mode 100644 index 00000000..cf529099 --- /dev/null +++ b/device/lib/pic16/libm/tancotf.c @@ -0,0 +1,89 @@ +/* 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 +#include + +#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); + } +} + diff --git a/device/lib/pic16/libm/tanf.c b/device/lib/pic16/libm/tanf.c new file mode 100644 index 00000000..bd247ef4 --- /dev/null +++ b/device/lib/pic16/libm/tanf.c @@ -0,0 +1,33 @@ +/* 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 + +float tancotf(const float x, const int iscot); + +float tanf(const float x) _MATH_REENTRANT +{ + return tancotf(x, 0); +} + diff --git a/device/lib/pic16/libm/tanhf.c b/device/lib/pic16/libm/tanhf.c new file mode 100644 index 00000000..27488d94 --- /dev/null +++ b/device/lib/pic16/libm/tanhf.c @@ -0,0 +1,64 @@ +/* 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 +#include + +#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