diff --git a/src/system/libroot/posix/glibc/arch/arm/Jamfile b/src/system/libroot/posix/glibc/arch/arm/Jamfile index b4c2f6dad5..deb19a0bb5 100644 --- a/src/system/libroot/posix/glibc/arch/arm/Jamfile +++ b/src/system/libroot/posix/glibc/arch/arm/Jamfile @@ -33,14 +33,41 @@ local genericSources = s_isinf.c s_isinff.c s_isnan.c s_isnanf.c s_signbit.c s_signbitf.c s_signbitl.c + + s_floor.c + s_ceil.c s_ceilf.c + s_modf.c + w_pow.c e_pow.c slowpow.c + e_exp.c slowexp.c + dosincos.c + doasin.c + sincos32.c + branred.c + halfulp.c + mpa.c mplog.c mpexp.c + s_sin.c + s_atan.c + s_tan.c + e_asin.c w_asin.c + e_log10.c w_log10.c + e_acos.c w_acos.c + e_atan2.c w_atan2.c mpatan2.c mpatan.c mptan.c mpsqrt.c + e_fmod.c w_fmod.c + e_log.c w_log.c + s_ldexp.c s_ldexpf.c ; MergeObject posix_gnu_arch_$(TARGET_ARCH)_generic.o : $(genericSources) ; +MergeObject posix_gnu_arch_$(TARGET_ARCH)_others.o : + e_sqrt.c +; + MergeObjectFromObjects posix_gnu_arch_$(TARGET_ARCH).o : : posix_gnu_arch_$(TARGET_ARCH)_generic.o + posix_gnu_arch_$(TARGET_ARCH)_others.o ; SEARCH on [ FGristFiles $(genericSources) ] diff --git a/src/system/libroot/posix/glibc/arch/arm/e_log.c b/src/system/libroot/posix/glibc/arch/arm/e_log.c new file mode 100644 index 0000000000..1a9967b546 --- /dev/null +++ b/src/system/libroot/posix/glibc/arch/arm/e_log.c @@ -0,0 +1,203 @@ +/* + * IBM Accurate Mathematical Library + * written by International Business Machines Corp. + * Copyright (C) 2001 Free Software Foundation + * + * This program 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 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 Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser 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. + */ +/*********************************************************************/ +/* */ +/* MODULE_NAME:ulog.c */ +/* */ +/* FUNCTION:ulog */ +/* */ +/* FILES NEEDED: dla.h endian.h mpa.h mydefs.h ulog.h */ +/* mpexp.c mplog.c mpa.c */ +/* ulog.tbl */ +/* */ +/* An ultimate log routine. Given an IEEE double machine number x */ +/* it computes the correctly rounded (to nearest) value of log(x). */ +/* Assumption: Machine arithmetic operations are performed in */ +/* round to nearest mode of IEEE 754 standard. */ +/* */ +/*********************************************************************/ + + +#include "endian.h" +#include "dla.h" +#include "mpa.h" +#include "MathLib.h" +#include "math_private.h" + +void __mplog(mp_no *, mp_no *, int); + +/*********************************************************************/ +/* An ultimate log routine. Given an IEEE double machine number x */ +/* it computes the correctly rounded (to nearest) value of log(x). */ +/*********************************************************************/ +double __ieee754_log(double x) { +#define M 4 + static const int pr[M]={8,10,18,32}; + int i,j,n,ux,dx,p; +#if 0 + int k; +#endif + double dbl_n,u,p0,q,r0,w,nln2a,luai,lubi,lvaj,lvbj, + sij,ssij,ttij,A,B,B0,y,y1,y2,polI,polII,sa,sb, + t1,t2,t3,t4,t5,t6,t7,t8,t,ra,rb,ww, + a0,aa0,s1,s2,ss2,s3,ss3,a1,aa1,a,aa,b,bb,c; + number num; + mp_no mpx,mpy,mpy1,mpy2,mperr; + +#include "ulog.tbl" +#include "ulog.h" + + /* Treating special values of x ( x<=0, x=INF, x=NaN etc.). */ + + num.d = x; ux = num.i[HIGH_HALF]; dx = num.i[LOW_HALF]; + n=0; + if (ux < 0x00100000) { + if (((ux & 0x7fffffff) | dx) == 0) return MHALF/ZERO; /* return -INF */ + if (ux < 0) return (x-x)/ZERO; /* return NaN */ + n -= 54; x *= two54.d; /* scale x */ + num.d = x; + } + if (ux >= 0x7ff00000) return x+x; /* INF or NaN */ + + /* Regular values of x */ + + w = x-ONE; + if (ABS(w) > U03) { goto case_03; } + + + /*--- Stage I, the case abs(x-1) < 0.03 */ + + t8 = MHALF*w; + EMULV(t8,w,a,aa,t1,t2,t3,t4,t5) + EADD(w,a,b,bb) + + /* Evaluate polynomial II */ + polII = (b0.d+w*(b1.d+w*(b2.d+w*(b3.d+w*(b4.d+ + w*(b5.d+w*(b6.d+w*(b7.d+w*b8.d))))))))*w*w*w; + c = (aa+bb)+polII; + + /* End stage I, case abs(x-1) < 0.03 */ + if ((y=b+(c+b*E2)) == b+(c-b*E2)) return y; + + /*--- Stage II, the case abs(x-1) < 0.03 */ + + a = d11.d+w*(d12.d+w*(d13.d+w*(d14.d+w*(d15.d+w*(d16.d+ + w*(d17.d+w*(d18.d+w*(d19.d+w*d20.d)))))))); + EMULV(w,a,s2,ss2,t1,t2,t3,t4,t5) + ADD2(d10.d,dd10.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d9.d,dd9.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d8.d,dd8.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d7.d,dd7.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d6.d,dd6.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d5.d,dd5.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d4.d,dd4.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d3.d,dd3.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(d2.d,dd2.d,s2,ss2,s3,ss3,t1,t2) + MUL2(w,ZERO,s3,ss3,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + MUL2(w,ZERO,s2,ss2,s3,ss3,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(w,ZERO, s3,ss3, b, bb,t1,t2) + + /* End stage II, case abs(x-1) < 0.03 */ + if ((y=b+(bb+b*E4)) == b+(bb-b*E4)) return y; + goto stage_n; + + /*--- Stage I, the case abs(x-1) > 0.03 */ + case_03: + + /* Find n,u such that x = u*2**n, 1/sqrt(2) < u < sqrt(2) */ + n += (num.i[HIGH_HALF] >> 20) - 1023; + num.i[HIGH_HALF] = (num.i[HIGH_HALF] & 0x000fffff) | 0x3ff00000; + if (num.d > SQRT_2) { num.d *= HALF; n++; } + u = num.d; dbl_n = (double) n; + + /* Find i such that ui=1+(i-75)/2**8 is closest to u (i= 0,1,2,...,181) */ + num.d += h1.d; + i = (num.i[HIGH_HALF] & 0x000fffff) >> 12; + + /* Find j such that vj=1+(j-180)/2**16 is closest to v=u/ui (j= 0,...,361) */ + num.d = u*Iu[i].d + h2.d; + j = (num.i[HIGH_HALF] & 0x000fffff) >> 4; + + /* Compute w=(u-ui*vj)/(ui*vj) */ + p0=(ONE+(i-75)*DEL_U)*(ONE+(j-180)*DEL_V); + q=u-p0; r0=Iu[i].d*Iv[j].d; w=q*r0; + + /* Evaluate polynomial I */ + polI = w+(a2.d+a3.d*w)*w*w; + + /* Add up everything */ + nln2a = dbl_n*LN2A; + luai = Lu[i][0].d; lubi = Lu[i][1].d; + lvaj = Lv[j][0].d; lvbj = Lv[j][1].d; + EADD(luai,lvaj,sij,ssij) + EADD(nln2a,sij,A ,ttij) + B0 = (((lubi+lvbj)+ssij)+ttij)+dbl_n*LN2B; + B = polI+B0; + + /* End stage I, case abs(x-1) >= 0.03 */ + if ((y=A+(B+E1)) == A+(B-E1)) return y; + + + /*--- Stage II, the case abs(x-1) > 0.03 */ + + /* Improve the accuracy of r0 */ + EMULV(p0,r0,sa,sb,t1,t2,t3,t4,t5) + t=r0*((ONE-sa)-sb); + EADD(r0,t,ra,rb) + + /* Compute w */ + MUL2(q,ZERO,ra,rb,w,ww,t1,t2,t3,t4,t5,t6,t7,t8) + + EADD(A,B0,a0,aa0) + + /* Evaluate polynomial III */ + s1 = (c3.d+(c4.d+c5.d*w)*w)*w; + EADD(c2.d,s1,s2,ss2) + MUL2(s2,ss2,w,ww,s3,ss3,t1,t2,t3,t4,t5,t6,t7,t8) + MUL2(s3,ss3,w,ww,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8) + ADD2(s2,ss2,w,ww,s3,ss3,t1,t2) + ADD2(s3,ss3,a0,aa0,a1,aa1,t1,t2) + + /* End stage II, case abs(x-1) >= 0.03 */ + if ((y=a1+(aa1+E3)) == a1+(aa1-E3)) return y; + + + /* Final stages. Use multi-precision arithmetic. */ + stage_n: + + for (i=0; i>14]; + s=a.x; + /*----------------- 2^-1022 <= | x |< 2^1024 -----------------*/ + if (k>0x000fffff && k<0x7ff00000) { + y=1.0-t*(t*s); + t=t*(rt0+y*(rt1+y*(rt2+y*rt3))); + c.i[HIGH_HALF]=0x20000000+((k&0x7fe00000)>>1); + y=t*s; + hy=(y+big)-big; + del=0.5*t*((s-hy*hy)-(y-hy)*(y+hy)); + res=y+del; + if (res == (res+1.002*((y-res)+del))) return res*c.x; + else { + res1=res+1.5*((y-res)+del); + EMULV(res,res1,z,zz,p,hx,tx,hy,ty); /* (z+zz)=res*res1 */ + return ((((z-s)+zz)<0)?max(res,res1):min(res,res1))*c.x; + } + } + else { + if ((k & 0x7ff00000) == 0x7ff00000) + return x*x+x; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf, sqrt(-inf)=sNaN */ + if (x==0) return x; /* sqrt(+0)=+0, sqrt(-0)=-0 */ + if (k<0) return (x-x)/(x-x); /* sqrt(-ve)=sNaN */ + return tm256.x*__ieee754_sqrt(x*t512.x); + } +} diff --git a/src/system/libroot/posix/glibc/arch/arm/uroot.h b/src/system/libroot/posix/glibc/arch/arm/uroot.h new file mode 100644 index 0000000000..02b74cb262 --- /dev/null +++ b/src/system/libroot/posix/glibc/arch/arm/uroot.h @@ -0,0 +1,44 @@ +/* + * IBM Accurate Mathematical Library + * Written by International Business Machines Corp. + * Copyright (C) 2001 Free Software Foundation, Inc. + * + * This program 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 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 Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser 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. + */ + +/******************************************************************/ +/* */ +/* MODULE_NAME:uroot.h */ +/* */ +/* common data and variables prototype and definition */ +/******************************************************************/ + +#ifndef UROOT_H +#define UROOT_H + +#ifdef BIG_ENDI + static const mynumber +/**/ t512 = {{0x5ff00000, 0x00000000 }}, /* 2^512 */ +/**/ tm256 = {{0x2ff00000, 0x00000000 }}; /* 2^-256 */ + +#else +#ifdef LITTLE_ENDI + static const mynumber +/**/ t512 = {{0x00000000, 0x5ff00000 }}, /* 2^512 */ +/**/ tm256 = {{0x00000000, 0x2ff00000 }}; /* 2^-256 */ +#endif +#endif + +#endif diff --git a/src/system/libroot/posix/glibc/include/arch/arm/bits/fenv.h b/src/system/libroot/posix/glibc/include/arch/arm/bits/fenv.h new file mode 100644 index 0000000000..7bd242385f --- /dev/null +++ b/src/system/libroot/posix/glibc/include/arch/arm/bits/fenv.h @@ -0,0 +1,58 @@ +/* Copyright (C) 1997, 1998, 1999 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C 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. + + The GNU C 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 the GNU C Library; if not, write to the Free + Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA + 02111-1307 USA. */ + +#ifndef _FENV_H +# error "Never use directly; include instead." +#endif + +/* Define bits representing exceptions in the FPU status word. */ +enum + { + FE_INVALID = 1, +#define FE_INVALID FE_INVALID + FE_DIVBYZERO = 2, +#define FE_DIVBYZERO FE_DIVBYZERO + FE_OVERFLOW = 4, +#define FE_OVERFLOW FE_OVERFLOW + FE_UNDERFLOW = 8, +#define FE_UNDERFLOW FE_UNDERFLOW + }; + +/* Amount to shift by to convert an exception to a mask bit. */ +#define FE_EXCEPT_SHIFT 16 + +/* All supported exceptions. */ +#define FE_ALL_EXCEPT \ + (FE_INVALID | FE_DIVBYZERO | FE_OVERFLOW | FE_UNDERFLOW) + +/* The ARM FPU basically only supports round-to-nearest. Other rounding + modes exist, but you have to encode them in the actual instruction. */ +#define FE_TONEAREST 0 + +/* Type representing exception flags. */ +typedef unsigned long int fexcept_t; + +/* Type representing floating-point environment. */ +typedef struct + { + unsigned long int __cw; + } +fenv_t; + +/* If the default argument is used we use this value. */ +#define FE_DFL_ENV ((fenv_t *) -1l) diff --git a/src/system/libroot/posix/glibc/include/arch/arm/bits/fenvinline.h b/src/system/libroot/posix/glibc/include/arch/arm/bits/fenvinline.h new file mode 100644 index 0000000000..e69de29bb2 diff --git a/src/system/libroot/posix/glibc/include/arch/arm/bits/mathdef.h b/src/system/libroot/posix/glibc/include/arch/arm/bits/mathdef.h index daca1102fc..5d9a3349cf 100644 --- a/src/system/libroot/posix/glibc/include/arch/arm/bits/mathdef.h +++ b/src/system/libroot/posix/glibc/include/arch/arm/bits/mathdef.h @@ -30,6 +30,12 @@ typedef float float_t; /* `float' expressions are evaluated as typedef double double_t; /* `double' expressions are evaluated as `double'. */ +/* Signal that both types are `long double'. */ +# define FLT_EVAL_METHOD 2 + +/* Define `INFINITY' as value of type `float'. */ +# define INFINITY HUGE_VALF + /* The values returned by `ilogb' for 0 and NaN respectively. */ # define FP_ILOGB0 (-2147483647) # define FP_ILOGBNAN (2147483647) @@ -41,3 +47,6 @@ typedef double double_t; /* `double' expressions are evaluated as declaration of all the `long double' function variants. */ # define __NO_LONG_DOUBLE_MATH 1 #endif + +/* Number of decimal digits for the `long double' type. */ +# define DECIMAL_DIG 21