Return to divegcd.c CVS log

Up to [local] / OpenXM_contrib / gmp / mpz

Annotation of OpenXM_contrib/gmp/mpz/divegcd.c, Revision 1.1

1.1     ! ohara       1: /* mpz_divexact_gcd -- exact division optimized for GCDs.
        !             2:
        !             3:    THE FUNCTIONS IN THIS FILE ARE FOR INTERNAL USE AND ARE ALMOST CERTAIN TO
        !             4:    BE SUBJECT TO INCOMPATIBLE CHANGES IN FUTURE GNU MP RELEASES.
        !             5:
        !             6: Copyright 2000 Free Software Foundation, Inc.
        !             7:
        !             8: This file is part of the GNU MP Library.
        !             9:
        !            10: The GNU MP Library is free software; you can redistribute it and/or modify
        !            11: it under the terms of the GNU Lesser General Public License as published by
        !            12: the Free Software Foundation; either version 2.1 of the License, or (at your
        !            13: option) any later version.
        !            14:
        !            15: The GNU MP Library is distributed in the hope that it will be useful, but
        !            16: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
        !            17: or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
        !            18: License for more details.
        !            19:
        !            20: You should have received a copy of the GNU Lesser General Public License
        !            21: along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
        !            22: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
        !            23: MA 02111-1307, USA.  */
        !            24:
        !            25: #include "gmp.h"
        !            26: #include "gmp-impl.h"
        !            27: #include "longlong.h"
        !            28:
        !            29:
        !            30: /* Set q to a/d, expecting d to be from a GCD and therefore usually small.
        !            31:
        !            32:    The distribution of GCDs of random numbers can be found in Knuth volume 2
        !            33:    section 4.5.2 theorem D.
        !            34:
        !            35:             GCD     chance
        !            36:              1       60.7%
        !            37:             2^k      20.2%     (1<=k<32)
        !            38:            3*2^k      9.0%     (1<=k<32)
        !            39:            other     10.1%
        !            40:
        !            41:    Only the low limb is examined for optimizations, since GCDs bigger than
        !            42:    2^32 (or 2^64) will occur very infrequently.
        !            43:
        !            44:    Future: This could change to an mpn_divexact_gcd, possibly partly
        !            45:    inlined, if/when the relevant mpq functions change to an mpn based
        !            46:    implementation.  */
        !            47:
        !            48:
        !            49: static void
        !            50: mpz_divexact_by3 (mpz_ptr q, mpz_srcptr a)
        !            51: {
        !            52:   mp_size_t  size = SIZ(a);
        !            53:   if (size == 0)
        !            54:     {
        !            55:       SIZ(q) = 0;
        !            56:       return;
        !            57:     }
        !            58:   else
        !            59:     {
        !            60:       mp_size_t  abs_size = ABS(size);
        !            61:       mp_ptr     qp;
        !            62:
        !            63:       MPZ_REALLOC (q, abs_size);
        !            64:
        !            65:       qp = PTR(q);
        !            66:       mpn_divexact_by3 (qp, PTR(a), abs_size);
        !            67:
        !            68:       abs_size -= (qp[abs_size-1] == 0);
        !            69:       SIZ(q) = (size>0 ? abs_size : -abs_size);
        !            70:     }
        !            71: }
        !            72:
        !            73: void
        !            74: mpz_divexact_gcd (mpz_ptr q, mpz_srcptr a, mpz_srcptr d)
        !            75: {
        !            76:   ASSERT (mpz_sgn (d) > 0);
        !            77:
        !            78:   if (SIZ(d) == 1)
        !            79:     {
        !            80:       mp_limb_t  dl = PTR(d)[0];
        !            81:       int        twos;
        !            82:
        !            83:       if (dl == 1)
        !            84:         {
        !            85:           if (q != a)
        !            86:             mpz_set (q, a);
        !            87:           return;
        !            88:         }
        !            89:       if (dl == 3)
        !            90:         {
        !            91:           mpz_divexact_by3 (q, a);
        !            92:           return;
        !            93:         }
        !            94:
        !            95:       count_trailing_zeros (twos, dl);
        !            96:       dl >>= twos;
        !            97:
        !            98:       if (dl == 1)
        !            99:         {
        !           100:           mpz_tdiv_q_2exp (q, a, twos);
        !           101:           return;
        !           102:         }
        !           103:       if (dl == 3)
        !           104:         {
        !           105:           mpz_tdiv_q_2exp (q, a, twos);
        !           106:           mpz_divexact_by3 (q, q);
        !           107:           return;
        !           108:         }
        !           109:     }
        !           110:
        !           111:   mpz_divexact (q, a, d);
        !           112: }