diff options
author | Jakub Jelinek <jakub@redhat.com> | 2006-11-09 20:12:02 +0000 |
---|---|---|
committer | Jakub Jelinek <jakub@redhat.com> | 2006-11-09 20:12:02 +0000 |
commit | bff4b406f7f47bb2aca83d90d8d886bfbc8ba447 (patch) | |
tree | 7c37b1de14721cc7466c0c3cfa631ff50b2d991a /powerpc-cpu | |
parent | dddf300c357a38a1d26a97ade2e905a0f85befe3 (diff) | |
download | glibc-bff4b406f7f47bb2aca83d90d8d886bfbc8ba447.tar glibc-bff4b406f7f47bb2aca83d90d8d886bfbc8ba447.tar.gz glibc-bff4b406f7f47bb2aca83d90d8d886bfbc8ba447.tar.bz2 glibc-bff4b406f7f47bb2aca83d90d8d886bfbc8ba447.zip |
powerpc-cpu add-on update to v0.4.
Diffstat (limited to 'powerpc-cpu')
44 files changed, 3756 insertions, 0 deletions
diff --git a/powerpc-cpu/ChangeLog b/powerpc-cpu/ChangeLog index a60cf322fd..79f9dc6b1f 100644 --- a/powerpc-cpu/ChangeLog +++ b/powerpc-cpu/ChangeLog @@ -1,3 +1,86 @@ +2006-10-20 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/power4/fpu/slowpow.c: New file. + * sysdeps/powerpc/powerpc64/power4/fpu/slowpow.c: New file. + +2006-10-03 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llround.S: New file. + * sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llroundf.S: New file. + + * sysdeps/powerpc/powerpc32/powerpc64/Makefile: Moved. + * sysdeps/powerpc/powerpc32/powerpc64/memcopy.h: Likewise. + * sysdeps/powerpc/powerpc32/powerpc64/wordcopy.c: Likewise. + * sysdeps/powerpc/powerpc32/power4/Makefile: To here. + * sysdeps/powerpc/powerpc32/power4/memcopy.h: Likewise. + * sysdeps/powerpc/powerpc32/power4/wordcopy.c: Likewise. + + * sysdeps/powerpc/powerpc32/powerpc64/fpu/Makefile: Moved. + * sysdeps/powerpc/powerpc32/powerpc64/fpu/mpa.c: Likewise. + * sysdeps/powerpc/powerpc32/power4/fpu/Makefile: To here. + * sysdeps/powerpc/powerpc32/power4/fpu/mpa.c: Likewise. + +2006-09-29 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/power6x/fpu/s_lrint.S: New file. + * sysdeps/powerpc/powerpc32/power6x/fpu/s_lround.S: New file. + * sysdeps/powerpc/powerpc64/power6x/fpu/s_llrint.S: New file. + * sysdeps/powerpc/powerpc64/power6x/fpu/s_llround.S: New file. + +2006-09-28 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/power5+/fpu/s_llround.S: New file. + * sysdeps/powerpc/powerpc32/power5+/fpu/s_llroundf.S: New file. + * sysdeps/powerpc/powerpc32/power5+/fpu/s_lround.S: New file. + * sysdeps/powerpc/powerpc32/power6x/Implies: New file. + * sysdeps/powerpc/powerpc32/power6x/fpu/Implies: New file. + * sysdeps/powerpc/powerpc64/power5+/fpu/s_llround.S: New file. + * sysdeps/powerpc/powerpc64/power6x/Implies: New file. + * sysdeps/powerpc/powerpc64/power6x/fpu/Implies: New file. + * sysdeps/unix/sysv/linux/powerpc/powerpc32/power6x/fpu/Implies: + New file. + * sysdeps/unix/sysv/linux/powerpc/powerpc64/power6x/fpu/Implies: + New file. + +2006-09-10 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/power6/memcpy.S: New file. + +2006-08-31 Steven Munroe <sjmunroe@us.ibm.com> + + * sysdeps/powerpc/powerpc32/970/Implies: Add + powerpc/powerpc32/powerpc64 to pick up memcopy.h and wordcopy.c + * sysdeps/powerpc/powerpc32/power4/Implies: New file. + * sysdeps/powerpc/powerpc32/power5/Implies: Add + powerpc/powerpc32/powerpc64 to pick up memcopy.h and wordcopy.c + * sysdeps/powerpc/powerpc32/power5+/Implies: Add + powerpc/powerpc32/powerpc64 to pick up memcopy.h and wordcopy.c + * sysdeps/powerpc/powerpc32/power6/Implies: Add + powerpc/powerpc32/powerpc64 to pick up memcopy.h and wordcopy.c + * sysdeps/powerpc/powerpc32/power6/wordcopy.c: New file. + * sysdeps/powerpc/powerpc32/powerpc64/Makefile: New file. + * sysdeps/powerpc/powerpc32/powerpc64/memcopy.h: New file. + * sysdeps/powerpc/powerpc32/powerpc64/wordcopy.c: New file. + * sysdeps/powerpc/powerpc32/powerpc64/fpu/Makefile: New file. + * sysdeps/powerpc/powerpc32/powerpc64/fpu/mpa.c: New file. + * sysdeps/powerpc/powerpc64/970/fpu/Implies: New file. + * sysdeps/powerpc/powerpc64/power4/Makefile: New file. + * sysdeps/powerpc/powerpc64/power4/fpu/Makefile: New file. + * sysdeps/powerpc/powerpc64/power4/fpu/mpa.c: New file. + * sysdeps/powerpc/powerpc64/power4/memcopy.h: New file. + * sysdeps/powerpc/powerpc64/power4/wordcopy.c: New file. + * sysdeps/powerpc/powerpc64/power5+/fpu/Implies: New file. + * sysdeps/powerpc/powerpc64/power6/fpu/Implies: Add + powerpc/powerpc64/power4/fpu to pick up optimized mpa.c. + * sysdeps/powerpc/powerpc64/power6/wordcopy.c: New file. + * sysdeps/unix/sysv/linux/powerpc/powerpc64/970/fpu/Implies: New file. + * sysdeps/unix/sysv/linux/powerpc/powerpc64/malloc-machine.h: New file. + * sysdeps/unix/sysv/linux/powerpc/powerpc64/power4/fpu/Implies: New + file. + * sysdeps/unix/sysv/linux/powerpc/powerpc64/power5/fpu/Implies: New + file. + + 2006-07-06 Steven Munroe <sjmunroe@us.ibm.com> * sysdeps/powerpc/powerpc64/power6/memcpy.S: New file. diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/970/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/970/Implies index 4e3a983426..a0ddab4863 100644 --- a/powerpc-cpu/sysdeps/powerpc/powerpc32/970/Implies +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/970/Implies @@ -1 +1,2 @@ powerpc/powerpc32/power4 +powerpc/powerpc32/powerpc64 diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Implies new file mode 100644 index 0000000000..b2ac1558f5 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Implies @@ -0,0 +1 @@ +powerpc/powerpc32/powerpc64
\ No newline at end of file diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Makefile b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Makefile new file mode 100644 index 0000000000..60aa508ba4 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/Makefile @@ -0,0 +1,6 @@ +# Makefile fragment for POWER4/5/5+. + +ifeq ($(subdir),string) +CFLAGS-wordcopy.c += --param max-variable-expansions-in-unroller=2 --param max-unroll-times=2 -funroll-loops -fpeel-loops -ftree-loop-linear +CFLAGS-memmove.c += --param max-variable-expansions-in-unroller=2 --param max-unroll-times=2 -funroll-loops -fpeel-loops -ftree-loop-linear +endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/Makefile b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/Makefile new file mode 100644 index 0000000000..a6fa75ecbc --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/Makefile @@ -0,0 +1,5 @@ +# Makefile fragment for POWER4/5/5+ with FPU. + +ifeq ($(subdir),math) +CFLAGS-mpa.c += --param max-unroll-times=4 -funroll-loops -fpeel-loops -ftree-loop-linear +endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/mpa.c b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/mpa.c new file mode 100644 index 0000000000..4a232e27bf --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/mpa.c @@ -0,0 +1,549 @@ + +/* + * IBM Accurate Mathematical Library + * written by International Business Machines Corp. + * Copyright (C) 2001, 2006 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: mpa.c */ +/* */ +/* FUNCTIONS: */ +/* mcr */ +/* acr */ +/* cr */ +/* cpy */ +/* cpymn */ +/* norm */ +/* denorm */ +/* mp_dbl */ +/* dbl_mp */ +/* add_magnitudes */ +/* sub_magnitudes */ +/* add */ +/* sub */ +/* mul */ +/* inv */ +/* dvd */ +/* */ +/* Arithmetic functions for multiple precision numbers. */ +/* Relative errors are bounded */ +/************************************************************************/ + + +#include "endian.h" +#include "mpa.h" +#include "mpa2.h" +#include <sys/param.h> /* For MIN() */ +/* mcr() compares the sizes of the mantissas of two multiple precision */ +/* numbers. Mantissas are compared regardless of the signs of the */ +/* numbers, even if x->d[0] or y->d[0] are zero. Exponents are also */ +/* disregarded. */ +static int mcr(const mp_no *x, const mp_no *y, int p) { + long i; + long p2 = p; + for (i=1; i<=p2; i++) { + if (X[i] == Y[i]) continue; + else if (X[i] > Y[i]) return 1; + else return -1; } + return 0; +} + + + +/* acr() compares the absolute values of two multiple precision numbers */ +int __acr(const mp_no *x, const mp_no *y, int p) { + long i; + + if (X[0] == ZERO) { + if (Y[0] == ZERO) i= 0; + else i=-1; + } + else if (Y[0] == ZERO) i= 1; + else { + if (EX > EY) i= 1; + else if (EX < EY) i=-1; + else i= mcr(x,y,p); + } + + return i; +} + + +/* cr90 compares the values of two multiple precision numbers */ +int __cr(const mp_no *x, const mp_no *y, int p) { + int i; + + if (X[0] > Y[0]) i= 1; + else if (X[0] < Y[0]) i=-1; + else if (X[0] < ZERO ) i= __acr(y,x,p); + else i= __acr(x,y,p); + + return i; +} + + +/* Copy a multiple precision number. Set *y=*x. x=y is permissible. */ +void __cpy(const mp_no *x, mp_no *y, int p) { + long i; + + EY = EX; + for (i=0; i <= p; i++) Y[i] = X[i]; + + return; +} + + +/* Copy a multiple precision number x of precision m into a */ +/* multiple precision number y of precision n. In case n>m, */ +/* the digits of y beyond the m'th are set to zero. In case */ +/* n<m, the digits of x beyond the n'th are ignored. */ +/* x=y is permissible. */ + +void __cpymn(const mp_no *x, int m, mp_no *y, int n) { + + long i,k; + long n2 = n; + long m2 = m; + + EY = EX; k=MIN(m2,n2); + for (i=0; i <= k; i++) Y[i] = X[i]; + for ( ; i <= n2; i++) Y[i] = ZERO; + + return; +} + +/* Convert a multiple precision number *x into a double precision */ +/* number *y, normalized case (|x| >= 2**(-1022))) */ +static void norm(const mp_no *x, double *y, int p) +{ + #define R radixi.d + long i; +#if 0 + int k; +#endif + double a,c,u,v,z[5]; + if (p<5) { + if (p==1) c = X[1]; + else if (p==2) c = X[1] + R* X[2]; + else if (p==3) c = X[1] + R*(X[2] + R* X[3]); + else if (p==4) c =(X[1] + R* X[2]) + R*R*(X[3] + R*X[4]); + } + else { + for (a=ONE, z[1]=X[1]; z[1] < TWO23; ) + {a *= TWO; z[1] *= TWO; } + + for (i=2; i<5; i++) { + z[i] = X[i]*a; + u = (z[i] + CUTTER)-CUTTER; + if (u > z[i]) u -= RADIX; + z[i] -= u; + z[i-1] += u*RADIXI; + } + + u = (z[3] + TWO71) - TWO71; + if (u > z[3]) u -= TWO19; + v = z[3]-u; + + if (v == TWO18) { + if (z[4] == ZERO) { + for (i=5; i <= p; i++) { + if (X[i] == ZERO) continue; + else {z[3] += ONE; break; } + } + } + else z[3] += ONE; + } + + c = (z[1] + R *(z[2] + R * z[3]))/a; + } + + c *= X[0]; + + for (i=1; i<EX; i++) c *= RADIX; + for (i=1; i>EX; i--) c *= RADIXI; + + *y = c; + return; +#undef R +} + +/* Convert a multiple precision number *x into a double precision */ +/* number *y, denormalized case (|x| < 2**(-1022))) */ +static void denorm(const mp_no *x, double *y, int p) +{ + long i,k; + long p2 = p; + double c,u,z[5]; +#if 0 + double a,v; +#endif + +#define R radixi.d + if (EX<-44 || (EX==-44 && X[1]<TWO5)) + { *y=ZERO; return; } + + if (p2==1) { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=ZERO; z[3]=ZERO; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; z[3]=ZERO; k=2;} + else {z[1]= TWO10; z[2]=ZERO; z[3]=X[1]; k=1;} + } + else if (p2==2) { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=X[2]; z[3]=ZERO; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; z[3]=X[2]; k=2;} + else {z[1]= TWO10; z[2]=ZERO; z[3]=X[1]; k=1;} + } + else { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=X[2]; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; k=2;} + else {z[1]= TWO10; z[2]=ZERO; k=1;} + z[3] = X[k]; + } + + u = (z[3] + TWO57) - TWO57; + if (u > z[3]) u -= TWO5; + + if (u==z[3]) { + for (i=k+1; i <= p2; i++) { + if (X[i] == ZERO) continue; + else {z[3] += ONE; break; } + } + } + + c = X[0]*((z[1] + R*(z[2] + R*z[3])) - TWO10); + + *y = c*TWOM1032; + return; + +#undef R +} + +/* Convert a multiple precision number *x into a double precision number *y. */ +/* The result is correctly rounded to the nearest/even. *x is left unchanged */ + +void __mp_dbl(const mp_no *x, double *y, int p) { +#if 0 + int i,k; + double a,c,u,v,z[5]; +#endif + + if (X[0] == ZERO) {*y = ZERO; return; } + + if (EX> -42) norm(x,y,p); + else if (EX==-42 && X[1]>=TWO10) norm(x,y,p); + else denorm(x,y,p); +} + + +/* dbl_mp() converts a double precision number x into a multiple precision */ +/* number *y. If the precision p is too small the result is truncated. x is */ +/* left unchanged. */ + +void __dbl_mp(double x, mp_no *y, int p) { + + long i,n; + long p2 = p; + double u; + + /* Sign */ + if (x == ZERO) {Y[0] = ZERO; return; } + else if (x > ZERO) Y[0] = ONE; + else {Y[0] = MONE; x=-x; } + + /* Exponent */ + for (EY=ONE; x >= RADIX; EY += ONE) x *= RADIXI; + for ( ; x < ONE; EY -= ONE) x *= RADIX; + + /* Digits */ + n=MIN(p2,4); + for (i=1; i<=n; i++) { + u = (x + TWO52) - TWO52; + if (u>x) u -= ONE; + Y[i] = u; x -= u; x *= RADIX; } + for ( ; i<=p2; i++) Y[i] = ZERO; + return; +} + + +/* add_magnitudes() adds the magnitudes of *x & *y assuming that */ +/* abs(*x) >= abs(*y) > 0. */ +/* The sign of the sum *z is undefined. x&y may overlap but not x&z or y&z. */ +/* No guard digit is used. The result equals the exact sum, truncated. */ +/* *x & *y are left unchanged. */ + +static void add_magnitudes(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i,j,k; + long p2 = p; + + EZ = EX; + + i=p2; j=p2+ EY - EX; k=p2+1; + + if (j<1) + {__cpy(x,z,p); return; } + else Z[k] = ZERO; + + for (; j>0; i--,j--) { + Z[k] += X[i] + Y[j]; + if (Z[k] >= RADIX) { + Z[k] -= RADIX; + Z[--k] = ONE; } + else + Z[--k] = ZERO; + } + + for (; i>0; i--) { + Z[k] += X[i]; + if (Z[k] >= RADIX) { + Z[k] -= RADIX; + Z[--k] = ONE; } + else + Z[--k] = ZERO; + } + + if (Z[1] == ZERO) { + for (i=1; i<=p2; i++) Z[i] = Z[i+1]; } + else EZ += ONE; +} + + +/* sub_magnitudes() subtracts the magnitudes of *x & *y assuming that */ +/* abs(*x) > abs(*y) > 0. */ +/* The sign of the difference *z is undefined. x&y may overlap but not x&z */ +/* or y&z. One guard digit is used. The error is less than one ulp. */ +/* *x & *y are left unchanged. */ + +static void sub_magnitudes(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i,j,k; + long p2 = p; + + EZ = EX; + + if (EX == EY) { + i=j=k=p2; + Z[k] = Z[k+1] = ZERO; } + else { + j= EX - EY; + if (j > p2) {__cpy(x,z,p); return; } + else { + i=p2; j=p2+1-j; k=p2; + if (Y[j] > ZERO) { + Z[k+1] = RADIX - Y[j--]; + Z[k] = MONE; } + else { + Z[k+1] = ZERO; + Z[k] = ZERO; j--;} + } + } + + for (; j>0; i--,j--) { + Z[k] += (X[i] - Y[j]); + if (Z[k] < ZERO) { + Z[k] += RADIX; + Z[--k] = MONE; } + else + Z[--k] = ZERO; + } + + for (; i>0; i--) { + Z[k] += X[i]; + if (Z[k] < ZERO) { + Z[k] += RADIX; + Z[--k] = MONE; } + else + Z[--k] = ZERO; + } + + for (i=1; Z[i] == ZERO; i++) ; + EZ = EZ - i + 1; + for (k=1; i <= p2+1; ) + Z[k++] = Z[i++]; + for (; k <= p2; ) + Z[k++] = ZERO; + + return; +} + + +/* Add two multiple precision numbers. Set *z = *x + *y. x&y may overlap */ +/* but not x&z or y&z. One guard digit is used. The error is less than */ +/* one ulp. *x & *y are left unchanged. */ + +void __add(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + int n; + + if (X[0] == ZERO) {__cpy(y,z,p); return; } + else if (Y[0] == ZERO) {__cpy(x,z,p); return; } + + if (X[0] == Y[0]) { + if (__acr(x,y,p) > 0) {add_magnitudes(x,y,z,p); Z[0] = X[0]; } + else {add_magnitudes(y,x,z,p); Z[0] = Y[0]; } + } + else { + if ((n=__acr(x,y,p)) == 1) {sub_magnitudes(x,y,z,p); Z[0] = X[0]; } + else if (n == -1) {sub_magnitudes(y,x,z,p); Z[0] = Y[0]; } + else Z[0] = ZERO; + } + return; +} + + +/* Subtract two multiple precision numbers. *z is set to *x - *y. x&y may */ +/* overlap but not x&z or y&z. One guard digit is used. The error is */ +/* less than one ulp. *x & *y are left unchanged. */ + +void __sub(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + int n; + + if (X[0] == ZERO) {__cpy(y,z,p); Z[0] = -Z[0]; return; } + else if (Y[0] == ZERO) {__cpy(x,z,p); return; } + + if (X[0] != Y[0]) { + if (__acr(x,y,p) > 0) {add_magnitudes(x,y,z,p); Z[0] = X[0]; } + else {add_magnitudes(y,x,z,p); Z[0] = -Y[0]; } + } + else { + if ((n=__acr(x,y,p)) == 1) {sub_magnitudes(x,y,z,p); Z[0] = X[0]; } + else if (n == -1) {sub_magnitudes(y,x,z,p); Z[0] = -Y[0]; } + else Z[0] = ZERO; + } + return; +} + + +/* Multiply two multiple precision numbers. *z is set to *x * *y. x&y */ +/* may overlap but not x&z or y&z. In case p=1,2,3 the exact result is */ +/* truncated to p digits. In case p>3 the error is bounded by 1.001 ulp. */ +/* *x & *y are left unchanged. */ + +void __mul(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i, i1, i2, j, k, k2; + long p2 = p; + double u, zk, zk2; + + /* Is z=0? */ + if (X[0]*Y[0]==ZERO) + { Z[0]=ZERO; return; } + + /* Multiply, add and carry */ + k2 = (p2<3) ? p2+p2 : p2+3; + zk = Z[k2]=ZERO; + for (k=k2; k>1; ) { + if (k > p2) {i1=k-p2; i2=p2+1; } + else {i1=1; i2=k; } +#if 1 + /* rearange this inner loop to allow the fmadd instructions to be + independent and execute in parallel on processors that have + dual symetrical FP pipelines. */ + if (i1 < (i2-1)) + { + /* make sure we have at least 2 iterations */ + if (((i2 - i1) & 1L) == 1L) + { + /* Handle the odd iterations case. */ + zk2 = x->d[i2-1]*y->d[i1]; + } + else + zk2 = zero.d; + /* Do two multiply/adds per loop iteration, using independent + accumulators; zk and zk2. */ + for (i=i1,j=i2-1; i<i2-1; i+=2,j-=2) + { + zk += x->d[i]*y->d[j]; + zk2 += x->d[i+1]*y->d[j-1]; + } + zk += zk2; /* final sum. */ + } + else + { + /* Special case when iterations is 1. */ + zk += x->d[i1]*y->d[i1]; + } +#else + /* The orginal code. */ + for (i=i1,j=i2-1; i<i2; i++,j--) zk += X[i]*Y[j]; +#endif + + u = (zk + CUTTER)-CUTTER; + if (u > zk) u -= RADIX; + Z[k] = zk - u; + zk = u*RADIXI; + --k; + } + Z[k] = zk; + + /* Is there a carry beyond the most significant digit? */ + if (Z[1] == ZERO) { + for (i=1; i<=p2; i++) Z[i]=Z[i+1]; + EZ = EX + EY - 1; } + else + EZ = EX + EY; + + Z[0] = X[0] * Y[0]; + return; +} + + +/* Invert a multiple precision number. Set *y = 1 / *x. */ +/* Relative error bound = 1.001*r**(1-p) for p=2, 1.063*r**(1-p) for p=3, */ +/* 2.001*r**(1-p) for p>3. */ +/* *x=0 is not permissible. *x is left unchanged. */ + +void __inv(const mp_no *x, mp_no *y, int p) { + long i; +#if 0 + int l; +#endif + double t; + mp_no z,w; + static const int np1[] = {0,0,0,0,1,2,2,2,2,3,3,3,3,3,3,3,3,3, + 4,4,4,4,4,4,4,4,4,4,4,4,4,4,4}; + const mp_no mptwo = {1,{1.0,2.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}}; + + __cpy(x,&z,p); z.e=0; __mp_dbl(&z,&t,p); + t=ONE/t; __dbl_mp(t,y,p); EY -= EX; + + for (i=0; i<np1[p]; i++) { + __cpy(y,&w,p); + __mul(x,&w,y,p); + __sub(&mptwo,y,&z,p); + __mul(&w,&z,y,p); + } + return; +} + + +/* Divide one multiple precision number by another.Set *z = *x / *y. *x & *y */ +/* are left unchanged. x&y may overlap but not x&z or y&z. */ +/* Relative error bound = 2.001*r**(1-p) for p=2, 2.063*r**(1-p) for p=3 */ +/* and 3.001*r**(1-p) for p>3. *y=0 is not permissible. */ + +void __dvd(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + mp_no w; + + if (X[0] == ZERO) Z[0] = ZERO; + else {__inv(y,&w,p); __mul(x,&w,z,p);} + return; +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/slowpow.c b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/slowpow.c new file mode 100644 index 0000000000..ad147a89a6 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/fpu/slowpow.c @@ -0,0 +1,94 @@ +/* + * IBM Accurate Mathematical Library + * written by International Business Machines Corp. + * Copyright (C) 2001, 2006 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:slowpow.c */ +/* */ +/* FUNCTION:slowpow */ +/* */ +/*FILES NEEDED:mpa.h */ +/* mpa.c mpexp.c mplog.c halfulp.c */ +/* */ +/* Given two IEEE double machine numbers y,x , routine computes the */ +/* correctly rounded (to nearest) value of x^y. Result calculated by */ +/* multiplication (in halfulp.c) or if result isn't accurate enough */ +/* then routine converts x and y into multi-precision doubles and */ +/* recompute. */ +/*************************************************************************/ + +#include "mpa.h" +#include "math_private.h" + +void __mpexp (mp_no * x, mp_no * y, int p); +void __mplog (mp_no * x, mp_no * y, int p); +double ulog (double); +double __halfulp (double x, double y); + +double +__slowpow (double x, double y, double z) +{ + double res, res1; + long double ldw, ldz, ldpp; + static const long double ldeps = 0x4.0p-96; + + res = __halfulp (x, y); /* halfulp() returns -10 or x^y */ + if (res >= 0) + return res; /* if result was really computed by halfulp */ + /* else, if result was not really computed by halfulp */ + + /* Compute pow as long double, 106 bits */ + ldz = __ieee754_logl ((long double) x); + ldw = (long double) y *ldz; + ldpp = __ieee754_expl (ldw); + res = (double) (ldpp + ldeps); + res1 = (double) (ldpp - ldeps); + + if (res != res1) /* if result still not accurate enough */ + { /* use mpa for higher persision. */ + mp_no mpx, mpy, mpz, mpw, mpp, mpr, mpr1; + static const mp_no eps = { -3, {1.0, 4.0} }; + int p; + + p = 10; /* p=precision 240 bits */ + __dbl_mp (x, &mpx, p); + __dbl_mp (y, &mpy, p); + __dbl_mp (z, &mpz, p); + __mplog (&mpx, &mpz, p); /* log(x) = z */ + __mul (&mpy, &mpz, &mpw, p); /* y * z =w */ + __mpexp (&mpw, &mpp, p); /* e^w =pp */ + __add (&mpp, &eps, &mpr, p); /* pp+eps =r */ + __mp_dbl (&mpr, &res, p); + __sub (&mpp, &eps, &mpr1, p); /* pp -eps =r1 */ + __mp_dbl (&mpr1, &res1, p); /* converting into double precision */ + if (res == res1) + return res; + + /* if we get here result wasn't calculated exactly, continue for + more exact calculation using 768 bits. */ + p = 32; + __dbl_mp (x, &mpx, p); + __dbl_mp (y, &mpy, p); + __dbl_mp (z, &mpz, p); + __mplog (&mpx, &mpz, p); /* log(c)=z */ + __mul (&mpy, &mpz, &mpw, p); /* y*z =w */ + __mpexp (&mpw, &mpp, p); /* e^w=pp */ + __mp_dbl (&mpp, &res, p); /* converting into double precision */ + } + return res; +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/memcopy.h b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/memcopy.h new file mode 100644 index 0000000000..c05208da55 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/memcopy.h @@ -0,0 +1,113 @@ +/* memcopy.h -- definitions for memory copy functions. Generic C version. + Copyright (C) 1991, 1992, 1993, 1997, 2004, 2006 Free Software Foundation, Inc. + This file is part of the GNU C Library. + Contributed by Torbjorn Granlund (tege@sics.se). + + 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. */ + +/* The strategy of the memory functions is: + + 1. Copy bytes until the destination pointer is aligned. + + 2. Copy words in unrolled loops. If the source and destination + are not aligned in the same way, use word memory operations, + but shift and merge two read words before writing. + + 3. Copy the few remaining bytes. + + This is fast on processors that have at least 10 registers for + allocation by GCC, and that can access memory at reg+const in one + instruction. + + I made an "exhaustive" test of this memmove when I wrote it, + exhaustive in the sense that I tried all alignment and length + combinations, with and without overlap. */ + +#include <sysdeps/generic/memcopy.h> + +/* The macros defined in this file are: + + BYTE_COPY_FWD(dst_beg_ptr, src_beg_ptr, nbytes_to_copy) + + BYTE_COPY_BWD(dst_end_ptr, src_end_ptr, nbytes_to_copy) + + WORD_COPY_FWD(dst_beg_ptr, src_beg_ptr, nbytes_remaining, nbytes_to_copy) + + WORD_COPY_BWD(dst_end_ptr, src_end_ptr, nbytes_remaining, nbytes_to_copy) + + MERGE(old_word, sh_1, new_word, sh_2) + [I fail to understand. I feel stupid. --roland] +*/ + + +/* Threshold value for when to enter the unrolled loops. */ +#undef OP_T_THRES +#define OP_T_THRES 16 + +/* Copy exactly NBYTES bytes from SRC_BP to DST_BP, + without any assumptions about alignment of the pointers. */ +#undef BYTE_COPY_FWD +#define BYTE_COPY_FWD(dst_bp, src_bp, nbytes) \ + do \ + { \ + size_t __nbytes = (nbytes); \ + if (__nbytes & 1) \ + { \ + ((byte *) dst_bp)[0] = ((byte *) src_bp)[0]; \ + src_bp += 1; \ + dst_bp += 1; \ + __nbytes -= 1; \ + } \ + while (__nbytes > 0) \ + { \ + byte __x = ((byte *) src_bp)[0]; \ + byte __y = ((byte *) src_bp)[1]; \ + src_bp += 2; \ + __nbytes -= 2; \ + ((byte *) dst_bp)[0] = __x; \ + ((byte *) dst_bp)[1] = __y; \ + dst_bp += 2; \ + } \ + } while (0) + +/* Copy exactly NBYTES_TO_COPY bytes from SRC_END_PTR to DST_END_PTR, + beginning at the bytes right before the pointers and continuing towards + smaller addresses. Don't assume anything about alignment of the + pointers. */ +#undef BYTE_COPY_BWD +#define BYTE_COPY_BWD(dst_ep, src_ep, nbytes) \ + do \ + { \ + size_t __nbytes = (nbytes); \ + if (__nbytes & 1) \ + { \ + src_ep -= 1; \ + dst_ep -= 1; \ + ((byte *) dst_ep)[0] = ((byte *) src_ep)[0]; \ + __nbytes -= 1; \ + } \ + while (__nbytes > 0) \ + { \ + byte __x, __y; \ + src_ep -= 2; \ + __y = ((byte *) src_ep)[1]; \ + __x = ((byte *) src_ep)[0]; \ + dst_ep -= 2; \ + __nbytes -= 2; \ + ((byte *) dst_ep)[1] = __y; \ + ((byte *) dst_ep)[0] = __x; \ + } \ + } while (0) diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/wordcopy.c b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/wordcopy.c new file mode 100644 index 0000000000..f71b41dc42 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power4/wordcopy.c @@ -0,0 +1,209 @@ +/* _memcopy.c -- subroutines for memory copy functions. + Copyright (C) 1991, 1996 Free Software Foundation, Inc. + This file is part of the GNU C Library. + Contributed by Torbjorn Granlund (tege@sics.se). + + 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. */ + +/* BE VERY CAREFUL IF YOU CHANGE THIS CODE...! */ + +#include <stddef.h> +#include <memcopy.h> + +/* _wordcopy_fwd_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + Both SRCP and DSTP should be aligned for memory operations on `op_t's. */ + +void +_wordcopy_fwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + do + { + a0 = ((op_t *) srcp)[0]; + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = a0; + ((op_t *) dstp)[1] = a1; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_fwd_dest_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + DSTP should be aligned for memory operations on `op_t's, but SRCP must + *not* be aligned. */ + +void +_wordcopy_fwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make SRCP aligned by rounding it down to the beginning of the `op_t' + it points in the middle of. */ + srcp &= -OPSIZ; + a0 = ((op_t *) srcp)[0]; + + if (len & 1) + { + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = MERGE (a0, sh_1, a1, sh_2); + + if (len == 1) + return; + + a0 = a1; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, sh_1, a1, sh_2); + ((op_t *) dstp)[1] = MERGE (a1, sh_1, a2, sh_2); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_bwd_aligned -- Copy block finishing right before + SRCP to block finishing right before DSTP with LEN `op_t' words + (not LEN bytes!). Both SRCP and DSTP should be aligned for memory + operations on `op_t's. */ + +void +_wordcopy_bwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + len -= 1; + } + + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = a1; + ((op_t *) dstp)[0] = a0; + + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_bwd_dest_aligned -- Copy block finishing right + before SRCP to block finishing right before DSTP with LEN `op_t' + words (not LEN bytes!). DSTP should be aligned for memory + operations on `op_t', but SRCP must *not* be aligned. */ + +void +_wordcopy_bwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make srcp aligned by rounding it down to the beginning of the op_t + it points in the middle of. */ + srcp &= -OPSIZ; + a2 = ((op_t *) srcp)[0]; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + a1 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[0] = MERGE (a1, sh_1, a2, sh_2); + + if (len == 1) + return; + + a2 = a1; + len -= 1; + } + + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, sh_1, a2, sh_2); + ((op_t *) dstp)[0] = MERGE (a0, sh_1, a1, sh_2); + a2 = a0; + + len -= 2; + } + while (len != 0); +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/Implies index 4e3a983426..a0ddab4863 100644 --- a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/Implies +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/Implies @@ -1 +1,2 @@ powerpc/powerpc32/power4 +powerpc/powerpc32/powerpc64 diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llround.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llround.S new file mode 100644 index 0000000000..a2171fe09a --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llround.S @@ -0,0 +1,60 @@ +/* lround function. POWER5+, PowerPC32 version. + Copyright (C) 2006 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., 1 Franklin Street, Fifth Floor, Boston MA + 02110-1301 USA. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + +/* long [r3] llround (float x [fp1]) + IEEE 1003.1 lround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we pre-round using the V2.02 Floating Round to Integer Nearest + instruction before we use the Floating Convert to Integer Word with + round to zero instruction. */ + + .machine "power5" +ENTRY (__llround) + stwu r1,-16(r1) + cfi_adjust_cfa_offset (16) + frin fp2,fp1 + fctidz fp3,fp2 /* Convert To Integer Word lround toward 0. */ + stfd fp3,8(r1) + nop /* Ensure the following load is in a different dispatch */ + nop /* group to avoid pipe stall on POWER4&5. */ + nop + lwz r4,12(r1) + lwz r3,8(r1) + addi r1,r1,16 + blr + END (__llround) + +weak_alias (__llround, llround) + +strong_alias (__llround, __llroundf) +weak_alias (__llround, llroundf) + +#ifdef NO_LONG_DOUBLE +weak_alias (__llround, llroundl) +strong_alias (__llround, __llroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __llround, llroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llroundf.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llroundf.S new file mode 100644 index 0000000000..ffe6b7eb37 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_llroundf.S @@ -0,0 +1,2 @@ +/* __llroundf is in s_llround.S */ +/* __llroundf is in s_llround.S */ diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_lround.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_lround.S new file mode 100644 index 0000000000..83107f6f97 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5+/fpu/s_lround.S @@ -0,0 +1,58 @@ +/* lround function. POWER5+, PowerPC32 version. + Copyright (C) 2006 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., 1 Franklin Street, Fifth Floor, Boston MA + 02110-1301 USA. */ +#include <sysdep.h> +#include <math_ldbl_opt.h> + +/* long [r3] lround (float x [fp1]) + IEEE 1003.1 lround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we pre-round using the V2.02 Floating Round to Integer Nearest + instruction before we use the Floating Convert to Integer Word with + round to zero instruction. */ + + .machine "power5" +ENTRY (__lround) + stwu r1,-16(r1) + cfi_adjust_cfa_offset (16) + frin fp2,fp1 + fctiwz fp3,fp2 /* Convert To Integer Word lround toward 0. */ + stfd fp3,8(r1) + nop /* Ensure the following load is in a different dispatch */ + nop /* group to avoid pipe stall on POWER4&5. */ + nop + lwz r3,12(r1) + addi r1,r1,16 + blr + END (__lround) + +weak_alias (__lround, lround) + +strong_alias (__lround, __lroundf) +weak_alias (__lround, lroundf) + +#ifdef NO_LONG_DOUBLE +weak_alias (__lround, lroundl) +strong_alias (__lround, __lroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __lround, lroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5/Implies index 4e3a983426..a0ddab4863 100644 --- a/powerpc-cpu/sysdeps/powerpc/powerpc32/power5/Implies +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power5/Implies @@ -1 +1,2 @@ powerpc/powerpc32/power4 +powerpc/powerpc32/powerpc64 diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/Implies index 59fcd3c09d..f488e227d3 100644 --- a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/Implies +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/Implies @@ -1 +1,2 @@ powerpc/powerpc32/power5+ +powerpc/powerpc32/powerpc64 diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/memcpy.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/memcpy.S new file mode 100644 index 0000000000..e8d56eb135 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/memcpy.S @@ -0,0 +1,842 @@ +/* Optimized memcpy implementation for PowerPC32 on POWER6. + Copyright (C) 2003, 2006 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., 51 Franklin Street, Fifth Floor, Boston MA + 02110-1301 USA. */ + +#include <sysdep.h> +#include <bp-sym.h> +#include <bp-asm.h> + +/* __ptr_t [r3] memcpy (__ptr_t dst [r3], __ptr_t src [r4], size_t len [r5]); + Returns 'dst'. + + Memcpy handles short copies (< 32-bytes) using a binary move blocks + (no loops) of lwz/stw. The tail (remaining 1-3) bytes is handled + with the appropriate combination of byte and halfword load/stores. + There is minimal effort to optimize the alignment of short moves. + + Longer moves (>= 32-bytes) justify the effort to get at least the + destination word (4-byte) aligned. Further optimization is + possible when both source and destination are word aligned. + Each case has an optimized unrolled loop. */ + +EALIGN (BP_SYM (memcpy), 5, 0) + CALL_MCOUNT + + stwu 1,-32(1) + cfi_adjust_cfa_offset(32) + cmplwi cr1,5,31 /* check for short move. */ + neg 0,3 + cmplwi cr1,5,31 + clrlwi 10,4,30 /* check alignment of src. */ + andi. 11,3,3 /* check alignment of dst. */ + clrlwi 0,0,30 /* Number of bytes until the 1st word of dst. */ + ble- cr1,L(word_unaligned_short) /* If move < 32 bytes. */ + cmplw cr6,10,11 + stw 31,24(1) + cfi_offset(31,(24-32)) + stw 30,20(1) + cfi_offset(30,(20-32)) + mr 30,3 + beq .L0 + mtcrf 0x01,0 + subf 31,0,5 /* Length after alignment. */ + add 12,4,0 /* Compute src addr after alignment. */ + /* Move 0-3 bytes as needed to get the destination word aligned. */ +1: bf 31,2f + lbz 6,0(4) + bf 30,3f + lhz 7,1(4) + stb 6,0(3) + sth 7,1(3) + addi 3,3,3 + b 0f +3: + stb 6,0(3) + addi 3,3,1 + b 0f +2: bf 30,0f + lhz 6,0(4) + sth 6,0(3) + addi 3,3,2 +0: + clrlwi 10,12,30 /* check alignment of src again. */ + srwi 9,31,2 /* Number of full words remaining. */ + bne- cr6,L(wdu) /* If source is not word aligned. .L6 */ + clrlwi 11,31,30 /* calculate the number of tail bytes */ + b L(word_aligned) + /* Copy words from source to destination, assuming the destination is + aligned on a word boundary. + + At this point we know there are at least 29 bytes left (32-3) to copy. + The next step is to determine if the source is also word aligned. + If not branch to the unaligned move code at .L6. which uses + a load, shift, store strategy. + + Otherwise source and destination are word aligned, and we can use + the optimized word copy loop. */ + .align 4 +.L0: + mr 31,5 + mr 12,4 + bne- cr6,L(wdu) /* If source is not word aligned. .L6 */ + srwi 9,5,2 /* Number of full words remaining. */ + clrlwi 11,5,30 /* calculate the number of tail bytes */ + + /* Move words where destination and source are word aligned. + Use an unrolled loop to copy 4 words (16-bytes) per iteration. + If the the copy is not an exact multiple of 16 bytes, 1-3 + words are copied as needed to set up the main loop. After + the main loop exits there may be a tail of 1-3 bytes. These bytes are + copied a halfword/byte at a time as needed to preserve alignment. */ +L(word_aligned): + mtcrf 0x01,9 + srwi 8,31,4 /* calculate the 16 byte loop count */ + cmplwi cr1,9,4 + cmplwi cr6,11,0 + mr 11,12 + + bf 30,1f + lwz 6,0(12) + lwz 7,4(12) + addi 11,12,8 + mtctr 8 + stw 6,0(3) + stw 7,4(3) + addi 10,3,8 + bf 31,4f + lwz 0,8(12) + stw 0,8(3) + blt cr1,3f + addi 11,12,12 + addi 10,3,12 + b 4f + .align 4 +1: + mr 10,3 + mtctr 8 + bf 31,4f + lwz 6,0(12) + addi 11,12,4 + stw 6,0(3) + addi 10,3,4 + + .align 4 +4: + lwz 6,0(11) + lwz 7,4(11) + lwz 8,8(11) + lwz 0,12(11) + stw 6,0(10) + stw 7,4(10) + stw 8,8(10) + stw 0,12(10) + addi 11,11,16 + addi 10,10,16 + bdnz 4b +3: + clrrwi 0,31,2 + mtcrf 0x01,31 + beq cr6,0f +.L9: + add 3,3,0 + add 12,12,0 + +/* At this point we have a tail of 0-3 bytes and we know that the + destination is word aligned. */ +2: bf 30,1f + lhz 6,0(12) + addi 12,12,2 + sth 6,0(3) + addi 3,3,2 +1: bf 31,0f + lbz 6,0(12) + stb 6,0(3) +0: + /* Return original dst pointer. */ + mr 3,30 + lwz 30,20(1) + lwz 31,24(1) + addi 1,1,32 + blr + +/* Copy up to 31 bytes. This divided into two cases 0-8 bytes and 9-31 + bytes. Each case is handled without loops, using binary (1,2,4,8) + tests. + + In the short (0-8 byte) case no attempt is made to force alignment + of either source or destination. The hardware will handle the + unaligned load/stores with small delays for crossing 32- 128-byte, + and 4096-byte boundaries. Since these short moves are unlikely to be + unaligned or cross these boundaries, the overhead to force + alignment is not justified. + + The longer (9-31 byte) move is more likely to cross 32- or 128-byte + boundaries. Since only loads are sensitive to the 32-/128-byte + boundaries it is more important to align the source then the + destination. If the source is not already word aligned, we first + move 1-3 bytes as needed. Since we are only word aligned we don't + use double word load/stores to insure that all loads are aligned. + While the destination and stores may still be unaligned, this + is only an issue for page (4096 byte boundary) crossing, which + should be rare for these short moves. The hardware handles this + case automatically with a small (~20 cycle) delay. */ + .align 4 + + cfi_same_value (31) + cfi_same_value (30) +L(word_unaligned_short): + mtcrf 0x01,5 + cmplwi cr6,5,8 + neg 8,4 + clrrwi 9,4,2 + andi. 0,8,3 + beq cr6,L(wus_8) /* Handle moves of 8 bytes. */ +/* At least 9 bytes left. Get the source word aligned. */ + cmpldi cr1,5,16 + mr 12,4 + ble cr6,L(wus_4) /* Handle moves of 0-8 bytes. */ + mr 11,3 + mr 10,5 + cmplwi cr6,0,2 + beq L(wus_tail) /* If the source is already word aligned skip this. */ +/* Copy 1-3 bytes to get source address word aligned. */ + lwz 6,0(9) + subf 10,0,5 + add 12,4,0 + blt cr6,5f + srdi 7,6,16 + bgt cr6,3f + sth 6,0(3) + b 7f + .align 4 +3: + stb 7,0(3) + sth 6,1(3) + b 7f + .align 4 +5: + stb 6,0(3) +7: + cmplwi cr1,10,16 + add 11,3,0 + mtcrf 0x01,10 + .align 4 +L(wus_tail): +/* At least 6 bytes left and the source is word aligned. This allows + some speculative loads up front. */ +/* We need to special case the fall-through because the biggest delays + are due to address computation not being ready in time for the + AGEN. */ + lwz 6,0(12) + lwz 7,4(12) + blt cr1,L(wus_tail8) + cmplwi cr0,10,24 +L(wus_tail16): /* Move 16 bytes. */ + stw 6,0(11) + stw 7,4(11) + lwz 6,8(12) + lwz 7,12(12) + stw 6,8(11) + stw 7,12(11) +/* Move 8 bytes more. */ + bf 28,L(wus_tail16p8) + cmplwi cr1,10,28 + lwz 6,16(12) + lwz 7,20(12) + stw 6,16(11) + stw 7,20(11) +/* Move 4 bytes more. */ + bf 29,L(wus_tail16p4) + lwz 6,24(12) + stw 6,24(11) + addi 12,12,28 + addi 11,11,28 + bgt cr1,L(wus_tail2) + /* exactly 28 bytes. Return original dst pointer and exit. */ + addi 1,1,32 + blr + .align 4 +L(wus_tail16p8): /* less then 8 bytes left. */ + beq cr1,L(wus_tailX) /* exactly 16 bytes, early exit. */ + cmplwi cr1,10,20 + bf 29,L(wus_tail16p2) +/* Move 4 bytes more. */ + lwz 6,16(12) + stw 6,16(11) + addi 12,12,20 + addi 11,11,20 + bgt cr1,L(wus_tail2) + /* exactly 20 bytes. Return original dst pointer and exit. */ + addi 1,1,32 + blr + .align 4 +L(wus_tail16p4): /* less then 4 bytes left. */ + addi 12,12,24 + addi 11,11,24 + bgt cr0,L(wus_tail2) + /* exactly 24 bytes. Return original dst pointer and exit. */ + addi 1,1,32 + blr + .align 4 +L(wus_tail16p2): /* 16 bytes moved, less then 4 bytes left. */ + addi 12,12,16 + addi 11,11,16 + b L(wus_tail2) + + .align 4 +L(wus_tail8): /* Move 8 bytes. */ +/* r6, r7 already loaded speculatively. */ + cmplwi cr1,10,8 + cmplwi cr0,10,12 + bf 28,L(wus_tail4) + stw 6,0(11) + stw 7,4(11) +/* Move 4 bytes more. */ + bf 29,L(wus_tail8p4) + lwz 6,8(12) + stw 6,8(11) + addi 12,12,12 + addi 11,11,12 + bgt cr0,L(wus_tail2) + /* exactly 12 bytes. Return original dst pointer and exit. */ + addi 1,1,32 + blr + .align 4 +L(wus_tail8p4): /* less then 4 bytes left. */ + addi 12,12,8 + addi 11,11,8 + bgt cr1,L(wus_tail2) + /* exactly 8 bytes. Return original dst pointer and exit. */ + addi 1,1,32 + blr + + .align 4 +L(wus_tail4): /* Move 4 bytes. */ +/* r6 already loaded speculatively. If we are here we know there is + more then 4 bytes left. So there is no need to test. */ + addi 12,12,4 + stw 6,0(11) + addi 11,11,4 +L(wus_tail2): /* Move 2-3 bytes. */ + bf 30,L(wus_tail1) + lhz 6,0(12) + sth 6,0(11) + bf 31,L(wus_tailX) + lbz 7,2(12) + stb 7,2(11) + addi 1,1,32 + blr +L(wus_tail1): /* Move 1 byte. */ + bf 31,L(wus_tailX) + lbz 6,0(12) + stb 6,0(11) +L(wus_tailX): + /* Return original dst pointer. */ + addi 1,1,32 + blr + +/* Special case to copy 0-8 bytes. */ + .align 4 +L(wus_8): + lwz 6,0(4) + lwz 7,4(4) + stw 6,0(3) + stw 7,4(3) + /* Return original dst pointer. */ + addi 1,1,32 + blr + .align 4 +L(wus_4): + bf 29,L(wus_2) + lwz 6,0(4) + stw 6,0(3) + bf 30,L(wus_5) + lhz 7,4(4) + sth 7,4(3) + bf 31,L(wus_0) + lbz 8,6(4) + stb 8,6(3) + addi 1,1,32 + blr + .align 4 +L(wus_5): + bf 31,L(wus_0) + lbz 6,4(4) + stb 6,4(3) + /* Return original dst pointer. */ + addi 1,1,32 + blr + .align 4 +L(wus_2): /* Move 2-3 bytes. */ + bf 30,L(wus_1) + lhz 6,0(4) + sth 6,0(3) + bf 31,L(wus_0) + lbz 7,2(4) + stb 7,2(3) + addi 1,1,32 + blr + .align 4 +L(wus_1): /* Move 1 byte. */ + bf 31,L(wus_0) + lbz 6,0(4) + stb 6,0(3) + .align 3 +L(wus_0): + /* Return original dst pointer. */ + addi 1,1,32 + blr + + .align 4 + cfi_offset(31,(24-32)) + cfi_offset(30,(20-32)) +L(wdu): + + /* Copy words where the destination is aligned but the source is + not. For power4, power5 and power6 machines there is penalty for + unaligned loads (src) that cross 32-byte, cacheline, or page + boundaries. So we want to use simple (unaligned) loads where + posible but avoid them where we know the load would span a 32-byte + boundary. + + At this point we know we have at least 29 (32-3) bytes to copy + the src is unaligned. and we may cross at least one 32-byte + boundary. Also we have the following regester values: + r3 == adjusted dst, word aligned + r4 == unadjusted src + r5 == unadjusted len + r9 == adjusted Word length + r10 == src alignment (1-3) + r12 == adjuested src, not aligned + r31 == adjusted len + + First we need to copy word upto but not crossing the next 32-byte + boundary. Then perform aligned loads just before and just after + the boundary and use shifts and or to gernerate the next aligned + word for dst. If more then 32 bytes remain we copy (unaligned src) + the next 7 words and repeat the loop until less then 32-bytes + remaim. + + Then if more then 4 bytes remain we again use aligned loads, + shifts and or to generate the next dst word. We then process the + remaining words using unaligned loads as needed. Finally we check + if there more then 0 bytes (1-3) bytes remainting and use + halfword and or byte load/stores to complete the copy. +*/ + mr 4,12 /* restore unaligned adjusted src ptr */ + clrlwi 0,12,27 /* Find dist from previous 32-byte boundary. */ + slwi 10,10,3 /* calculate number of bits to shift 1st word left */ + cmplwi cr5,0,16 + subfic 8,0,32 /* Number of bytes to next 32-byte boundary. */ + + mtcrf 0x01,8 + cmplwi cr1,10,16 + subfic 9,10,32 /* number of bits to shift 2nd word right */ +/* This test is reversed because the timing to compare the bytes to + 32-byte boundary could not be meet. So we compare the bytes from + previous 32-byte boundary and invert the test. */ + bge cr5,L(wdu_h32_8) + .align 4 + lwz 6,0(4) + lwz 7,4(4) + addi 12,4,16 /* generate alternate pointers to avoid agen */ + addi 11,3,16 /* timing issues downstream. */ + stw 6,0(3) + stw 7,4(3) + subi 31,31,16 + lwz 6,8(4) + lwz 7,12(4) + addi 4,4,16 + stw 6,8(3) + stw 7,12(3) + addi 3,3,16 + bf 28,L(wdu_h32_4) + lwz 6,0(12) + lwz 7,4(12) + subi 31,31,8 + addi 4,4,8 + stw 6,0(11) + stw 7,4(11) + addi 3,3,8 + bf 29,L(wdu_h32_0) + lwz 6,8(12) + addi 4,4,4 + subi 31,31,4 + stw 6,8(11) + addi 3,3,4 + b L(wdu_h32_0) + .align 4 +L(wdu_h32_8): + bf 28,L(wdu_h32_4) + lwz 6,0(4) + lwz 7,4(4) + subi 31,31,8 + bf 29,L(wdu_h32_8x) + stw 6,0(3) + stw 7,4(3) + lwz 6,8(4) + addi 4,4,12 + subi 31,31,4 + stw 6,8(3) + addi 3,3,12 + b L(wdu_h32_0) + .align 4 +L(wdu_h32_8x): + addi 4,4,8 + stw 6,0(3) + stw 7,4(3) + addi 3,3,8 + b L(wdu_h32_0) + .align 4 +L(wdu_h32_4): + bf 29,L(wdu_h32_0) + lwz 6,0(4) + subi 31,31,4 + addi 4,4,4 + stw 6,0(3) + addi 3,3,4 + .align 4 +L(wdu_h32_0): +/* set up for 32-byte boundry crossing word move and possibly 32-byte + move loop. */ + clrrwi 12,4,2 + cmplwi cr5,31,32 + bge cr1,L(wdu2_32) +#if 0 + b L(wdu1_32) +/* + cmplwi cr1,10,8 + beq cr1,L(wdu1_32) + cmplwi cr1,10,16 + beq cr1,L(wdu2_32) + cmplwi cr1,10,24 + beq cr1,L(wdu3_32) +*/ +L(wdu_32): + lwz 6,0(12) + cmplwi cr6,31,4 + srwi 8,31,5 /* calculate the 32 byte loop count */ + slw 0,6,10 + clrlwi 31,31,27 /* The remaining bytes, < 32. */ + blt cr5,L(wdu_32tail) + mtctr 8 + cmplwi cr6,31,4 + .align 4 +L(wdu_loop32): + /* copy 32 bytes at a time */ + lwz 8,4(12) + addi 12,12,32 + lwz 7,4(4) + srw 8,8,9 + or 0,0,8 + stw 0,0(3) + stw 7,4(3) + lwz 6,8(4) + lwz 7,12(4) + stw 6,8(3) + stw 7,12(3) + lwz 6,16(4) + lwz 7,20(4) + stw 6,16(3) + stw 7,20(3) + lwz 6,24(4) + lwz 7,28(4) + lwz 8,0(12) + addi 4,4,32 + stw 6,24(3) + stw 7,28(3) + addi 3,3,32 + slw 0,8,10 + bdnz+ L(wdu_loop32) + +L(wdu_32tail): + mtcrf 0x01,31 + cmplwi cr5,31,16 + blt cr6,L(wdu_4tail) + /* calculate and store the final word */ + lwz 8,4(12) + srw 8,8,9 + or 6,0,8 + b L(wdu_32tailx) +#endif + .align 4 +L(wdu1_32): + lwz 6,-1(4) + cmplwi cr6,31,4 + srwi 8,31,5 /* calculate the 32 byte loop count */ + slwi 6,6,8 + clrlwi 31,31,27 /* The remaining bytes, < 32. */ + blt cr5,L(wdu1_32tail) + mtctr 8 + cmplwi cr6,31,4 + + lwz 8,3(4) + lwz 7,4(4) +/* Equivalent to: srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,8,(32-8),31 + b L(wdu1_loop32x) + .align 4 +L(wdu1_loop32): + /* copy 32 bytes at a time */ + lwz 8,3(4) + lwz 7,4(4) + stw 10,-8(3) + stw 11,-4(3) +/* Equivalent to srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,8,(32-8),31 +L(wdu1_loop32x): + lwz 10,8(4) + lwz 11,12(4) + stw 6,0(3) + stw 7,4(3) + lwz 6,16(4) + lwz 7,20(4) + stw 10,8(3) + stw 11,12(3) + lwz 10,24(4) + lwz 11,28(4) + lwz 8,32-1(4) + addi 4,4,32 + stw 6,16(3) + stw 7,20(3) + addi 3,3,32 + slwi 6,8,8 + bdnz+ L(wdu1_loop32) + stw 10,-8(3) + stw 11,-4(3) + +L(wdu1_32tail): + mtcrf 0x01,31 + cmplwi cr5,31,16 + blt cr6,L(wdu_4tail) + /* calculate and store the final word */ + lwz 8,3(4) +/* Equivalent to: srwi 8,8,32-9; or 6,6,8 */ + rlwimi 6,8,8,(32-8),31 + b L(wdu_32tailx) + +L(wdu2_32): + bgt cr1,L(wdu3_32) + lwz 6,-2(4) + cmplwi cr6,31,4 + srwi 8,31,5 /* calculate the 32 byte loop count */ + slwi 6,6,16 + clrlwi 31,31,27 /* The remaining bytes, < 32. */ + blt cr5,L(wdu2_32tail) + mtctr 8 + cmplwi cr6,31,4 + + lwz 8,2(4) + lwz 7,4(4) +/* Equivalent to: srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,16,(32-16),31 + b L(wdu2_loop32x) + .align 4 +L(wdu2_loop32): + /* copy 32 bytes at a time */ + lwz 8,2(4) + lwz 7,4(4) + stw 10,-8(3) + stw 11,-4(3) +/* Equivalent to srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,16,(32-16),31 +L(wdu2_loop32x): + lwz 10,8(4) + lwz 11,12(4) + stw 6,0(3) + stw 7,4(3) + lwz 6,16(4) + lwz 7,20(4) + stw 10,8(3) + stw 11,12(3) + lwz 10,24(4) + lwz 11,28(4) +/* lwz 8,0(12) */ + lwz 8,32-2(4) + addi 4,4,32 + stw 6,16(3) + stw 7,20(3) + addi 3,3,32 + slwi 6,8,16 + bdnz+ L(wdu2_loop32) + stw 10,-8(3) + stw 11,-4(3) + +L(wdu2_32tail): + mtcrf 0x01,31 + cmplwi cr5,31,16 + blt cr6,L(wdu_4tail) + /* calculate and store the final word */ + lwz 8,2(4) +/* Equivalent to: srwi 8,8,32-9; or 6,6,8 */ + rlwimi 6,8,16,(32-16),31 + b L(wdu_32tailx) + +L(wdu3_32): +/* lwz 6,0(12) */ + lwz 6,-3(4) + cmplwi cr6,31,4 + srwi 8,31,5 /* calculate the 32 byte loop count */ + slwi 6,6,24 + clrlwi 31,31,27 /* The remaining bytes, < 32. */ + blt cr5,L(wdu3_32tail) + mtctr 8 + cmplwi cr6,31,4 + + lwz 8,1(4) + lwz 7,4(4) +/* Equivalent to: srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,24,(32-24),31 + b L(wdu3_loop32x) + .align 4 +L(wdu3_loop32): + /* copy 32 bytes at a time */ + lwz 8,1(4) + lwz 7,4(4) + stw 10,-8(3) + stw 11,-4(3) +/* Equivalent to srwi 8,8,32-8; or 6,6,8 */ + rlwimi 6,8,24,(32-24),31 +L(wdu3_loop32x): + lwz 10,8(4) + lwz 11,12(4) + stw 6,0(3) + stw 7,4(3) + lwz 6,16(4) + lwz 7,20(4) + stw 10,8(3) + stw 11,12(3) + lwz 10,24(4) + lwz 11,28(4) + lwz 8,32-3(4) + addi 4,4,32 + stw 6,16(3) + stw 7,20(3) + addi 3,3,32 + slwi 6,8,24 + bdnz+ L(wdu3_loop32) + stw 10,-8(3) + stw 11,-4(3) + +L(wdu3_32tail): + mtcrf 0x01,31 + cmplwi cr5,31,16 + blt cr6,L(wdu_4tail) + /* calculate and store the final word */ + lwz 8,1(4) +/* Equivalent to: srwi 8,8,32-9; or 6,6,8 */ + rlwimi 6,8,24,(32-24),31 + b L(wdu_32tailx) + .align 4 +L(wdu_32tailx): + blt cr5,L(wdu_t32_8) + lwz 7,4(4) + addi 12,4,16 /* generate alternate pointers to avoid agen */ + addi 11,3,16 /* timing issues downstream. */ + stw 6,0(3) + stw 7,4(3) + subi 31,31,16 + lwz 6,8(4) + lwz 7,12(4) + addi 4,4,16 + stw 6,8(3) + stw 7,12(3) + addi 3,3,16 + bf 28,L(wdu_t32_4x) + lwz 6,0(12) + lwz 7,4(12) + addi 4,4,8 + subi 31,31,8 + stw 6,0(11) + stw 7,4(11) + addi 3,3,8 + bf 29,L(wdu_t32_0) + lwz 6,8(12) + addi 4,4,4 + subi 31,31,4 + stw 6,8(11) + addi 3,3,4 + b L(wdu_t32_0) + .align 4 +L(wdu_t32_4x): + bf 29,L(wdu_t32_0) + lwz 6,0(4) + addi 4,4,4 + subi 31,31,4 + stw 6,0(3) + addi 3,3,4 + b L(wdu_t32_0) + .align 4 +L(wdu_t32_8): + bf 28,L(wdu_t32_4) + lwz 7,4(4) + subi 31,31,8 + bf 29,L(wdu_t32_8x) + stw 6,0(3) + stw 7,4(3) + lwz 6,8(4) + subi 31,31,4 + addi 4,4,12 + stw 6,8(3) + addi 3,3,12 + b L(wdu_t32_0) + .align 4 +L(wdu_t32_8x): + addi 4,4,8 + stw 6,0(3) + stw 7,4(3) + addi 3,3,8 + b L(wdu_t32_0) + .align 4 +L(wdu_t32_4): + subi 31,31,4 + stw 6,0(3) + addi 4,4,4 + addi 3,3,4 + .align 4 +L(wdu_t32_0): +L(wdu_4tail): + cmplwi cr6,31,0 + beq cr6,L(wdus_0) /* If the tail is 0 bytes we are done! */ + bf 30,L(wdus_3) + lhz 7,0(4) + sth 7,0(3) + bf 31,L(wdus_0) + lbz 8,2(4) + stb 8,2(3) + mr 3,30 + lwz 30,20(1) + lwz 31,24(1) + addi 1,1,32 + blr + .align 4 +L(wdus_3): + bf 31,L(wus_0) + lbz 6,0(4) + stb 6,0(3) + .align 4 +L(wdus_0): + /* Return original dst pointer. */ + mr 3,30 + lwz 30,20(1) + lwz 31,24(1) + addi 1,1,32 + blr +END (BP_SYM (memcpy)) + +libc_hidden_builtin_def (memcpy) diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/wordcopy.c b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/wordcopy.c new file mode 100644 index 0000000000..ddf28659f9 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6/wordcopy.c @@ -0,0 +1,287 @@ +/* _memcopy.c -- subroutines for memory copy functions. + Copyright (C) 1991, 1996, 2006 Free Software Foundation, Inc. + This file is part of the GNU C Library. + Contributed by Torbjorn Granlund (tege@sics.se). + Updated for POWER6 by Steven Munroe (sjmunroe@us.ibm.com). + + 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. */ + +/* BE VERY CAREFUL IF YOU CHANGE THIS CODE...! */ + +#include <stddef.h> +#include <memcopy.h> + +/* _wordcopy_fwd_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + Both SRCP and DSTP should be aligned for memory operations on `op_t's. */ + +void +_wordcopy_fwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + do + { + a0 = ((op_t *) srcp)[0]; + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = a0; + ((op_t *) dstp)[1] = a1; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_fwd_dest_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + DSTP should be aligned for memory operations on `op_t's, but SRCP must + *not* be aligned. */ + +void +_wordcopy_fwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + int align; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + align = srcp % OPSIZ; + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make SRCP aligned by rounding it down to the beginning of the `op_t' + it points in the middle of. */ + srcp &= -OPSIZ; + a0 = ((op_t *) srcp)[0]; + + if (len & 1) + { + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = MERGE (a0, sh_1, a1, sh_2); + + if (len == 1) + return; + + a0 = a1; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + switch (align) + { + case 1: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 8, a1, (32-8)); + ((op_t *) dstp)[1] = MERGE (a1, 8, a2, (32-8)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 2: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 16, a1, (32-16)); + ((op_t *) dstp)[1] = MERGE (a1, 16, a2, (32-16)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 3: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 24, a1, (32-24)); + ((op_t *) dstp)[1] = MERGE (a1, 24, a2, (32-24)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + } + +} + +/* _wordcopy_bwd_aligned -- Copy block finishing right before + SRCP to block finishing right before DSTP with LEN `op_t' words + (not LEN bytes!). Both SRCP and DSTP should be aligned for memory + operations on `op_t's. */ + +void +_wordcopy_bwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + len -= 1; + } + + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = a1; + ((op_t *) dstp)[0] = a0; + + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_bwd_dest_aligned -- Copy block finishing right + before SRCP to block finishing right before DSTP with LEN `op_t' + words (not LEN bytes!). DSTP should be aligned for memory + operations on `op_t', but SRCP must *not* be aligned. */ + +void +_wordcopy_bwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + int align; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + align = srcp % OPSIZ; + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make srcp aligned by rounding it down to the beginning of the op_t + it points in the middle of. */ + srcp &= -OPSIZ; + a2 = ((op_t *) srcp)[0]; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + a1 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[0] = MERGE (a1, sh_1, a2, sh_2); + + if (len == 1) + return; + + a2 = a1; + len -= 1; + } + + switch (align) + { + case 1: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 8, a2, (32-8)); + ((op_t *) dstp)[0] = MERGE (a0, 8, a1, (32-8)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 2: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 16, a2, (32-16)); + ((op_t *) dstp)[0] = MERGE (a0, 16, a1, (32-16)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 3: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 24, a2, (32-24)); + ((op_t *) dstp)[0] = MERGE (a0, 24, a1, (32-24)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + } +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/Implies new file mode 100644 index 0000000000..d041ecab15 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/Implies @@ -0,0 +1,3 @@ +powerpc/powerpc32/power6 +powerpc/powerpc32/power5+ +powerpc/powerpc32/powerpc64 diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/Implies new file mode 100644 index 0000000000..da50634d7a --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/Implies @@ -0,0 +1,3 @@ +powerpc/powerpc32/power6/fpu +powerpc/powerpc32/power5+/fpu +powerpc/powerpc32/powerpc64/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lrint.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lrint.S new file mode 100644 index 0000000000..76d8432920 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lrint.S @@ -0,0 +1,42 @@ +/* Round double to long int. POWER6x PowerPC32 version. + Copyright (C) 2006 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. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + + .machine "power6" +/* long int[r3] __lrint (double x[fp1]) */ +ENTRY (__lrint) + fctiw fp13,fp1 + mftgpr r3,fp13 + blr + END (__lrint) + +weak_alias (__lrint, lrint) + +strong_alias (__lrint, __lrintf) +weak_alias (__lrint, lrintf) + +#ifdef NO_LONG_DOUBLE +strong_alias (__lrint, __lrintl) +weak_alias (__lrint, lrintl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __lrint, lrintl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lround.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lround.S new file mode 100644 index 0000000000..dc0ace821a --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/power6x/fpu/s_lround.S @@ -0,0 +1,52 @@ +/* lround function. POWER6x, PowerPC32 version. + Copyright (C) 2006 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., 1 Franklin Street, Fifth Floor, Boston MA + 02110-1301 USA. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + +/* long [r3] lround (float x [fp1]) + IEEE 1003.1 lround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we pre-round using the V2.02 Floating Round to Integer Nearest + instruction before we use the Floating Convert to Integer Word with + round to zero instruction. */ + + .machine "power6" +ENTRY (__lround) + frin fp2,fp1 /* Pre-round +-0.5. */ + fctiwz fp3,fp2 /* Convert To Integer Word lround toward 0. */ + mftgpr r3,fp3 /* Transfer fpr3 to r3. */ + blr + END (__lround) + +weak_alias (__lround, lround) + +strong_alias (__lround, __lroundf) +weak_alias (__lround, lroundf) + +#ifdef NO_LONG_DOUBLE +weak_alias (__lround, lroundl) +strong_alias (__lround, __lroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __lround, lroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llround.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llround.S new file mode 100644 index 0000000000..952d2aa6a5 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llround.S @@ -0,0 +1,97 @@ +/* llround function. PowerPC32 on PowerPC64 version. + Copyright (C) 2004, 2006 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., 1 Franklin Street, Fifth Floor, Boston MA + 02110-1301 USA. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + + .section .rodata.cst8,"aM",@progbits,8 + .align 2 +.LC0: /* 0.0 */ + .long 0x00000000 +.LC1: /* 0.5 */ + .long 0x3f000000 + + .section ".text" + +/* long [r3] lround (float x [fp1]) + IEEE 1003.1 lround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we can't use the PowerPC "round to Nearest" mode. Instead we set + "round toward Zero" mode and round by adding +-0.5 before rounding + to the integer value. */ + +ENTRY (__llround) + stwu r1,-16(r1) + cfi_adjust_cfa_offset (16) +#ifdef SHARED + mflr r11 + cfi_register(lr,r11) +# ifdef HAVE_ASM_PPC_REL16 + bcl 20,31,1f +1: mflr r9 + addis r9,r9,.LC0-1b@ha + addi r9,r9,.LC0-1b@l +# else + bl _GLOBAL_OFFSET_TABLE_@local-4 + mflr r10 + lwz r9,.LC0@got(10) +# endif + mtlr r11 + cfi_same_value (lr) + lfs fp12,0(r9) + lfs fp10,.LC1-.LC0(r9) +#else + lis r9,.LC0@ha + lis r10,.LC1@ha + lfs fp12,.LC0@l(r9) + lfs fp10,.LC1@l(r10) +#endif + fcmpu cr6,fp1,fp12 /* if (x > 0.0) */ + ble- cr6,.L4 + fadd fp1,fp1,fp10 /* x+= 0.5; */ +.L9: + fctidz fp2,fp1 /* Convert To Integer DW round toward 0. */ + stfd fp2,8(r1) + nop /* Ensure the following load is in a different dispatch */ + nop /* group to avoid pipe stall on POWER4&5. */ + nop + lwz r4,12(r1) + lwz r3,8(r1) + addi r1,r1,16 + blr +.L4: + fsub fp1,fp1,fp10 /* x-= 0.5; */ + b .L9 + END (__llround) + +weak_alias (__llround, llround) + +strong_alias (__llround, __llroundf) +weak_alias (__llround, llroundf) + +#ifdef NO_LONG_DOUBLE +weak_alias (__llround, llroundl) +strong_alias (__llround, __llroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __llround, llroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llroundf.S b/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llroundf.S new file mode 100644 index 0000000000..72d6181541 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc32/powerpc64/fpu/s_llroundf.S @@ -0,0 +1 @@ +/* __llroundf is in s_llround.S */ diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/970/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/970/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/970/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/Makefile b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/Makefile new file mode 100644 index 0000000000..60aa508ba4 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/Makefile @@ -0,0 +1,6 @@ +# Makefile fragment for POWER4/5/5+. + +ifeq ($(subdir),string) +CFLAGS-wordcopy.c += --param max-variable-expansions-in-unroller=2 --param max-unroll-times=2 -funroll-loops -fpeel-loops -ftree-loop-linear +CFLAGS-memmove.c += --param max-variable-expansions-in-unroller=2 --param max-unroll-times=2 -funroll-loops -fpeel-loops -ftree-loop-linear +endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/Makefile b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/Makefile new file mode 100644 index 0000000000..89dfa5ef35 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/Makefile @@ -0,0 +1,5 @@ +# Makefile fragment for POWER4/5/5+ platforms with FPU. + +ifeq ($(subdir),math) +CFLAGS-mpa.c += --param max-unroll-times=4 -funroll-loops -fpeel-loops -ftree-loop-linear +endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/mpa.c b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/mpa.c new file mode 100644 index 0000000000..4a232e27bf --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/mpa.c @@ -0,0 +1,549 @@ + +/* + * IBM Accurate Mathematical Library + * written by International Business Machines Corp. + * Copyright (C) 2001, 2006 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: mpa.c */ +/* */ +/* FUNCTIONS: */ +/* mcr */ +/* acr */ +/* cr */ +/* cpy */ +/* cpymn */ +/* norm */ +/* denorm */ +/* mp_dbl */ +/* dbl_mp */ +/* add_magnitudes */ +/* sub_magnitudes */ +/* add */ +/* sub */ +/* mul */ +/* inv */ +/* dvd */ +/* */ +/* Arithmetic functions for multiple precision numbers. */ +/* Relative errors are bounded */ +/************************************************************************/ + + +#include "endian.h" +#include "mpa.h" +#include "mpa2.h" +#include <sys/param.h> /* For MIN() */ +/* mcr() compares the sizes of the mantissas of two multiple precision */ +/* numbers. Mantissas are compared regardless of the signs of the */ +/* numbers, even if x->d[0] or y->d[0] are zero. Exponents are also */ +/* disregarded. */ +static int mcr(const mp_no *x, const mp_no *y, int p) { + long i; + long p2 = p; + for (i=1; i<=p2; i++) { + if (X[i] == Y[i]) continue; + else if (X[i] > Y[i]) return 1; + else return -1; } + return 0; +} + + + +/* acr() compares the absolute values of two multiple precision numbers */ +int __acr(const mp_no *x, const mp_no *y, int p) { + long i; + + if (X[0] == ZERO) { + if (Y[0] == ZERO) i= 0; + else i=-1; + } + else if (Y[0] == ZERO) i= 1; + else { + if (EX > EY) i= 1; + else if (EX < EY) i=-1; + else i= mcr(x,y,p); + } + + return i; +} + + +/* cr90 compares the values of two multiple precision numbers */ +int __cr(const mp_no *x, const mp_no *y, int p) { + int i; + + if (X[0] > Y[0]) i= 1; + else if (X[0] < Y[0]) i=-1; + else if (X[0] < ZERO ) i= __acr(y,x,p); + else i= __acr(x,y,p); + + return i; +} + + +/* Copy a multiple precision number. Set *y=*x. x=y is permissible. */ +void __cpy(const mp_no *x, mp_no *y, int p) { + long i; + + EY = EX; + for (i=0; i <= p; i++) Y[i] = X[i]; + + return; +} + + +/* Copy a multiple precision number x of precision m into a */ +/* multiple precision number y of precision n. In case n>m, */ +/* the digits of y beyond the m'th are set to zero. In case */ +/* n<m, the digits of x beyond the n'th are ignored. */ +/* x=y is permissible. */ + +void __cpymn(const mp_no *x, int m, mp_no *y, int n) { + + long i,k; + long n2 = n; + long m2 = m; + + EY = EX; k=MIN(m2,n2); + for (i=0; i <= k; i++) Y[i] = X[i]; + for ( ; i <= n2; i++) Y[i] = ZERO; + + return; +} + +/* Convert a multiple precision number *x into a double precision */ +/* number *y, normalized case (|x| >= 2**(-1022))) */ +static void norm(const mp_no *x, double *y, int p) +{ + #define R radixi.d + long i; +#if 0 + int k; +#endif + double a,c,u,v,z[5]; + if (p<5) { + if (p==1) c = X[1]; + else if (p==2) c = X[1] + R* X[2]; + else if (p==3) c = X[1] + R*(X[2] + R* X[3]); + else if (p==4) c =(X[1] + R* X[2]) + R*R*(X[3] + R*X[4]); + } + else { + for (a=ONE, z[1]=X[1]; z[1] < TWO23; ) + {a *= TWO; z[1] *= TWO; } + + for (i=2; i<5; i++) { + z[i] = X[i]*a; + u = (z[i] + CUTTER)-CUTTER; + if (u > z[i]) u -= RADIX; + z[i] -= u; + z[i-1] += u*RADIXI; + } + + u = (z[3] + TWO71) - TWO71; + if (u > z[3]) u -= TWO19; + v = z[3]-u; + + if (v == TWO18) { + if (z[4] == ZERO) { + for (i=5; i <= p; i++) { + if (X[i] == ZERO) continue; + else {z[3] += ONE; break; } + } + } + else z[3] += ONE; + } + + c = (z[1] + R *(z[2] + R * z[3]))/a; + } + + c *= X[0]; + + for (i=1; i<EX; i++) c *= RADIX; + for (i=1; i>EX; i--) c *= RADIXI; + + *y = c; + return; +#undef R +} + +/* Convert a multiple precision number *x into a double precision */ +/* number *y, denormalized case (|x| < 2**(-1022))) */ +static void denorm(const mp_no *x, double *y, int p) +{ + long i,k; + long p2 = p; + double c,u,z[5]; +#if 0 + double a,v; +#endif + +#define R radixi.d + if (EX<-44 || (EX==-44 && X[1]<TWO5)) + { *y=ZERO; return; } + + if (p2==1) { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=ZERO; z[3]=ZERO; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; z[3]=ZERO; k=2;} + else {z[1]= TWO10; z[2]=ZERO; z[3]=X[1]; k=1;} + } + else if (p2==2) { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=X[2]; z[3]=ZERO; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; z[3]=X[2]; k=2;} + else {z[1]= TWO10; z[2]=ZERO; z[3]=X[1]; k=1;} + } + else { + if (EX==-42) {z[1]=X[1]+TWO10; z[2]=X[2]; k=3;} + else if (EX==-43) {z[1]= TWO10; z[2]=X[1]; k=2;} + else {z[1]= TWO10; z[2]=ZERO; k=1;} + z[3] = X[k]; + } + + u = (z[3] + TWO57) - TWO57; + if (u > z[3]) u -= TWO5; + + if (u==z[3]) { + for (i=k+1; i <= p2; i++) { + if (X[i] == ZERO) continue; + else {z[3] += ONE; break; } + } + } + + c = X[0]*((z[1] + R*(z[2] + R*z[3])) - TWO10); + + *y = c*TWOM1032; + return; + +#undef R +} + +/* Convert a multiple precision number *x into a double precision number *y. */ +/* The result is correctly rounded to the nearest/even. *x is left unchanged */ + +void __mp_dbl(const mp_no *x, double *y, int p) { +#if 0 + int i,k; + double a,c,u,v,z[5]; +#endif + + if (X[0] == ZERO) {*y = ZERO; return; } + + if (EX> -42) norm(x,y,p); + else if (EX==-42 && X[1]>=TWO10) norm(x,y,p); + else denorm(x,y,p); +} + + +/* dbl_mp() converts a double precision number x into a multiple precision */ +/* number *y. If the precision p is too small the result is truncated. x is */ +/* left unchanged. */ + +void __dbl_mp(double x, mp_no *y, int p) { + + long i,n; + long p2 = p; + double u; + + /* Sign */ + if (x == ZERO) {Y[0] = ZERO; return; } + else if (x > ZERO) Y[0] = ONE; + else {Y[0] = MONE; x=-x; } + + /* Exponent */ + for (EY=ONE; x >= RADIX; EY += ONE) x *= RADIXI; + for ( ; x < ONE; EY -= ONE) x *= RADIX; + + /* Digits */ + n=MIN(p2,4); + for (i=1; i<=n; i++) { + u = (x + TWO52) - TWO52; + if (u>x) u -= ONE; + Y[i] = u; x -= u; x *= RADIX; } + for ( ; i<=p2; i++) Y[i] = ZERO; + return; +} + + +/* add_magnitudes() adds the magnitudes of *x & *y assuming that */ +/* abs(*x) >= abs(*y) > 0. */ +/* The sign of the sum *z is undefined. x&y may overlap but not x&z or y&z. */ +/* No guard digit is used. The result equals the exact sum, truncated. */ +/* *x & *y are left unchanged. */ + +static void add_magnitudes(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i,j,k; + long p2 = p; + + EZ = EX; + + i=p2; j=p2+ EY - EX; k=p2+1; + + if (j<1) + {__cpy(x,z,p); return; } + else Z[k] = ZERO; + + for (; j>0; i--,j--) { + Z[k] += X[i] + Y[j]; + if (Z[k] >= RADIX) { + Z[k] -= RADIX; + Z[--k] = ONE; } + else + Z[--k] = ZERO; + } + + for (; i>0; i--) { + Z[k] += X[i]; + if (Z[k] >= RADIX) { + Z[k] -= RADIX; + Z[--k] = ONE; } + else + Z[--k] = ZERO; + } + + if (Z[1] == ZERO) { + for (i=1; i<=p2; i++) Z[i] = Z[i+1]; } + else EZ += ONE; +} + + +/* sub_magnitudes() subtracts the magnitudes of *x & *y assuming that */ +/* abs(*x) > abs(*y) > 0. */ +/* The sign of the difference *z is undefined. x&y may overlap but not x&z */ +/* or y&z. One guard digit is used. The error is less than one ulp. */ +/* *x & *y are left unchanged. */ + +static void sub_magnitudes(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i,j,k; + long p2 = p; + + EZ = EX; + + if (EX == EY) { + i=j=k=p2; + Z[k] = Z[k+1] = ZERO; } + else { + j= EX - EY; + if (j > p2) {__cpy(x,z,p); return; } + else { + i=p2; j=p2+1-j; k=p2; + if (Y[j] > ZERO) { + Z[k+1] = RADIX - Y[j--]; + Z[k] = MONE; } + else { + Z[k+1] = ZERO; + Z[k] = ZERO; j--;} + } + } + + for (; j>0; i--,j--) { + Z[k] += (X[i] - Y[j]); + if (Z[k] < ZERO) { + Z[k] += RADIX; + Z[--k] = MONE; } + else + Z[--k] = ZERO; + } + + for (; i>0; i--) { + Z[k] += X[i]; + if (Z[k] < ZERO) { + Z[k] += RADIX; + Z[--k] = MONE; } + else + Z[--k] = ZERO; + } + + for (i=1; Z[i] == ZERO; i++) ; + EZ = EZ - i + 1; + for (k=1; i <= p2+1; ) + Z[k++] = Z[i++]; + for (; k <= p2; ) + Z[k++] = ZERO; + + return; +} + + +/* Add two multiple precision numbers. Set *z = *x + *y. x&y may overlap */ +/* but not x&z or y&z. One guard digit is used. The error is less than */ +/* one ulp. *x & *y are left unchanged. */ + +void __add(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + int n; + + if (X[0] == ZERO) {__cpy(y,z,p); return; } + else if (Y[0] == ZERO) {__cpy(x,z,p); return; } + + if (X[0] == Y[0]) { + if (__acr(x,y,p) > 0) {add_magnitudes(x,y,z,p); Z[0] = X[0]; } + else {add_magnitudes(y,x,z,p); Z[0] = Y[0]; } + } + else { + if ((n=__acr(x,y,p)) == 1) {sub_magnitudes(x,y,z,p); Z[0] = X[0]; } + else if (n == -1) {sub_magnitudes(y,x,z,p); Z[0] = Y[0]; } + else Z[0] = ZERO; + } + return; +} + + +/* Subtract two multiple precision numbers. *z is set to *x - *y. x&y may */ +/* overlap but not x&z or y&z. One guard digit is used. The error is */ +/* less than one ulp. *x & *y are left unchanged. */ + +void __sub(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + int n; + + if (X[0] == ZERO) {__cpy(y,z,p); Z[0] = -Z[0]; return; } + else if (Y[0] == ZERO) {__cpy(x,z,p); return; } + + if (X[0] != Y[0]) { + if (__acr(x,y,p) > 0) {add_magnitudes(x,y,z,p); Z[0] = X[0]; } + else {add_magnitudes(y,x,z,p); Z[0] = -Y[0]; } + } + else { + if ((n=__acr(x,y,p)) == 1) {sub_magnitudes(x,y,z,p); Z[0] = X[0]; } + else if (n == -1) {sub_magnitudes(y,x,z,p); Z[0] = -Y[0]; } + else Z[0] = ZERO; + } + return; +} + + +/* Multiply two multiple precision numbers. *z is set to *x * *y. x&y */ +/* may overlap but not x&z or y&z. In case p=1,2,3 the exact result is */ +/* truncated to p digits. In case p>3 the error is bounded by 1.001 ulp. */ +/* *x & *y are left unchanged. */ + +void __mul(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + long i, i1, i2, j, k, k2; + long p2 = p; + double u, zk, zk2; + + /* Is z=0? */ + if (X[0]*Y[0]==ZERO) + { Z[0]=ZERO; return; } + + /* Multiply, add and carry */ + k2 = (p2<3) ? p2+p2 : p2+3; + zk = Z[k2]=ZERO; + for (k=k2; k>1; ) { + if (k > p2) {i1=k-p2; i2=p2+1; } + else {i1=1; i2=k; } +#if 1 + /* rearange this inner loop to allow the fmadd instructions to be + independent and execute in parallel on processors that have + dual symetrical FP pipelines. */ + if (i1 < (i2-1)) + { + /* make sure we have at least 2 iterations */ + if (((i2 - i1) & 1L) == 1L) + { + /* Handle the odd iterations case. */ + zk2 = x->d[i2-1]*y->d[i1]; + } + else + zk2 = zero.d; + /* Do two multiply/adds per loop iteration, using independent + accumulators; zk and zk2. */ + for (i=i1,j=i2-1; i<i2-1; i+=2,j-=2) + { + zk += x->d[i]*y->d[j]; + zk2 += x->d[i+1]*y->d[j-1]; + } + zk += zk2; /* final sum. */ + } + else + { + /* Special case when iterations is 1. */ + zk += x->d[i1]*y->d[i1]; + } +#else + /* The orginal code. */ + for (i=i1,j=i2-1; i<i2; i++,j--) zk += X[i]*Y[j]; +#endif + + u = (zk + CUTTER)-CUTTER; + if (u > zk) u -= RADIX; + Z[k] = zk - u; + zk = u*RADIXI; + --k; + } + Z[k] = zk; + + /* Is there a carry beyond the most significant digit? */ + if (Z[1] == ZERO) { + for (i=1; i<=p2; i++) Z[i]=Z[i+1]; + EZ = EX + EY - 1; } + else + EZ = EX + EY; + + Z[0] = X[0] * Y[0]; + return; +} + + +/* Invert a multiple precision number. Set *y = 1 / *x. */ +/* Relative error bound = 1.001*r**(1-p) for p=2, 1.063*r**(1-p) for p=3, */ +/* 2.001*r**(1-p) for p>3. */ +/* *x=0 is not permissible. *x is left unchanged. */ + +void __inv(const mp_no *x, mp_no *y, int p) { + long i; +#if 0 + int l; +#endif + double t; + mp_no z,w; + static const int np1[] = {0,0,0,0,1,2,2,2,2,3,3,3,3,3,3,3,3,3, + 4,4,4,4,4,4,4,4,4,4,4,4,4,4,4}; + const mp_no mptwo = {1,{1.0,2.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0, + 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}}; + + __cpy(x,&z,p); z.e=0; __mp_dbl(&z,&t,p); + t=ONE/t; __dbl_mp(t,y,p); EY -= EX; + + for (i=0; i<np1[p]; i++) { + __cpy(y,&w,p); + __mul(x,&w,y,p); + __sub(&mptwo,y,&z,p); + __mul(&w,&z,y,p); + } + return; +} + + +/* Divide one multiple precision number by another.Set *z = *x / *y. *x & *y */ +/* are left unchanged. x&y may overlap but not x&z or y&z. */ +/* Relative error bound = 2.001*r**(1-p) for p=2, 2.063*r**(1-p) for p=3 */ +/* and 3.001*r**(1-p) for p>3. *y=0 is not permissible. */ + +void __dvd(const mp_no *x, const mp_no *y, mp_no *z, int p) { + + mp_no w; + + if (X[0] == ZERO) Z[0] = ZERO; + else {__inv(y,&w,p); __mul(x,&w,z,p);} + return; +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/slowpow.c b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/slowpow.c new file mode 100644 index 0000000000..ad147a89a6 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/fpu/slowpow.c @@ -0,0 +1,94 @@ +/* + * IBM Accurate Mathematical Library + * written by International Business Machines Corp. + * Copyright (C) 2001, 2006 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:slowpow.c */ +/* */ +/* FUNCTION:slowpow */ +/* */ +/*FILES NEEDED:mpa.h */ +/* mpa.c mpexp.c mplog.c halfulp.c */ +/* */ +/* Given two IEEE double machine numbers y,x , routine computes the */ +/* correctly rounded (to nearest) value of x^y. Result calculated by */ +/* multiplication (in halfulp.c) or if result isn't accurate enough */ +/* then routine converts x and y into multi-precision doubles and */ +/* recompute. */ +/*************************************************************************/ + +#include "mpa.h" +#include "math_private.h" + +void __mpexp (mp_no * x, mp_no * y, int p); +void __mplog (mp_no * x, mp_no * y, int p); +double ulog (double); +double __halfulp (double x, double y); + +double +__slowpow (double x, double y, double z) +{ + double res, res1; + long double ldw, ldz, ldpp; + static const long double ldeps = 0x4.0p-96; + + res = __halfulp (x, y); /* halfulp() returns -10 or x^y */ + if (res >= 0) + return res; /* if result was really computed by halfulp */ + /* else, if result was not really computed by halfulp */ + + /* Compute pow as long double, 106 bits */ + ldz = __ieee754_logl ((long double) x); + ldw = (long double) y *ldz; + ldpp = __ieee754_expl (ldw); + res = (double) (ldpp + ldeps); + res1 = (double) (ldpp - ldeps); + + if (res != res1) /* if result still not accurate enough */ + { /* use mpa for higher persision. */ + mp_no mpx, mpy, mpz, mpw, mpp, mpr, mpr1; + static const mp_no eps = { -3, {1.0, 4.0} }; + int p; + + p = 10; /* p=precision 240 bits */ + __dbl_mp (x, &mpx, p); + __dbl_mp (y, &mpy, p); + __dbl_mp (z, &mpz, p); + __mplog (&mpx, &mpz, p); /* log(x) = z */ + __mul (&mpy, &mpz, &mpw, p); /* y * z =w */ + __mpexp (&mpw, &mpp, p); /* e^w =pp */ + __add (&mpp, &eps, &mpr, p); /* pp+eps =r */ + __mp_dbl (&mpr, &res, p); + __sub (&mpp, &eps, &mpr1, p); /* pp -eps =r1 */ + __mp_dbl (&mpr1, &res1, p); /* converting into double precision */ + if (res == res1) + return res; + + /* if we get here result wasn't calculated exactly, continue for + more exact calculation using 768 bits. */ + p = 32; + __dbl_mp (x, &mpx, p); + __dbl_mp (y, &mpy, p); + __dbl_mp (z, &mpz, p); + __mplog (&mpx, &mpz, p); /* log(c)=z */ + __mul (&mpy, &mpz, &mpw, p); /* y*z =w */ + __mpexp (&mpw, &mpp, p); /* e^w=pp */ + __mp_dbl (&mpp, &res, p); /* converting into double precision */ + } + return res; +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/memcopy.h b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/memcopy.h new file mode 100644 index 0000000000..9a4ff79f4a --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/memcopy.h @@ -0,0 +1 @@ +#include "../../powerpc32/power4/memcopy.h" diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/wordcopy.c b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/wordcopy.c new file mode 100644 index 0000000000..f427b48e7a --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power4/wordcopy.c @@ -0,0 +1 @@ +#include "../../powerpc32/power4/wordcopy.c" diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/s_llround.S b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/s_llround.S new file mode 100644 index 0000000000..6a022769f4 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5+/fpu/s_llround.S @@ -0,0 +1,59 @@ +/* llround function. POWER5+, PowerPC64 version. + Copyright (C) 2006 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. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + +/* long long [r3] llround (float x [fp1]) + IEEE 1003.1 llround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we pre-round using the V2.02 Floating Round to Integer Nearest + instruction before we use Floating Convert to Integer Word with + round to zero instruction. */ + + .machine "power5" +EALIGN (__llround, 4, 0) + CALL_MCOUNT 0 + frin fp2, fp1 /* Round to nearest +-0.5. */ + fctidz fp3, fp2 /* Convert To Integer DW round toward 0. */ + stfd fp3, -16(r1) + nop /* Insure the following load is in a different dispatch group */ + nop /* to avoid pipe stall on POWER4&5. */ + nop + ld r3, -16(r1) + blr + END (__llround) + +strong_alias (__llround, __lround) +weak_alias (__llround, llround) +weak_alias (__lround, lround) + +#ifdef NO_LONG_DOUBLE +weak_alias (__llround, llroundl) +strong_alias (__llround, __llroundl) +weak_alias (__lround, lroundl) +strong_alias (__lround, __lroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __llround, llroundl, GLIBC_2_1) +compat_symbol (libm, __lround, lroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power5/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power5/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/fpu/Implies index f09854edb6..885e39c003 100644 --- a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/fpu/Implies +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/fpu/Implies @@ -1 +1,2 @@ powerpc/powerpc64/power5+/fpu +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/wordcopy.c b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/wordcopy.c new file mode 100644 index 0000000000..faddd945b3 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6/wordcopy.c @@ -0,0 +1,410 @@ +/* _memcopy.c -- subroutines for memory copy functions. + Copyright (C) 1991, 1996 Free Software Foundation, Inc. + This file is part of the GNU C Library. + Contributed by Torbjorn Granlund (tege@sics.se). + + 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. */ + +/* BE VERY CAREFUL IF YOU CHANGE THIS CODE...! */ + +#include <stddef.h> +#include <memcopy.h> + +/* _wordcopy_fwd_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + Both SRCP and DSTP should be aligned for memory operations on `op_t's. */ + +void +_wordcopy_fwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + do + { + a0 = ((op_t *) srcp)[0]; + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = a0; + ((op_t *) dstp)[1] = a1; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_fwd_dest_aligned -- Copy block beginning at SRCP to + block beginning at DSTP with LEN `op_t' words (not LEN bytes!). + DSTP should be aligned for memory operations on `op_t's, but SRCP must + *not* be aligned. */ + +void +_wordcopy_fwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + int align; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + align = srcp % OPSIZ; + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make SRCP aligned by rounding it down to the beginning of the `op_t' + it points in the middle of. */ + srcp &= -OPSIZ; + a0 = ((op_t *) srcp)[0]; + + if (len & 1) + { + a1 = ((op_t *) srcp)[1]; + ((op_t *) dstp)[0] = MERGE (a0, sh_1, a1, sh_2); + + if (len == 1) + return; + + a0 = a1; + srcp += OPSIZ; + dstp += OPSIZ; + len -= 1; + } + + switch (align) + { + case 1: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 8, a1, (64-8)); + ((op_t *) dstp)[1] = MERGE (a1, 8, a2, (64-8)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 2: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 16, a1, (64-16)); + ((op_t *) dstp)[1] = MERGE (a1, 16, a2, (64-16)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 3: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 24, a1, (64-24)); + ((op_t *) dstp)[1] = MERGE (a1, 24, a2, (64-24)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 4: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 32, a1, (64-32)); + ((op_t *) dstp)[1] = MERGE (a1, 32, a2, (64-32)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 5: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 40, a1, (64-40)); + ((op_t *) dstp)[1] = MERGE (a1, 40, a2, (64-40)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 6: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 48, a1, (64-48)); + ((op_t *) dstp)[1] = MERGE (a1, 48, a2, (64-48)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + case 7: + do + { + a1 = ((op_t *) srcp)[1]; + a2 = ((op_t *) srcp)[2]; + ((op_t *) dstp)[0] = MERGE (a0, 56, a1, (64-56)); + ((op_t *) dstp)[1] = MERGE (a1, 56, a2, (64-56)); + a0 = a2; + + srcp += 2 * OPSIZ; + dstp += 2 * OPSIZ; + len -= 2; + } + while (len != 0); + break; + } + +} + +/* _wordcopy_bwd_aligned -- Copy block finishing right before + SRCP to block finishing right before DSTP with LEN `op_t' words + (not LEN bytes!). Both SRCP and DSTP should be aligned for memory + operations on `op_t's. */ + +void +_wordcopy_bwd_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + ((op_t *) dstp)[0] = ((op_t *) srcp)[0]; + + if (len == 1) + return; + len -= 1; + } + + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = a1; + ((op_t *) dstp)[0] = a0; + + len -= 2; + } + while (len != 0); +} + +/* _wordcopy_bwd_dest_aligned -- Copy block finishing right + before SRCP to block finishing right before DSTP with LEN `op_t' + words (not LEN bytes!). DSTP should be aligned for memory + operations on `op_t', but SRCP must *not* be aligned. */ + +void +_wordcopy_bwd_dest_aligned (dstp, srcp, len) + long int dstp; + long int srcp; + size_t len; +{ + op_t a0, a1, a2; + int sh_1, sh_2; + int align; + + /* Calculate how to shift a word read at the memory operation + aligned srcp to make it aligned for copy. */ + + align = srcp % OPSIZ; + sh_1 = 8 * (srcp % OPSIZ); + sh_2 = 8 * OPSIZ - sh_1; + + /* Make srcp aligned by rounding it down to the beginning of the op_t + it points in the middle of. */ + srcp &= -OPSIZ; + a2 = ((op_t *) srcp)[0]; + + if (len & 1) + { + srcp -= OPSIZ; + dstp -= OPSIZ; + a1 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[0] = MERGE (a1, sh_1, a2, sh_2); + + if (len == 1) + return; + + a2 = a1; + len -= 1; + } + + switch (align) + { + case 1: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 8, a2, (64-8)); + ((op_t *) dstp)[0] = MERGE (a0, 8, a1, (64-8)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 2: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 16, a2, (64-16)); + ((op_t *) dstp)[0] = MERGE (a0, 16, a1, (64-16)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 3: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 24, a2, (64-24)); + ((op_t *) dstp)[0] = MERGE (a0, 24, a1, (64-24)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 4: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 32, a2, (64-32)); + ((op_t *) dstp)[0] = MERGE (a0, 32, a1, (64-32)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 5: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 40, a2, (64-40)); + ((op_t *) dstp)[0] = MERGE (a0, 40, a1, (64-40)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 6: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 48, a2, (64-48)); + ((op_t *) dstp)[0] = MERGE (a0, 48, a1, (64-48)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + case 7: + do + { + srcp -= 2 * OPSIZ; + dstp -= 2 * OPSIZ; + + a1 = ((op_t *) srcp)[1]; + a0 = ((op_t *) srcp)[0]; + ((op_t *) dstp)[1] = MERGE (a1, 56, a2, (64-56)); + ((op_t *) dstp)[0] = MERGE (a0, 56, a1, (64-56)); + a2 = a0; + + len -= 2; + } + while (len != 0); + break; + } +} diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/Implies new file mode 100644 index 0000000000..d993b0cd73 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/Implies @@ -0,0 +1,2 @@ +powerpc/powerpc64/power6 +powerpc/powerpc64/power5+ diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/Implies b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/Implies new file mode 100644 index 0000000000..84510991b2 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/Implies @@ -0,0 +1,3 @@ +powerpc/powerpc64/power6/fpu +powerpc/powerpc64/power5+/fpu +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llrint.S b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llrint.S new file mode 100644 index 0000000000..4a4a5f8855 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llrint.S @@ -0,0 +1,45 @@ +/* Round double to long int. POWER6x PowerPC64 version. + Copyright (C) 2006 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. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + + .machine "power6" +/* long long int[r3] __llrint (double x[fp1]) */ +ENTRY (__llrint) + CALL_MCOUNT 0 + fctid fp13,fp1 + mftgpr r3,fp13 + blr + END (__llrint) + +strong_alias (__llrint, __lrint) +weak_alias (__llrint, llrint) +weak_alias (__lrint, lrint) + +#ifdef NO_LONG_DOUBLE +strong_alias (__llrint, __llrintl) +weak_alias (__llrint, llrintl) +strong_alias (__lrint, __lrintl) +weak_alias (__lrint, lrintl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __llrint, llrintl, GLIBC_2_1) +compat_symbol (libm, __lrint, lrintl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llround.S b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llround.S new file mode 100644 index 0000000000..149feaf965 --- /dev/null +++ b/powerpc-cpu/sysdeps/powerpc/powerpc64/power6x/fpu/s_llround.S @@ -0,0 +1,55 @@ +/* llround function. POWER6x PowerPC64 version. + Copyright (C) 2006 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. */ + +#include <sysdep.h> +#include <math_ldbl_opt.h> + +/* long long [r3] llround (float x [fp1]) + IEEE 1003.1 llround function. IEEE specifies "round to the nearest + integer value, rounding halfway cases away from zero, regardless of + the current rounding mode." However PowerPC Architecture defines + "round to Nearest" as "Choose the best approximation. In case of a + tie, choose the one that is even (least significant bit o).". + So we pre-round using the V2.02 Floating Round to Integer Nearest + instruction before we use Floating Convert to Integer Word with + round to zero instruction. */ + + .machine "power6" +ENTRY (__llround) + CALL_MCOUNT 0 + frin fp2,fp1 /* Round to nearest +-0.5. */ + fctidz fp3,fp2 /* Convert To Integer DW round toward 0. */ + mftgpr r3,fp3 /* Transfer integer to R3. */ + blr + END (__llround) + +strong_alias (__llround, __lround) +weak_alias (__llround, llround) +weak_alias (__lround, lround) + +#ifdef NO_LONG_DOUBLE +weak_alias (__llround, llroundl) +strong_alias (__llround, __llroundl) +weak_alias (__lround, lroundl) +strong_alias (__lround, __lroundl) +#endif +#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1) +compat_symbol (libm, __llround, llroundl, GLIBC_2_1) +compat_symbol (libm, __lround, lroundl, GLIBC_2_1) +#endif diff --git a/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc32/power6x/fpu/Implies b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc32/power6x/fpu/Implies new file mode 100644 index 0000000000..ffca13c0ef --- /dev/null +++ b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc32/power6x/fpu/Implies @@ -0,0 +1,4 @@ +# Make sure this comes before the powerpc/powerpc32/fpu that's +# listed in unix/sysv/linux/powerpc/powerpc32/fpu/Implies. +powerpc/powerpc32/power6x/fpu +powerpc/powerpc32/power6/fpu diff --git a/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/970/fpu/Implies b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/970/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/970/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power4/fpu/Implies b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power4/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power4/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power5/fpu/Implies b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power5/fpu/Implies new file mode 100644 index 0000000000..558a5fb174 --- /dev/null +++ b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power5/fpu/Implies @@ -0,0 +1 @@ +powerpc/powerpc64/power4/fpu diff --git a/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power6x/fpu/Implies b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power6x/fpu/Implies new file mode 100644 index 0000000000..fd74b8341c --- /dev/null +++ b/powerpc-cpu/sysdeps/unix/sysv/linux/powerpc/powerpc64/power6x/fpu/Implies @@ -0,0 +1,4 @@ +# Make sure this comes before the powerpc/powerpc64/fpu that's +# listed in unix/sysv/linux/powerpc/powerpc64/fpu/Implies. +powerpc/powerpc64/power6x/fpu +powerpc/powerpc64/power6/fpu |