1 1.1 mrg /* mpz_divexact_gcd -- exact division optimized for GCDs. 2 1.1 mrg 3 1.1 mrg THE FUNCTIONS IN THIS FILE ARE FOR INTERNAL USE AND ARE ALMOST CERTAIN TO 4 1.1 mrg BE SUBJECT TO INCOMPATIBLE CHANGES IN FUTURE GNU MP RELEASES. 5 1.1 mrg 6 1.1.1.2 mrg Copyright 2000, 2005, 2011, 2012 Free Software Foundation, Inc. 7 1.1 mrg 8 1.1 mrg This file is part of the GNU MP Library. 9 1.1 mrg 10 1.1 mrg The GNU MP Library is free software; you can redistribute it and/or modify 11 1.1.1.3 mrg it under the terms of either: 12 1.1.1.3 mrg 13 1.1.1.3 mrg * the GNU Lesser General Public License as published by the Free 14 1.1.1.3 mrg Software Foundation; either version 3 of the License, or (at your 15 1.1.1.3 mrg option) any later version. 16 1.1.1.3 mrg 17 1.1.1.3 mrg or 18 1.1.1.3 mrg 19 1.1.1.3 mrg * the GNU General Public License as published by the Free Software 20 1.1.1.3 mrg Foundation; either version 2 of the License, or (at your option) any 21 1.1.1.3 mrg later version. 22 1.1.1.3 mrg 23 1.1.1.3 mrg or both in parallel, as here. 24 1.1 mrg 25 1.1 mrg The GNU MP Library is distributed in the hope that it will be useful, but 26 1.1 mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 27 1.1.1.3 mrg or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 28 1.1.1.3 mrg for more details. 29 1.1 mrg 30 1.1.1.3 mrg You should have received copies of the GNU General Public License and the 31 1.1.1.3 mrg GNU Lesser General Public License along with the GNU MP Library. If not, 32 1.1.1.3 mrg see https://www.gnu.org/licenses/. */ 33 1.1 mrg 34 1.1 mrg #include "gmp-impl.h" 35 1.1 mrg #include "longlong.h" 36 1.1 mrg 37 1.1 mrg 38 1.1 mrg /* Set q to a/d, expecting d to be from a GCD and therefore usually small. 39 1.1 mrg 40 1.1 mrg The distribution of GCDs of random numbers can be found in Knuth volume 2 41 1.1 mrg section 4.5.2 theorem D. 42 1.1 mrg 43 1.1 mrg GCD chance 44 1.1 mrg 1 60.8% 45 1.1 mrg 2^k 20.2% (1<=k<32) 46 1.1 mrg 3*2^k 9.0% (1<=k<32) 47 1.1 mrg other 10.1% 48 1.1 mrg 49 1.1 mrg Only the low limb is examined for optimizations, since GCDs bigger than 50 1.1 mrg 2^32 (or 2^64) will occur very infrequently. 51 1.1 mrg 52 1.1 mrg Future: This could change to an mpn_divexact_gcd, possibly partly 53 1.1 mrg inlined, if/when the relevant mpq functions change to an mpn based 54 1.1 mrg implementation. */ 55 1.1 mrg 56 1.1 mrg 57 1.1.1.2 mrg #if GMP_NUMB_BITS % 2 == 0 58 1.1 mrg static void 59 1.1 mrg mpz_divexact_by3 (mpz_ptr q, mpz_srcptr a) 60 1.1 mrg { 61 1.1 mrg mp_size_t size = SIZ(a); 62 1.1.1.2 mrg mp_size_t abs_size = ABS(size); 63 1.1.1.2 mrg mp_ptr qp; 64 1.1 mrg 65 1.1.1.2 mrg qp = MPZ_REALLOC (q, abs_size); 66 1.1 mrg 67 1.1.1.2 mrg mpn_bdiv_dbm1 (qp, PTR(a), abs_size, GMP_NUMB_MASK / 3); 68 1.1 mrg 69 1.1.1.2 mrg abs_size -= (qp[abs_size-1] == 0); 70 1.1.1.2 mrg SIZ(q) = (size>0 ? abs_size : -abs_size); 71 1.1.1.2 mrg } 72 1.1.1.2 mrg #endif 73 1.1.1.2 mrg 74 1.1.1.2 mrg #if GMP_NUMB_BITS % 4 == 0 75 1.1.1.2 mrg static void 76 1.1.1.2 mrg mpz_divexact_by5 (mpz_ptr q, mpz_srcptr a) 77 1.1.1.2 mrg { 78 1.1.1.2 mrg mp_size_t size = SIZ(a); 79 1.1.1.2 mrg mp_size_t abs_size = ABS(size); 80 1.1.1.2 mrg mp_ptr qp; 81 1.1.1.2 mrg 82 1.1.1.2 mrg qp = MPZ_REALLOC (q, abs_size); 83 1.1.1.2 mrg 84 1.1.1.2 mrg mpn_bdiv_dbm1 (qp, PTR(a), abs_size, GMP_NUMB_MASK / 5); 85 1.1.1.2 mrg 86 1.1.1.2 mrg abs_size -= (qp[abs_size-1] == 0); 87 1.1.1.2 mrg SIZ(q) = (size>0 ? abs_size : -abs_size); 88 1.1.1.2 mrg } 89 1.1.1.2 mrg #endif 90 1.1.1.2 mrg 91 1.1.1.2 mrg static void 92 1.1.1.2 mrg mpz_divexact_limb (mpz_ptr q, mpz_srcptr a, mp_limb_t d) 93 1.1.1.2 mrg { 94 1.1.1.2 mrg mp_size_t size = SIZ(a); 95 1.1.1.2 mrg mp_size_t abs_size = ABS(size); 96 1.1.1.2 mrg mp_ptr qp; 97 1.1.1.2 mrg 98 1.1.1.2 mrg qp = MPZ_REALLOC (q, abs_size); 99 1.1.1.2 mrg 100 1.1.1.4 mrg MPN_DIVREM_OR_DIVEXACT_1 (qp, PTR(a), abs_size, d); 101 1.1.1.2 mrg 102 1.1.1.2 mrg abs_size -= (qp[abs_size-1] == 0); 103 1.1.1.4 mrg SIZ(q) = (size > 0 ? abs_size : -abs_size); 104 1.1 mrg } 105 1.1 mrg 106 1.1 mrg void 107 1.1 mrg mpz_divexact_gcd (mpz_ptr q, mpz_srcptr a, mpz_srcptr d) 108 1.1 mrg { 109 1.1 mrg ASSERT (mpz_sgn (d) > 0); 110 1.1 mrg 111 1.1.1.2 mrg if (SIZ(a) == 0) 112 1.1.1.2 mrg { 113 1.1.1.2 mrg SIZ(q) = 0; 114 1.1.1.2 mrg return; 115 1.1.1.2 mrg } 116 1.1.1.2 mrg 117 1.1 mrg if (SIZ(d) == 1) 118 1.1 mrg { 119 1.1 mrg mp_limb_t dl = PTR(d)[0]; 120 1.1 mrg int twos; 121 1.1 mrg 122 1.1.1.2 mrg if ((dl & 1) == 0) 123 1.1.1.2 mrg { 124 1.1.1.2 mrg count_trailing_zeros (twos, dl); 125 1.1.1.2 mrg dl >>= twos; 126 1.1.1.2 mrg mpz_tdiv_q_2exp (q, a, twos); 127 1.1.1.2 mrg a = q; 128 1.1.1.2 mrg } 129 1.1 mrg 130 1.1 mrg if (dl == 1) 131 1.1.1.2 mrg { 132 1.1.1.2 mrg if (q != a) 133 1.1.1.2 mrg mpz_set (q, a); 134 1.1.1.2 mrg return; 135 1.1.1.2 mrg } 136 1.1.1.2 mrg #if GMP_NUMB_BITS % 2 == 0 137 1.1 mrg if (dl == 3) 138 1.1.1.2 mrg { 139 1.1.1.2 mrg mpz_divexact_by3 (q, a); 140 1.1.1.2 mrg return; 141 1.1.1.2 mrg } 142 1.1.1.2 mrg #endif 143 1.1.1.2 mrg #if GMP_NUMB_BITS % 4 == 0 144 1.1.1.2 mrg if (dl == 5) 145 1.1.1.2 mrg { 146 1.1.1.2 mrg mpz_divexact_by5 (q, a); 147 1.1.1.2 mrg return; 148 1.1.1.2 mrg } 149 1.1.1.2 mrg #endif 150 1.1.1.2 mrg 151 1.1.1.2 mrg mpz_divexact_limb (q, a, dl); 152 1.1.1.2 mrg return; 153 1.1 mrg } 154 1.1 mrg 155 1.1 mrg mpz_divexact (q, a, d); 156 1.1 mrg } 157
Indexes created Sat Feb 28 05:31:39 UTC 2026