我正在使用Python v3.1中的分数模块来计算最大公约数。我想知道使用了什么算法。我猜是欧几里得方法,但我想确定一下。文档(http://docs.python.org/py3k/library/fracts.html?highlight=fracts.gcd#fracts.gcd)没有帮助。有人能给我提供线索吗?
quhf5bfb1#
根据the 3.1.2 source code online,gcd的定义如下:
gcd
def gcd(a, b): """Calculate the Greatest Common Divisor of a and b. Unless b==0, the result will have the same sign as b (so that when b is divided by it, the result comes out positive). """ while b: a, b = b, a%b return a
没错,这就是欧几里德算法,用纯Python编写的。
gzszwxb42#
来自分数python“自版本3.5起弃用:请改用math.gcd()。”我也在找算法,希望有帮助。
deikduxw3#
从Python 3.5开始,GCD代码被移到了math.gcd,从Python 3.9开始,math.gcd接受任意数量的参数。实际的GCD代码现在是用C实现的(对于CPython来说),这使得它比原来的纯Python实现要快得多。样板文件:
math.gcd
static PyObject * math_gcd(PyObject *module, PyObject * const *args, Py_ssize_t nargs) { PyObject *res, *x; Py_ssize_t i; if (nargs == 0) { return PyLong_FromLong(0); } res = PyNumber_Index(args[0]); if (res == NULL) { return NULL; } if (nargs == 1) { Py_SETREF(res, PyNumber_Absolute(res)); return res; } PyObject *one = _PyLong_GetOne(); // borrowed ref for (i = 1; i < nargs; i++) { x = _PyNumber_Index(args[i]); if (x == NULL) { Py_DECREF(res); return NULL; } if (res == one) { /* Fast path: just check arguments. It is okay to use identity comparison here. */ Py_DECREF(x); continue; } Py_SETREF(res, _PyLong_GCD(res, x)); Py_DECREF(x); if (res == NULL) { return NULL; } } return res; }
实际计算使用Lehmer's GCD algorithm:
PyObject * _PyLong_GCD(PyObject *aarg, PyObject *barg) { PyLongObject *a, *b, *c = NULL, *d = NULL, *r; stwodigits x, y, q, s, t, c_carry, d_carry; stwodigits A, B, C, D, T; int nbits, k; Py_ssize_t size_a, size_b, alloc_a, alloc_b; digit *a_digit, *b_digit, *c_digit, *d_digit, *a_end, *b_end; a = (PyLongObject *)aarg; b = (PyLongObject *)barg; size_a = Py_SIZE(a); size_b = Py_SIZE(b); if (-2 <= size_a && size_a <= 2 && -2 <= size_b && size_b <= 2) { Py_INCREF(a); Py_INCREF(b); goto simple; } /* Initial reduction: make sure that 0 <= b <= a. */ a = (PyLongObject *)long_abs(a); if (a == NULL) return NULL; b = (PyLongObject *)long_abs(b); if (b == NULL) { Py_DECREF(a); return NULL; } if (long_compare(a, b) < 0) { r = a; a = b; b = r; } /* We now own references to a and b */ alloc_a = Py_SIZE(a); alloc_b = Py_SIZE(b); /* reduce until a fits into 2 digits */ while ((size_a = Py_SIZE(a)) > 2) { nbits = bit_length_digit(a->ob_digit[size_a-1]); /* extract top 2*PyLong_SHIFT bits of a into x, along with corresponding bits of b into y */ size_b = Py_SIZE(b); assert(size_b <= size_a); if (size_b == 0) { if (size_a < alloc_a) { r = (PyLongObject *)_PyLong_Copy(a); Py_DECREF(a); } else r = a; Py_DECREF(b); Py_XDECREF(c); Py_XDECREF(d); return (PyObject *)r; } x = (((twodigits)a->ob_digit[size_a-1] << (2*PyLong_SHIFT-nbits)) | ((twodigits)a->ob_digit[size_a-2] << (PyLong_SHIFT-nbits)) | (a->ob_digit[size_a-3] >> nbits)); y = ((size_b >= size_a - 2 ? b->ob_digit[size_a-3] >> nbits : 0) | (size_b >= size_a - 1 ? (twodigits)b->ob_digit[size_a-2] << (PyLong_SHIFT-nbits) : 0) | (size_b >= size_a ? (twodigits)b->ob_digit[size_a-1] << (2*PyLong_SHIFT-nbits) : 0)); /* inner loop of Lehmer's algorithm; A, B, C, D never grow larger than PyLong_MASK during the algorithm. */ A = 1; B = 0; C = 0; D = 1; for (k=0;; k++) { if (y-C == 0) break; q = (x+(A-1))/(y-C); s = B+q*D; t = x-q*y; if (s > t) break; x = y; y = t; t = A+q*C; A = D; B = C; C = s; D = t; } if (k == 0) { /* no progress; do a Euclidean step */ if (l_mod(a, b, &r) < 0) goto error; Py_SETREF(a, b); b = r; alloc_a = alloc_b; alloc_b = Py_SIZE(b); continue; } /* a, b = A*b-B*a, D*a-C*b if k is odd a, b = A*a-B*b, D*b-C*a if k is even */ if (k&1) { T = -A; A = -B; B = T; T = -C; C = -D; D = T; } if (c != NULL) { Py_SET_SIZE(c, size_a); } else if (Py_REFCNT(a) == 1) { c = (PyLongObject*)Py_NewRef(a); } else { alloc_a = size_a; c = _PyLong_New(size_a); if (c == NULL) goto error; } if (d != NULL) { Py_SET_SIZE(d, size_a); } else if (Py_REFCNT(b) == 1 && size_a <= alloc_b) { d = (PyLongObject*)Py_NewRef(b); Py_SET_SIZE(d, size_a); } else { alloc_b = size_a; d = _PyLong_New(size_a); if (d == NULL) goto error; } a_end = a->ob_digit + size_a; b_end = b->ob_digit + size_b; /* compute new a and new b in parallel */ a_digit = a->ob_digit; b_digit = b->ob_digit; c_digit = c->ob_digit; d_digit = d->ob_digit; c_carry = 0; d_carry = 0; while (b_digit < b_end) { c_carry += (A * *a_digit) - (B * *b_digit); d_carry += (D * *b_digit++) - (C * *a_digit++); *c_digit++ = (digit)(c_carry & PyLong_MASK); *d_digit++ = (digit)(d_carry & PyLong_MASK); c_carry >>= PyLong_SHIFT; d_carry >>= PyLong_SHIFT; } while (a_digit < a_end) { c_carry += A * *a_digit; d_carry -= C * *a_digit++; *c_digit++ = (digit)(c_carry & PyLong_MASK); *d_digit++ = (digit)(d_carry & PyLong_MASK); c_carry >>= PyLong_SHIFT; d_carry >>= PyLong_SHIFT; } assert(c_carry == 0); assert(d_carry == 0); Py_INCREF(c); Py_INCREF(d); Py_DECREF(a); Py_DECREF(b); a = long_normalize(c); b = long_normalize(d); } Py_XDECREF(c); Py_XDECREF(d); simple: assert(Py_REFCNT(a) > 0); assert(Py_REFCNT(b) > 0); /* Issue #24999: use two shifts instead of ">> 2*PyLong_SHIFT" to avoid undefined behaviour when LONG_MAX type is smaller than 60 bits */ #if LONG_MAX >> PyLong_SHIFT >> PyLong_SHIFT /* a fits into a long, so b must too */ x = PyLong_AsLong((PyObject *)a); y = PyLong_AsLong((PyObject *)b); #elif LLONG_MAX >> PyLong_SHIFT >> PyLong_SHIFT x = PyLong_AsLongLong((PyObject *)a); y = PyLong_AsLongLong((PyObject *)b); #else # error "_PyLong_GCD" #endif x = Py_ABS(x); y = Py_ABS(y); Py_DECREF(a); Py_DECREF(b); /* usual Euclidean algorithm for longs */ while (y != 0) { t = y; y = x % y; x = t; } #if LONG_MAX >> PyLong_SHIFT >> PyLong_SHIFT return PyLong_FromLong(x); #elif LLONG_MAX >> PyLong_SHIFT >> PyLong_SHIFT return PyLong_FromLongLong(x); #else # error "_PyLong_GCD" #endif error: Py_DECREF(a); Py_DECREF(b); Py_XDECREF(c); Py_XDECREF(d); return NULL; }
3条答案
按热度按时间quhf5bfb1#
根据the 3.1.2 source code online,
gcd
的定义如下:没错,这就是欧几里德算法,用纯Python编写的。
gzszwxb42#
来自分数python
“自版本3.5起弃用:请改用math.gcd()。”
我也在找算法,希望有帮助。
deikduxw3#
从Python 3.5开始,GCD代码被移到了
math.gcd
,从Python 3.9开始,math.gcd
接受任意数量的参数。实际的GCD代码现在是用C实现的(对于CPython来说),这使得它比原来的纯Python实现要快得多。
样板文件:
实际计算使用Lehmer's GCD algorithm: