8000 Merge pull request #8774 from eric-wieser/lcm-gcd · numpy/numpy@d407f24 · GitHub
[go: up one dir, main page]
Skip to content

Commit d407f24

Browse files
mhvkmhvk
authored
Merge pull request #8774 from eric-wieser/lcm-gcd
ENH: Add gcd and lcm ufuncs
2 parents 963c704 + 58998a8 commit d407f24

File tree

10 files changed

+334
-0
lines changed

10 files changed

+334
-0
lines changed

‎doc/release/1.15.0-notes.rst

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,9 @@ Highlights
1010
New functions
1111
=============
1212

13+
* `np.gcd` and `np.lcm`, to compute the greatest common divisor and least
14+
common multiple.
15+
1316

1417
Deprecations
1518
============
@@ -40,6 +43,12 @@ C API changes
4043
New Features
4144
============
4245

46+
``np.gcd`` and ``np.lcm`` ufuncs added for integer and objects types
47+
--------------------------------------------------------------------
48+
These compute the greatest common divisor, and lowest common multiple,
49+
respectively. These work on all the numpy integer types, as well as the
50+
builtin arbitrary-precision `Decimal` and `long` types.
51+
4352

4453
Improvements
4554
============

‎doc/source/reference/routines.math.rst

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -101,6 +101,14 @@ Floating point routines
101101
nextafter
102102
spacing
103103

104+
Rational routines
105+
-----------------
106+
.. autosummary::
107+
:toctree: generated/
108+
109+
lcm
110+
gcd
111+
104112
Arithmetic operations
105113
---------------------
106114
.. autosummary::

‎doc/source/reference/ufuncs.rst

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -550,6 +550,8 @@ Math operations
550550
square
551551
cbrt
552552
reciprocal
553+
gcd
554+
lcm
553555

554556
.. tip::
555557

‎numpy/core/code_generators/generate_umath.py

Lines changed: 14 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -875,6 +875,20 @@ def english_upper(s):
875875
TypeDescription('d', None, 'd', 'di'),
876876
TypeDescription('g', None, 'g', 'gi'),
877877
],
878+
),
879+
'gcd' :
880+
Ufunc(2, 1, Zero,
881+
docstrings.get('numpy.core.umath.gcd'),
882+
"PyUFunc_SimpleBinaryOperationTypeResolver",
883+
TD(ints),
884+
TD('O', f='npy_ObjectGCD'),
885+
),
886+
'lcm' :
887+
Ufunc(2, 1, None,
888+
docstrings.get('numpy.core.umath.lcm'),
889+
"PyUFunc_SimpleBinaryOperationTypeResolver",
890+
TD(ints),
891+
TD('O', f='npy_ObjectLCM'),
878892
)
879893
}
880894

‎numpy/core/code_generators/ufunc_docstrings.py

Lines changed: 60 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3679,3 +3679,63 @@ def add_newdoc(place, name, doc):
36793679
array([ 0., 1., 2., 3., 4., 5.])
36803680
36813681
""")
3682+
3683+
add_newdoc('numpy.core.umath', 'gcd',
3684+
"""
3685+
Returns the grea F438 test common divisor of |x1| and |x2|
3686+
3687+
Parameters
3688+
----------
3689+
x1, x2 : array_like, int
3690+
Arrays of values
3691+
3692+
Returns
3693+
-------
3694+
y : ndarray or scalar
3695+
The greatest common divisor of the absolute value of the inputs
3696+
3697+
See Also
3698+
--------
3699+
lcm : The lowest common multiple
3700+
3701+
Examples
3702+
--------
3703+
>>> np.gcd(12, 20)
3704+
4
3705+
>>> np.gcd.reduce([15, 25, 35])
3706+
5
3707+
>>> np.gcd(np.arange(6), 20)
3708+
array([20, 1, 2, 1, 4, 5])
3709+
3710+
""")
3711+
3712+
add_newdoc('numpy.core.umath', 'lcm',
3713+
"""
3714+
Returns the lowest common multiple of |x1| and |x2|
3715+
3716+
Parameters
3717+
----------
3718+
x1, x2 : array_like, int
3719+
Arrays of values
3720+
3721+
Returns
3722+
-------
3723+
y : ndarray or scalar
3724+
The lowest common multiple of the absolute value of the inputs
3725+
3726+
See Also
3727+
--------
3728+
gcd : The greatest common divisor
3729+
3730+
Examples
3731+
--------
3732+
>>> np.lcm(12, 20)
3733+
60
3734+
>>> np.lcm.reduce([3, 12, 20])
3735+
60
3736+
>>> np.lcm.reduce([40, 12, 20])
3737+
120
3738+
>>> np.lcm(np.arange(6), 20)
3739+
array([ 0, 20, 20, 60, 20, 20])
3740+
3741+
""")

‎numpy/core/src/npymath/npy_math_internal.h.src

Lines changed: 38 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -678,3 +678,41 @@ npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus)
678678
#undef DEG2RAD
679679

680680
/**end repeat**/
681+
682+
/**begin repeat
683+
*
684+
* #type = npy_uint, npy_ulong, npy_ulonglong#
685+
* #c = u,ul,ull#
686+
*/
687+
NPY_INPLACE @type@
688+
npy_gcd@c@(@type@ a, @type@ b)
689+
{
690+
@type@ c;
691+
while (a != 0) {
692+
c = a;
693+
a = b%a;
694+
b = c;
695+
}
696+
return b;
697+
}
698+
699+
NPY_INPLACE @type@
700+
npy_lcm@c@(@type@ a, @type@ b)
701+
{
702+
@type@ gcd = npy_gcd@c@(a, b);
703+
return gcd == 0 ? 0 : a / gcd * b;
704+
}
705+
/**end repeat**/
706+
707+
/**begin repeat
708+
*
709+
* #type = (npy_int, npy_long, npy_longlong)*2#
710+
* #c = (,l,ll)*2#
711+
* #func=gcd*3,lcm*3#
712+
*/
713+
NPY_INPLACE @type@
714+
npy_@func@@c@(@type@ a, @type@ b)
715+
{
716+
return npy_@func@u@c@(a < 0 ? -a : a, b < 0 ? -b : b);
717+
}
718+
/**end repeat**/

‎numpy/core/src/umath/funcs.inc.src

Lines changed: 68 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,7 @@
88

99
#define NPY_NO_DEPRECATED_API NPY_API_VERSION
1010
#include "npy_pycompat.h"
11+
#include "npy_import.h"
1112

1213

1314
/*
@@ -158,6 +159,73 @@ npy_ObjectLogicalNot(PyObject *i1)
158159
}
159160
}
160161

162+
static PyObject *
163+
npy_ObjectGCD(PyObject *i1, PyObject *i2)
164+
{
165+
PyObject *gcd = NULL;
166+
167+
/* use math.gcd if available, and valid on the provided types */
168+
#if PY_VERSION_HEX >= 0x03050000
169+
{
170+
static PyObject *math_gcd_func = NULL;
171+
172+
npy_cache_import("math", "gcd", &math_gcd_func);
173+
if (math_gcd_func == NULL) {
174+
return NULL;
175+
}
176+
gcd = PyObject_CallFunction(math_gcd_func, "OO", i1, i2);
177+
if (gcd != NULL) {
178+
return gcd;
179+
}
180+
/* silence errors, and fall back on pure-python gcd */
181+
PyErr_Clear();
182+
}
183+
#endif
184+
185+
/* otherwise, use our internal one, written in python */
186+
{
187+
static PyObject *internal_gcd_func = NULL;
188+
189+
npy_cache_import("numpy.core._internal", "_gcd", &internal_gcd_func);
190+
if (internal_gcd_func == NULL) {
191+
return NULL;
192+
}
193+
gcd = PyObject_CallFunction(internal_gcd_func, "OO", i1, i2);
194+
if (gcd == NULL) {
195+
return NULL;
196+
}
197+
/* _gcd has some unusual behaviour regarding sign */
198+
return PyNumber_Absolute(gcd);
199+
}
200+
}
201+
202+
static PyObject *
203+
npy_ObjectLCM(PyObject *i1, PyObject *i2)
204+
{
205+
/* lcm(a, b) = abs(a // gcd(a, b) * b) */
206+
207+
PyObject *gcd = npy_ObjectGCD(i1, i2);
208+
PyObject *tmp;
209+
if(gcd == NULL) {
210+
return NULL;
211+
}
212+
/* Floor divide preserves integer types - we know the division will have
213+
* no remainder
214+
*/
215+
tmp = PyNumber_FloorDivide(i1, gcd);
216+
if(tmp == NULL) {
217+
return NULL;
218+
}
219+
220+
tmp = PyNumber_Multiply(tmp, i2);
221+
if(tmp == NULL) {
222+
return NULL;
223+
}
224+
225+
/* even though we fix gcd to be positive, we need to do it again here */
226+
return PyNumber_Absolute(tmp);
227+
}
228+
161229
/*
162230
*****************************************************************************
163231
** COMPLEX FUNCTIONS **

‎numpy/core/src/umath/loops.c.src

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1041,6 +1041,7 @@ NPY_NO_EXPORT void
10411041
/**begin repeat
10421042
* #TYPE = BYTE, SHORT, INT, LONG, LONGLONG#
10431043
* #type = npy_byte, npy_short, npy_int, npy_long, npy_longlong#
1044+
* #c = ,,,l,ll#
10441045
*/
10451046

10461047
NPY_NO_EXPORT NPY_GCC_OPT_3 void
@@ -1132,11 +1133,26 @@ NPY_NO_EXPORT void
11321133
}
11331134
}
11341135

1136+
/**begin repeat1
1137+
* #kind = gcd, lcm#
1138+
**/
1139+
NPY_NO_EXPORT void
1140+
@TYPE@_@kind@(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func))
1141+
{
1142+
BINARY_LOOP {
1143+
const @type@ in1 = *(@type@ *)ip1;
1144+
const @type@ in2 = *(@type@ *)ip2;
1145+
*((@type@ *)op1) = npy_@kind@@c@(in1, in2);
1146+
}
1147+
}
1148+
/**end repeat1**/
1149+
11351150
/**end repeat**/
11361151

11371152
/**begin repeat
11381153
* #TYPE = UBYTE, USHORT, UINT, ULONG, ULONGLONG#
11391154
* #type = npy_ubyte, npy_ushort, npy_uint, npy_ulong, npy_ulonglong#
1155+
* #c = u,u,u,ul,ull#
11401156
*/
11411157

11421158
NPY_NO_EXPORT void
@@ -1204,6 +1220,20 @@ NPY_NO_EXPORT void
12041220
}
12051221
}
12061222

1223+
/**begin repeat1
1224+
* #kind = gcd, lcm#
1225+
**/
1226+
NPY_NO_EXPORT void
1227+
@TYPE@_@kind@(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func))
1228+
{
1229+
BINARY_LOOP {
1230+
const @type@ in1 = *(@type@ *)ip1;
1231+
const @type@ in2 = *(@type@ *)ip2;
1232+
*((@type@ *)op1) = npy_@kind@@c@(in1, in2);
1233+
}
1234+
}
1235+
/**end repeat1**/
1236+
12071237
/**end repeat**/
12081238

12091239
/*

‎numpy/core/src/umath/loops.h.src

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -140,6 +140,12 @@ NPY_NO_EXPORT void
140140
NPY_NO_EXPORT void
141141
@S@@TYPE@_divmod(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func));
142142

143+
NPY_NO_EXPORT void
144+
@S@@TYPE@_gcd(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func));
145+
146+
NPY_NO_EXPORT void
147+
@S@@TYPE@_lcm(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func));
148+
143149
/**end repeat1**/
144150

145151
/**end repeat**/

0 commit comments

Comments
 (0)
0