[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [Gcl-devel] gcd proposal
From: |
Camm Maguire |
Subject: |
Re: [Gcl-devel] gcd proposal |
Date: |
12 Jan 2004 13:47:29 -0500 |
User-agent: |
Gnus/5.09 (Gnus v5.9.0) Emacs/21.2 |
Greetings, and thanks for your suggestion! Tested and installed.
About a 10% speedup on a pIII, ~16% on a P4/xeon.
Take care,
Andrei Zorine <address@hidden> writes:
> hello,
> i have a small proposal about gcd computation. I propose to use a
> binary gcd algorithm in the get-gcd function (num_arithm.c). It gives
> speed improvement over the existing algorithm. Maxima source tree
> contains a commac.lisp file with gcd redefined and commented out. This
> way it works.
> --
> Andrei Zorine
> object
> get_gcd(object x, object y)
> {
> int i, j, k;
> long k1,t;
> object q, r;
>
> if (number_minusp(x))
> x = number_negate(x);
> if (number_minusp(y))
> y = number_negate(y);
>
> L:
> if (type_of(x) == t_fixnum && type_of(y) == t_fixnum) {
> i = fix(x);
> j = fix(y);
> /* binary gcd starts here */
> k1=0;
> while(i%2==0 && j%2==0)
> {
> k1++;
> i=i>>1;
> j=j>>1;
> };
> if(i%2==1) t=-j;
> else {
> t=i;
> t=t>>1;
> };
> B3:
> while(t%2==0) t=t>>1;
> if(t>0) i=t; else j=-t;
> t=i-j;
> if(t!=0) goto B3;
> return(make_fixnum(i<<k1));
> } /*binary gcd ends here */
>
> if (number_compare(x, y) < 0) {
> r = x;
> x = y;
> y = r;
> }
> if (type_of(y) == t_fixnum && fix(y) == 0) {
> return(x);
> }
> integer_quotient_remainder_1(x, y, &q, &r);
> x = y;
> y = r;
> goto L;
> }
> #include <stdio.h>
> #include <stdlib.h>
> #include <time.h>
> #include <sys/times.h>
>
> long gcl_gcd(long i, long j)
> {
> long k;
> LL:
> if (i < j) {
> k = i;
> i = j;
> j = k;
> }
> if (j == 0) {
> return(i);
> }
> k = i % j;
> i = j;
> j = k;
> goto LL;
> }
>
> long binary_gcd(long u, long v)
> {
> long k,t;
> k=0;
> while(u%2==0 && v%2==0)
> {
> k++;
> u=u/2;
> v=v/2;
> };
> if(u%2==1) t=-v;
> else {
> t=u;
> t=t/2;
> };
> B3:
> while(t%2==0) t=t/2;
> if(t>0) u=t; else v=-t;
> t=u-v;
> if(t!=0) goto B3;
> return u<<k;
> };
>
> int main(int argc, char** argv)
> {
> long i,j,c,cmax=45000,d,modulus=0;
> clock_t start, end;
> struct tms beg, fin;
> double used1, used2;
>
> if(argc>1) modulus = atoi(argv[1]);
>
> srand( time(0) );
> for(;;)
> {
> i = rand();
> j = rand();
> if(modulus>0)
> {
> i=i%modulus;
> j=j%modulus;
> };
> printf("gcd(%d, %d) = 1) %d, 2) %d\t",i,j,gcl_gcd(i,j),
> binary_gcd(i,j));
> start = times(&beg);
> for(c=0;c<cmax;c++)
> d=gcl_gcd(i,j);
> end = times(&fin);
> used1 = ((double) fin.tms_utime - (double) beg.tms_utime);
> start = times(&beg);
> for(c=0;c<cmax;c++)
> d=binary_gcd(i,j);
> end = times(&fin);
> used2 = ((double) fin.tms_utime - (double) beg.tms_utime);
> printf("gcl: %10.8f, binary: %10.8f, gcl-binary = %10.8f\n",
> used1/(double)cmax, used2/(double)cmax, (used1-used2)/(double)cmax);
> }
> return 0;
> }
> _______________________________________________
> Gcl-devel mailing list
> address@hidden
> http://mail.gnu.org/mailman/listinfo/gcl-devel
--
Camm Maguire address@hidden
==========================================================================
"The earth is but one country, and mankind its citizens." -- Baha'u'llah