1
#define PY_SSIZE_T_CLEAN
2
#include <Python.h>
3
#include "structmember.h"
4

5
#define NPY_NO_DEPRECATED_API NPY_API_VERSION
6
#define _MULTIARRAYMODULE
7
#include "numpy/arrayobject.h"
8
#include "numpy/arrayscalars.h"
9

10
#include "numpy/npy_math.h"
11

12
#include "npy_config.h"
13

14
#include "npy_pycompat.h"
15

16
#include "ctors.h"
17

18
#include "shape.h"
19

20
#include "multiarraymodule.h" /* for interned strings */
21
#include "templ_common.h" /* for npy_mul_with_overflow_intp */
22
#include "common.h" /* for convert_shape_to_string */
23
#include "alloc.h"
24

25
static int
26
_fix_unknown_dimension(PyArray_Dims *newshape, PyArrayObject *arr);
27

28
static int
29
_attempt_nocopy_reshape(PyArrayObject *self, int newnd, const npy_intp *newdims,
30
                        npy_intp *newstrides, int is_f_order);
31

32
static void
33
_putzero(char *optr, PyObject *zero, PyArray_Descr *dtype);
34

35
/*NUMPY_API
36
 * Resize (reallocate data).  Only works if nothing else is referencing this
37
 * array and it is contiguous.  If refcheck is 0, then the reference count is
38
 * not checked and assumed to be 1.  You still must own this data and have no
39
 * weak-references and no base object.
40
 */
41
NPY_NO_EXPORT PyObject *
42 1
PyArray_Resize(PyArrayObject *self, PyArray_Dims *newshape, int refcheck,
43
               NPY_ORDER NPY_UNUSED(order))
44
{
45
    npy_intp oldnbytes, newnbytes;
46
    npy_intp oldsize, newsize;
47 1
    int new_nd=newshape->len, k, elsize;
48
    int refcnt;
49 1
    npy_intp* new_dimensions=newshape->ptr;
50
    npy_intp new_strides[NPY_MAXDIMS];
51
    npy_intp *dimptr;
52
    char *new_data;
53

54 1
    if (!PyArray_ISONESEGMENT(self)) {
55 0
        PyErr_SetString(PyExc_ValueError,
56
                "resize only works on single-segment arrays");
57 0
        return NULL;
58
    }
59

60
    /* Compute total size of old and new arrays. The new size might overflow */
61 1
    oldsize = PyArray_SIZE(self);
62 1
    newsize = 1;
63 1
    for(k = 0; k < new_nd; k++) {
64 1
        if (new_dimensions[k] == 0) {
65
            newsize = 0;
66
            break;
67
        }
68 1
        if (new_dimensions[k] < 0) {
69 1
            PyErr_SetString(PyExc_ValueError,
70
                    "negative dimensions not allowed");
71 1
            return NULL;
72
        }
73 1
        if (npy_mul_with_overflow_intp(&newsize, newsize, new_dimensions[k])) {
74 0
            return PyErr_NoMemory();
75
        }
76
    }
77

78
    /* Convert to number of bytes. The new count might overflow */
79 1
    elsize = PyArray_DESCR(self)->elsize;
80 1
    oldnbytes = oldsize * elsize;
81 1
    if (npy_mul_with_overflow_intp(&newnbytes, newsize, elsize)) {
82 0
        return PyErr_NoMemory();
83
    }
84

85 1
    if (oldnbytes != newnbytes) {
86 1
        if (!(PyArray_FLAGS(self) & NPY_ARRAY_OWNDATA)) {
87 0
            PyErr_SetString(PyExc_ValueError,
88
                    "cannot resize this array: it does not own its data");
89 0
            return NULL;
90
        }
91

92 1
        if (PyArray_BASE(self) != NULL
93 1
              || (((PyArrayObject_fields *)self)->weakreflist != NULL)) {
94 1
            PyErr_SetString(PyExc_ValueError,
95
                    "cannot resize an array that "
96
                    "references or is referenced\n"
97
                    "by another array in this way. Use the np.resize function.");
98 1
            return NULL;
99
        }
100 1
        if (refcheck) {
101
#ifdef PYPY_VERSION
102
            PyErr_SetString(PyExc_ValueError,
103
                    "cannot resize an array with refcheck=True on PyPy.\n"
104
                    "Use the np.resize function or refcheck=False");
105
            return NULL;
106
#else
107 1
            refcnt = PyArray_REFCOUNT(self);
108
#endif /* PYPY_VERSION */
109
        }
110
        else {
111
            refcnt = 1;
112
        }
113 1
        if (refcnt > 2) {
114 1
            PyErr_SetString(PyExc_ValueError,
115
                    "cannot resize an array that "
116
                    "references or is referenced\n"
117
                    "by another array in this way.\n"
118
                    "Use the np.resize function or refcheck=False");
119 1
            return NULL;
120
        }
121

122
        /* Reallocate space if needed - allocating 0 is forbidden */
123 1
        new_data = PyDataMem_RENEW(
124
            PyArray_DATA(self), newnbytes == 0 ? elsize : newnbytes);
125 1
        if (new_data == NULL) {
126 0
            PyErr_SetString(PyExc_MemoryError,
127
                    "cannot allocate memory for array");
128 0
            return NULL;
129
        }
130 1
        ((PyArrayObject_fields *)self)->data = new_data;
131
    }
132

133 1
    if (newnbytes > oldnbytes && PyArray_ISWRITEABLE(self)) {
134
        /* Fill new memory with zeros */
135 1
        if (PyDataType_FLAGCHK(PyArray_DESCR(self), NPY_ITEM_REFCOUNT)) {
136 1
            PyObject *zero = PyLong_FromLong(0);
137
            char *optr;
138 1
            optr = PyArray_BYTES(self) + oldnbytes;
139 1
            npy_intp n_new = newsize - oldsize;
140 1
            for (npy_intp i = 0; i < n_new; i++) {
141 1
                _putzero((char *)optr, zero, PyArray_DESCR(self));
142 1
                optr += elsize;
143
            }
144 1
            Py_DECREF(zero);
145
        }
146
        else{
147 1
            memset(PyArray_BYTES(self) + oldnbytes, 0, newnbytes - oldnbytes);
148
        }
149
    }
150

151 1
    if (new_nd > 0) {
152 1
        if (PyArray_NDIM(self) != new_nd) {
153
            /* Different number of dimensions. */
154 1
            ((PyArrayObject_fields *)self)->nd = new_nd;
155
            /* Need new dimensions and strides arrays */
156 1
            dimptr = PyDimMem_RENEW(PyArray_DIMS(self), 3*new_nd);
157 1
            if (dimptr == NULL) {
158 0
                PyErr_SetString(PyExc_MemoryError,
159
                                "cannot allocate memory for array");
160 0
                return NULL;
161
            }
162 1
            ((PyArrayObject_fields *)self)->dimensions = dimptr;
163 1
            ((PyArrayObject_fields *)self)->strides = dimptr + new_nd;
164
        }
165
        /* make new_strides variable */
166 1
        _array_fill_strides(new_strides, new_dimensions, new_nd,
167 1
                            PyArray_DESCR(self)->elsize, PyArray_FLAGS(self),
168
                            &(((PyArrayObject_fields *)self)->flags));
169 1
        memmove(PyArray_DIMS(self), new_dimensions, new_nd*sizeof(npy_intp));
170 1
        memmove(PyArray_STRIDES(self), new_strides, new_nd*sizeof(npy_intp));
171
    }
172
    else {
173 1
        PyDimMem_FREE(((PyArrayObject_fields *)self)->dimensions);
174 1
        ((PyArrayObject_fields *)self)->nd = 0;
175 1
        ((PyArrayObject_fields *)self)->dimensions = NULL;
176 1
        ((PyArrayObject_fields *)self)->strides = NULL;
177
    }
178 1
    Py_RETURN_NONE;
179
}
180

181
/*
182
 * Returns a new array
183
 * with the new shape from the data
184
 * in the old array --- order-perspective depends on order argument.
185
 * copy-only-if-necessary
186
 */
187

188
/*NUMPY_API
189
 * New shape for an array
190
 */
191
NPY_NO_EXPORT PyObject *
192 1
PyArray_Newshape(PyArrayObject *self, PyArray_Dims *newdims,
193
                 NPY_ORDER order)
194
{
195
    npy_intp i;
196 1
    npy_intp *dimensions = newdims->ptr;
197
    PyArrayObject *ret;
198 1
    int ndim = newdims->len;
199
    npy_bool same;
200 1
    npy_intp *strides = NULL;
201
    npy_intp newstrides[NPY_MAXDIMS];
202
    int flags;
203

204 1
    if (order == NPY_ANYORDER) {
205 1
        order = PyArray_ISFORTRAN(self);
206
    }
207 1
    else if (order == NPY_KEEPORDER) {
208 1
        PyErr_SetString(PyExc_ValueError,
209
                "order 'K' is not permitted for reshaping");
210 1
        return NULL;
211
    }
212
    /*  Quick check to make sure anything actually needs to be done */
213 1
    if (ndim == PyArray_NDIM(self)) {
214
        same = NPY_TRUE;
215
        i = 0;
216 1
        while (same && i < ndim) {
217 1
            if (PyArray_DIM(self,i) != dimensions[i]) {
218 1
                same=NPY_FALSE;
219
            }
220 1
            i++;
221
        }
222 1
        if (same) {
223 1
            return PyArray_View(self, NULL, NULL);
224
        }
225
    }
226

227
    /*
228
     * fix any -1 dimensions and check new-dimensions against old size
229
     */
230 1
    if (_fix_unknown_dimension(newdims, self) < 0) {
231
        return NULL;
232
    }
233
    /*
234
     * sometimes we have to create a new copy of the array
235
     * in order to get the right orientation and
236
     * because we can't just re-use the buffer with the
237
     * data in the order it is in.
238
     * NPY_RELAXED_STRIDES_CHECKING: size check is unnecessary when set.
239
     */
240 1
    Py_INCREF(self);
241 1
    if ((PyArray_SIZE(self) > 1) &&
242 1
        ((order == NPY_CORDER && !PyArray_IS_C_CONTIGUOUS(self)) ||
243 1
         (order == NPY_FORTRANORDER && !PyArray_IS_F_CONTIGUOUS(self)))) {
244 1
        int success = 0;
245 1
        success = _attempt_nocopy_reshape(self, ndim, dimensions,
246
                                          newstrides, order);
247 1
        if (success) {
248
            /* no need to copy the array after all */
249
            strides = newstrides;
250
        }
251
        else {
252
            PyObject *newcopy;
253 1
            newcopy = PyArray_NewCopy(self, order);
254 1
            Py_DECREF(self);
255 1
            if (newcopy == NULL) {
256
                return NULL;
257
            }
258
            self = (PyArrayObject *)newcopy;
259
        }
260
    }
261
    /* We always have to interpret the contiguous buffer correctly */
262

263
    /* Make sure the flags argument is set. */
264 1
    flags = PyArray_FLAGS(self);
265 1
    if (ndim > 1) {
266 1
        if (order == NPY_FORTRANORDER) {
267 1
            flags &= ~NPY_ARRAY_C_CONTIGUOUS;
268 1
            flags |= NPY_ARRAY_F_CONTIGUOUS;
269
        }
270
        else {
271 1
            flags &= ~NPY_ARRAY_F_CONTIGUOUS;
272 1
            flags |= NPY_ARRAY_C_CONTIGUOUS;
273
        }
274
    }
275

276 1
    Py_INCREF(PyArray_DESCR(self));
277 1
    ret = (PyArrayObject *)PyArray_NewFromDescr_int(
278 1
            Py_TYPE(self), PyArray_DESCR(self),
279
            ndim, dimensions, strides, PyArray_DATA(self),
280
            flags, (PyObject *)self, (PyObject *)self,
281
            0, 1);
282 1
    Py_DECREF(self);
283
    return (PyObject *)ret;
284
}
285

286

287

288
/* For backward compatibility -- Not recommended */
289

290
/*NUMPY_API
291
 * Reshape
292
 */
293
NPY_NO_EXPORT PyObject *
294 1
PyArray_Reshape(PyArrayObject *self, PyObject *shape)
295
{
296
    PyObject *ret;
297
    PyArray_Dims newdims;
298

299 1
    if (!PyArray_IntpConverter(shape, &newdims)) {
300
        return NULL;
301
    }
302 1
    ret = PyArray_Newshape(self, &newdims, NPY_CORDER);
303 1
    npy_free_cache_dim_obj(newdims);
304 1
    return ret;
305
}
306

307

308
static void
309 1
_putzero(char *optr, PyObject *zero, PyArray_Descr *dtype)
310
{
311 1
    if (!PyDataType_FLAGCHK(dtype, NPY_ITEM_REFCOUNT)) {
312 0
        memset(optr, 0, dtype->elsize);
313
    }
314 1
    else if (PyDataType_HASFIELDS(dtype)) {
315 1
        PyObject *key, *value, *title = NULL;
316
        PyArray_Descr *new;
317
        int offset;
318 1
        Py_ssize_t pos = 0;
319 1
        while (PyDict_Next(dtype->fields, &pos, &key, &value)) {
320 1
            if (NPY_TITLE_KEY(key, value)) {
321 0
                continue;
322
            }
323 1
            if (!PyArg_ParseTuple(value, "Oi|O", &new, &offset, &title)) {
324 0
                return;
325
            }
326 1
            _putzero(optr + offset, zero, new);
327
        }
328
    }
329
    else {
330
        npy_intp i;
331 1
        npy_intp nsize = dtype->elsize / sizeof(zero);
332

333 1
        for (i = 0; i < nsize; i++) {
334 1
            Py_INCREF(zero);
335 1
            memcpy(optr, &zero, sizeof(zero));
336 1
            optr += sizeof(zero);
337
        }
338
    }
339
    return;
340
}
341

342

343
/*
344
 * attempt to reshape an array without copying data
345
 *
346
 * The requested newdims are not checked, but must be compatible with
347
 * the size of self, which must be non-zero. Other than that this
348
 * function should correctly handle all reshapes, including axes of
349
 * length 1. Zero strides should work but are untested.
350
 *
351
 * If a copy is needed, returns 0
352
 * If no copy is needed, returns 1 and fills newstrides
353
 *     with appropriate strides
354
 *
355
 * The "is_f_order" argument describes how the array should be viewed
356
 * during the reshape, not how it is stored in memory (that
357
 * information is in PyArray_STRIDES(self)).
358
 *
359
 * If some output dimensions have length 1, the strides assigned to
360
 * them are arbitrary. In the current implementation, they are the
361
 * stride of the next-fastest index.
362
 */
363
static int
364 1
_attempt_nocopy_reshape(PyArrayObject *self, int newnd, const npy_intp *newdims,
365
                        npy_intp *newstrides, int is_f_order)
366
{
367
    int oldnd;
368
    npy_intp olddims[NPY_MAXDIMS];
369
    npy_intp oldstrides[NPY_MAXDIMS];
370
    npy_intp last_stride;
371
    int oi, oj, ok, ni, nj, nk;
372

373 1
    oldnd = 0;
374
    /*
375
     * Remove axes with dimension 1 from the old array. They have no effect
376
     * but would need special cases since their strides do not matter.
377
     */
378 1
    for (oi = 0; oi < PyArray_NDIM(self); oi++) {
379 1
        if (PyArray_DIMS(self)[oi]!= 1) {
380 1
            olddims[oldnd] = PyArray_DIMS(self)[oi];
381 1
            oldstrides[oldnd] = PyArray_STRIDES(self)[oi];
382 1
            oldnd++;
383
        }
384
    }
385

386
    /* oi to oj and ni to nj give the axis ranges currently worked with */
387
    oi = 0;
388
    oj = 1;
389
    ni = 0;
390
    nj = 1;
391 1
    while (ni < newnd && oi < oldnd) {
392 1
        npy_intp np = newdims[ni];
393 1
        npy_intp op = olddims[oi];
394

395 1
        while (np != op) {
396 1
            if (np < op) {
397
                /* Misses trailing 1s, these are handled later */
398 1
                np *= newdims[nj++];
399
            } else {
400 1
                op *= olddims[oj++];
401
            }
402
        }
403

404
        /* Check whether the original axes can be combined */
405 1
        for (ok = oi; ok < oj - 1; ok++) {
406 1
            if (is_f_order) {
407 1
                if (oldstrides[ok+1] != olddims[ok]*oldstrides[ok]) {
408
                     /* not contiguous enough */
409
                    return 0;
410
                }
411
            }
412
            else {
413
                /* C order */
414 1
                if (oldstrides[ok] != olddims[ok+1]*oldstrides[ok+1]) {
415
                    /* not contiguous enough */
416
                    return 0;
417
                }
418
            }
419
        }
420

421
        /* Calculate new strides for all axes currently worked with */
422 1
        if (is_f_order) {
423 1
            newstrides[ni] = oldstrides[oi];
424 1
            for (nk = ni + 1; nk < nj; nk++) {
425 1
                newstrides[nk] = newstrides[nk - 1]*newdims[nk - 1];
426
            }
427
        }
428
        else {
429
            /* C order */
430 1
            newstrides[nj - 1] = oldstrides[oj - 1];
431 1
            for (nk = nj - 1; nk > ni; nk--) {
432 1
                newstrides[nk - 1] = newstrides[nk]*newdims[nk];
433
            }
434
        }
435 1
        ni = nj++;
436 1
        oi = oj++;
437
    }
438

439
    /*
440
     * Set strides corresponding to trailing 1s of the new shape.
441
     */
442 1
    if (ni >= 1) {
443 1
        last_stride = newstrides[ni - 1];
444
    }
445
    else {
446 0
        last_stride = PyArray_ITEMSIZE(self);
447
    }
448 1
    if (is_f_order) {
449 1
        last_stride *= newdims[ni - 1];
450
    }
451 1
    for (nk = ni; nk < newnd; nk++) {
452 1
        newstrides[nk] = last_stride;
453
    }
454

455
    return 1;
456
}
457

458
static void
459 1
raise_reshape_size_mismatch(PyArray_Dims *newshape, PyArrayObject *arr)
460
{
461 1
    PyObject *msg = PyUnicode_FromFormat("cannot reshape array of size %zd "
462 1
                                         "into shape ", PyArray_SIZE(arr));
463 1
    PyObject *tmp = convert_shape_to_string(newshape->len, newshape->ptr, "");
464

465 1
    PyUString_ConcatAndDel(&msg, tmp);
466 1
    if (msg != NULL) {
467 1
        PyErr_SetObject(PyExc_ValueError, msg);
468 1
        Py_DECREF(msg);
469
    }
470
}
471

472
static int
473 1
_fix_unknown_dimension(PyArray_Dims *newshape, PyArrayObject *arr)
474
{
475
    npy_intp *dimensions;
476 1
    npy_intp s_original = PyArray_SIZE(arr);
477
    npy_intp i_unknown, s_known;
478
    int i, n;
479

480 1
    dimensions = newshape->ptr;
481 1
    n = newshape->len;
482 1
    s_known = 1;
483 1
    i_unknown = -1;
484

485 1
    for (i = 0; i < n; i++) {
486 1
        if (dimensions[i] < 0) {
487 1
            if (i_unknown == -1) {
488 1
                i_unknown = i;
489
            }
490
            else {
491 0
                PyErr_SetString(PyExc_ValueError,
492
                                "can only specify one unknown dimension");
493 0
                return -1;
494
            }
495
        }
496 1
        else if (npy_mul_with_overflow_intp(&s_known, s_known,
497
                                            dimensions[i])) {
498 1
            raise_reshape_size_mismatch(newshape, arr);
499 1
            return -1;
500
        }
501
    }
502

503 1
    if (i_unknown >= 0) {
504 1
        if (s_known == 0 || s_original % s_known != 0) {
505 0
            raise_reshape_size_mismatch(newshape, arr);
506 0
            return -1;
507
        }
508 1
        dimensions[i_unknown] = s_original / s_known;
509
    }
510
    else {
511 1
        if (s_original != s_known) {
512 1
            raise_reshape_size_mismatch(newshape, arr);
513 1
            return -1;
514
        }
515
    }
516
    return 0;
517
}
518

519
/*NUMPY_API
520
 *
521
 * return a new view of the array object with all of its unit-length
522
 * dimensions squeezed out if needed, otherwise
523
 * return the same array.
524
 */
525
NPY_NO_EXPORT PyObject *
526 1
PyArray_Squeeze(PyArrayObject *self)
527
{
528
    PyArrayObject *ret;
529
    npy_bool unit_dims[NPY_MAXDIMS];
530
    int idim, ndim, any_ones;
531
    npy_intp *shape;
532

533 1
    ndim = PyArray_NDIM(self);
534 1
    shape = PyArray_SHAPE(self);
535

536 1
    any_ones = 0;
537 1
    for (idim = 0; idim < ndim; ++idim) {
538 1
        if (shape[idim] == 1) {
539 1
            unit_dims[idim] = 1;
540 1
            any_ones = 1;
541
        }
542
        else {
543 1
            unit_dims[idim] = 0;
544
        }
545
    }
546

547
    /* If there were no ones to squeeze out, return the same array */
548 1
    if (!any_ones) {
549 1
        Py_INCREF(self);
550 1
        return (PyObject *)self;
551
    }
552

553 1
    ret = (PyArrayObject *)PyArray_View(self, NULL, &PyArray_Type);
554 1
    if (ret == NULL) {
555
        return NULL;
556
    }
557

558 1
    PyArray_RemoveAxesInPlace(ret, unit_dims);
559

560
    /*
561
     * If self isn't not a base class ndarray, call its
562
     * __array_wrap__ method
563
     */
564 1
    if (Py_TYPE(self) != &PyArray_Type) {
565 1
        PyArrayObject *tmp = PyArray_SubclassWrap(self, ret);
566 1
        Py_DECREF(ret);
567
        ret = tmp;
568
    }
569

570
    return (PyObject *)ret;
571
}
572

573
/*
574
 * Just like PyArray_Squeeze, but allows the caller to select
575
 * a subset of the size-one dimensions to squeeze out.
576
 */
577
NPY_NO_EXPORT PyObject *
578 1
PyArray_SqueezeSelected(PyArrayObject *self, npy_bool *axis_flags)
579
{
580
    PyArrayObject *ret;
581
    int idim, ndim, any_ones;
582
    npy_intp *shape;
583

584 1
    ndim = PyArray_NDIM(self);
585 1
    shape = PyArray_SHAPE(self);
586

587
    /* Verify that the axes requested are all of size one */
588 1
    any_ones = 0;
589 1
    for (idim = 0; idim < ndim; ++idim) {
590 1
        if (axis_flags[idim] != 0) {
591 1
            if (shape[idim] == 1) {
592
                any_ones = 1;
593
            }
594
            else {
595 1
                PyErr_SetString(PyExc_ValueError,
596
                        "cannot select an axis to squeeze out "
597
                        "which has size not equal to one");
598 1
                return NULL;
599
            }
600
        }
601
    }
602

603
    /* If there were no axes to squeeze out, return the same array */
604 1
    if (!any_ones) {
605 1
        Py_INCREF(self);
606 1
        return (PyObject *)self;
607
    }
608

609 1
    ret = (PyArrayObject *)PyArray_View(self, NULL, &PyArray_Type);
610 1
    if (ret == NULL) {
611
        return NULL;
612
    }
613

614 1
    PyArray_RemoveAxesInPlace(ret, axis_flags);
615

616
    /*
617
     * If self isn't not a base class ndarray, call its
618
     * __array_wrap__ method
619
     */
620 1
    if (Py_TYPE(self) != &PyArray_Type) {
621 1
        PyArrayObject *tmp = PyArray_SubclassWrap(self, ret);
622 1
        Py_DECREF(ret);
623
        ret = tmp;
624
    }
625

626
    return (PyObject *)ret;
627
}
628

629
/*NUMPY_API
630
 * SwapAxes
631
 */
632
NPY_NO_EXPORT PyObject *
633 1
PyArray_SwapAxes(PyArrayObject *ap, int a1, int a2)
634
{
635
    PyArray_Dims new_axes;
636
    npy_intp dims[NPY_MAXDIMS];
637 1
    int n = PyArray_NDIM(ap);
638
    int i;
639

640 1
    if (check_and_adjust_axis_msg(&a1, n, npy_ma_str_axis1) < 0) {
641
        return NULL;
642
    }
643 1
    if (check_and_adjust_axis_msg(&a2, n, npy_ma_str_axis2) < 0) {
644
        return NULL;
645
    }
646

647 1
    for (i = 0; i < n; ++i) {
648 1
        dims[i] = i;
649
    }
650 1
    dims[a1] = a2;
651 1
    dims[a2] = a1;
652

653 1
    new_axes.ptr = dims;
654 1
    new_axes.len = n;
655

656 1
    return PyArray_Transpose(ap, &new_axes);
657
}
658

659

660
/*NUMPY_API
661
 * Return Transpose.
662
 */
663
NPY_NO_EXPORT PyObject *
664 1
PyArray_Transpose(PyArrayObject *ap, PyArray_Dims *permute)
665
{
666
    npy_intp *axes;
667
    int i, n;
668
    int permutation[NPY_MAXDIMS], reverse_permutation[NPY_MAXDIMS];
669 1
    PyArrayObject *ret = NULL;
670
    int flags;
671

672 1
    if (permute == NULL) {
673 1
        n = PyArray_NDIM(ap);
674 1
        for (i = 0; i < n; i++) {
675 1
            permutation[i] = n-1-i;
676
        }
677
    }
678
    else {
679 1
        n = permute->len;
680 1
        axes = permute->ptr;
681 1
        if (n != PyArray_NDIM(ap)) {
682 1
            PyErr_SetString(PyExc_ValueError,
683
                            "axes don't match array");
684 1
            return NULL;
685
        }
686 1
        for (i = 0; i < n; i++) {
687 1
            reverse_permutation[i] = -1;
688
        }
689 1
        for (i = 0; i < n; i++) {
690 1
            int axis = axes[i];
691 1
            if (check_and_adjust_axis(&axis, PyArray_NDIM(ap)) < 0) {
692 1
                return NULL;
693
            }
694 1
            if (reverse_permutation[axis] != -1) {
695 1
                PyErr_SetString(PyExc_ValueError,
696
                                "repeated axis in transpose");
697 1
                return NULL;
698
            }
699 1
            reverse_permutation[axis] = i;
700 1
            permutation[i] = axis;
701
        }
702
    }
703

704 1
    flags = PyArray_FLAGS(ap);
705

706
    /*
707
     * this allocates memory for dimensions and strides (but fills them
708
     * incorrectly), sets up descr, and points data at PyArray_DATA(ap).
709
     */
710 1
    Py_INCREF(PyArray_DESCR(ap));
711 1
    ret = (PyArrayObject *) PyArray_NewFromDescrAndBase(
712 1
            Py_TYPE(ap), PyArray_DESCR(ap),
713 1
            n, PyArray_DIMS(ap), NULL, PyArray_DATA(ap),
714
            flags, (PyObject *)ap, (PyObject *)ap);
715 1
    if (ret == NULL) {
716
        return NULL;
717
    }
718

719
    /* fix the dimensions and strides of the return-array */
720 1
    for (i = 0; i < n; i++) {
721 1
        PyArray_DIMS(ret)[i] = PyArray_DIMS(ap)[permutation[i]];
722 1
        PyArray_STRIDES(ret)[i] = PyArray_STRIDES(ap)[permutation[i]];
723
    }
724 1
    PyArray_UpdateFlags(ret, NPY_ARRAY_C_CONTIGUOUS | NPY_ARRAY_F_CONTIGUOUS |
725
                        NPY_ARRAY_ALIGNED);
726 1
    return (PyObject *)ret;
727
}
728

729
/*
730
 * Sorts items so stride is descending, because C-order
731
 * is the default in the face of ambiguity.
732
 */
733 1
static int _npy_stride_sort_item_comparator(const void *a, const void *b)
734
{
735 1
    npy_intp astride = ((const npy_stride_sort_item *)a)->stride,
736 1
            bstride = ((const npy_stride_sort_item *)b)->stride;
737

738
    /* Sort the absolute value of the strides */
739 1
    if (astride < 0) {
740 1
        astride = -astride;
741
    }
742 1
    if (bstride < 0) {
743 1
        bstride = -bstride;
744
    }
745

746 1
    if (astride == bstride) {
747
        /*
748
         * Make the qsort stable by next comparing the perm order.
749
         * (Note that two perm entries will never be equal)
750
         */
751 1
        npy_intp aperm = ((const npy_stride_sort_item *)a)->perm,
752 1
                bperm = ((const npy_stride_sort_item *)b)->perm;
753 1
        return (aperm < bperm) ? -1 : 1;
754
    }
755 1
    if (astride > bstride) {
756
        return -1;
757
    }
758 1
    return 1;
759
}
760

761
/*NUMPY_API
762
 *
763
 * This function populates the first ndim elements
764
 * of strideperm with sorted descending by their absolute values.
765
 * For example, the stride array (4, -2, 12) becomes
766
 * [(2, 12), (0, 4), (1, -2)].
767
 */
768
NPY_NO_EXPORT void
769 1
PyArray_CreateSortedStridePerm(int ndim, npy_intp const *strides,
770
                        npy_stride_sort_item *out_strideperm)
771
{
772
    int i;
773

774
    /* Set up the strideperm values */
775 1
    for (i = 0; i < ndim; ++i) {
776 1
        out_strideperm[i].perm = i;
777 1
        out_strideperm[i].stride = strides[i];
778
    }
779

780
    /* Sort them */
781 1
    qsort(out_strideperm, ndim, sizeof(npy_stride_sort_item),
782
                                    &_npy_stride_sort_item_comparator);
783
}
784

785
static NPY_INLINE npy_intp
786
s_intp_abs(npy_intp x)
787
{
788 1
    return (x < 0) ? -x : x;
789
}
790

791

792

793
/*
794
 * Creates a sorted stride perm matching the KEEPORDER behavior
795
 * of the NpyIter object. Because this operates based on multiple
796
 * input strides, the 'stride' member of the npy_stride_sort_item
797
 * would be useless and we simply argsort a list of indices instead.
798
 *
799
 * The caller should have already validated that 'ndim' matches for
800
 * every array in the arrays list.
801
 */
802
NPY_NO_EXPORT void
803 1
PyArray_CreateMultiSortedStridePerm(int narrays, PyArrayObject **arrays,
804
                        int ndim, int *out_strideperm)
805
{
806
    int i0, i1, ipos, ax_j0, ax_j1, iarrays;
807

808
    /* Initialize the strideperm values to the identity. */
809 1
    for (i0 = 0; i0 < ndim; ++i0) {
810 1
        out_strideperm[i0] = i0;
811
    }
812

813
    /*
814
     * This is the same as the custom stable insertion sort in
815
     * the NpyIter object, but sorting in the reverse order as
816
     * in the iterator. The iterator sorts from smallest stride
817
     * to biggest stride (Fortran order), whereas here we sort
818
     * from biggest stride to smallest stride (C order).
819
     */
820 1
    for (i0 = 1; i0 < ndim; ++i0) {
821

822 1
        ipos = i0;
823 1
        ax_j0 = out_strideperm[i0];
824

825 1
        for (i1 = i0 - 1; i1 >= 0; --i1) {
826 1
            int ambig = 1, shouldswap = 0;
827

828 1
            ax_j1 = out_strideperm[i1];
829

830 1
            for (iarrays = 0; iarrays < narrays; ++iarrays) {
831 1
                if (PyArray_SHAPE(arrays[iarrays])[ax_j0] != 1 &&
832 1
                            PyArray_SHAPE(arrays[iarrays])[ax_j1] != 1) {
833 1
                    if (s_intp_abs(PyArray_STRIDES(arrays[iarrays])[ax_j0]) <=
834 1
                            s_intp_abs(PyArray_STRIDES(arrays[iarrays])[ax_j1])) {
835
                        /*
836
                         * Set swap even if it's not ambiguous already,
837
                         * because in the case of conflicts between
838
                         * different operands, C-order wins.
839
                         */
840
                        shouldswap = 0;
841
                    }
842
                    else {
843
                        /* Only set swap if it's still ambiguous */
844 1
                        if (ambig) {
845 1
                            shouldswap = 1;
846
                        }
847
                    }
848

849
                    /*
850
                     * A comparison has been done, so it's
851
                     * no longer ambiguous
852
                     */
853
                    ambig = 0;
854
                }
855
            }
856
            /*
857
             * If the comparison was unambiguous, either shift
858
             * 'ipos' to 'i1' or stop looking for an insertion point
859
             */
860 1
            if (!ambig) {
861 1
                if (shouldswap) {
862
                    ipos = i1;
863
                }
864
                else {
865
                    break;
866
                }
867
            }
868
        }
869

870
        /* Insert out_strideperm[i0] into the right place */
871 1
        if (ipos != i0) {
872 1
            for (i1 = i0; i1 > ipos; --i1) {
873 1
                out_strideperm[i1] = out_strideperm[i1-1];
874
            }
875 1
            out_strideperm[ipos] = ax_j0;
876
        }
877
    }
878
}
879

880
/*NUMPY_API
881
 * Ravel
882
 * Returns a contiguous array
883
 */
884
NPY_NO_EXPORT PyObject *
885 1
PyArray_Ravel(PyArrayObject *arr, NPY_ORDER order)
886
{
887 1
    PyArray_Dims newdim = {NULL,1};
888 1
    npy_intp val[1] = {-1};
889

890 1
    newdim.ptr = val;
891

892 1
    if (order == NPY_KEEPORDER) {
893
        /* This handles some corner cases, such as 0-d arrays as well */
894 1
        if (PyArray_IS_C_CONTIGUOUS(arr)) {
895
            order = NPY_CORDER;
896
        }
897 1
        else if (PyArray_IS_F_CONTIGUOUS(arr)) {
898
            order = NPY_FORTRANORDER;
899
        }
900
    }
901 1
    else if (order == NPY_ANYORDER) {
902 1
        order = PyArray_ISFORTRAN(arr) ? NPY_FORTRANORDER : NPY_CORDER;
903
    }
904

905 1
    if (order == NPY_CORDER && PyArray_IS_C_CONTIGUOUS(arr)) {
906 1
        return PyArray_Newshape(arr, &newdim, NPY_CORDER);
907
    }
908 1
    else if (order == NPY_FORTRANORDER && PyArray_IS_F_CONTIGUOUS(arr)) {
909 1
        return PyArray_Newshape(arr, &newdim, NPY_FORTRANORDER);
910
    }
911
    /* For KEEPORDER, check if we can make a flattened view */
912 1
    else if (order == NPY_KEEPORDER) {
913
        npy_stride_sort_item strideperm[NPY_MAXDIMS];
914
        npy_intp stride;
915 1
        int i, ndim = PyArray_NDIM(arr);
916

917 1
        PyArray_CreateSortedStridePerm(PyArray_NDIM(arr),
918 1
                                PyArray_STRIDES(arr), strideperm);
919

920
        /* The output array must be contiguous, so the first stride is fixed */
921 1
        stride = PyArray_ITEMSIZE(arr);
922

923 1
        for (i = ndim-1; i >= 0; --i) {
924 1
            if (PyArray_DIM(arr, strideperm[i].perm) == 1) {
925
                /* A size one dimension does not matter */
926 1
                continue;
927
            }
928 1
            if (strideperm[i].stride != stride) {
929
                break;
930
            }
931 1
            stride *= PyArray_DIM(arr, strideperm[i].perm);
932
        }
933

934
        /* If all the strides matched a contiguous layout, return a view */
935 1
        if (i < 0) {
936 1
            stride = PyArray_ITEMSIZE(arr);
937 1
            val[0] = PyArray_SIZE(arr);
938

939 1
            Py_INCREF(PyArray_DESCR(arr));
940 1
            return PyArray_NewFromDescrAndBase(
941 1
                    Py_TYPE(arr), PyArray_DESCR(arr),
942 1
                    1, val, &stride, PyArray_BYTES(arr),
943
                    PyArray_FLAGS(arr), (PyObject *)arr, (PyObject *)arr);
944
        }
945
    }
946

947 1
    return PyArray_Flatten(arr, order);
948
}
949

950
/*NUMPY_API
951
 * Flatten
952
 */
953
NPY_NO_EXPORT PyObject *
954 1
PyArray_Flatten(PyArrayObject *a, NPY_ORDER order)
955
{
956
    PyArrayObject *ret;
957
    npy_intp size;
958

959 1
    if (order == NPY_ANYORDER) {
960 1
        order = PyArray_ISFORTRAN(a) ? NPY_FORTRANORDER : NPY_CORDER;
961
    }
962

963 1
    size = PyArray_SIZE(a);
964 1
    Py_INCREF(PyArray_DESCR(a));
965 1
    ret = (PyArrayObject *)PyArray_NewFromDescr(Py_TYPE(a),
966
                               PyArray_DESCR(a),
967
                               1, &size,
968
                               NULL,
969
                               NULL,
970
                               0, (PyObject *)a);
971 1
    if (ret == NULL) {
972
        return NULL;
973
    }
974

975 1
    if (PyArray_CopyAsFlat(ret, a, order) < 0) {
976 0
        Py_DECREF(ret);
977
        return NULL;
978
    }
979
    return (PyObject *)ret;
980
}
981

982
/* See shape.h for parameters documentation */
983
NPY_NO_EXPORT PyObject *
984 1
build_shape_string(npy_intp n, npy_intp const *vals)
985
{
986
    npy_intp i;
987
    PyObject *ret, *tmp;
988

989
    /*
990
     * Negative dimension indicates "newaxis", which can
991
     * be discarded for printing if it's a leading dimension.
992
     * Find the first non-"newaxis" dimension.
993
     */
994 1
    i = 0;
995 1
    while (i < n && vals[i] < 0) {
996 0
        ++i;
997
    }
998

999 1
    if (i == n) {
1000 0
        return PyUnicode_FromFormat("()");
1001
    }
1002
    else {
1003 1
        ret = PyUnicode_FromFormat("(%" NPY_INTP_FMT, vals[i++]);
1004 1
        if (ret == NULL) {
1005
            return NULL;
1006
        }
1007
    }
1008

1009 1
    for (; i < n; ++i) {
1010 1
        if (vals[i] < 0) {
1011 0
            tmp = PyUnicode_FromString(",newaxis");
1012
        }
1013
        else {
1014 1
            tmp = PyUnicode_FromFormat(",%" NPY_INTP_FMT, vals[i]);
1015
        }
1016 1
        if (tmp == NULL) {
1017 0
            Py_DECREF(ret);
1018
            return NULL;
1019
        }
1020

1021 1
        PyUString_ConcatAndDel(&ret, tmp);
1022 1
        if (ret == NULL) {
1023
            return NULL;
1024
        }
1025
    }
1026

1027 1
    tmp = PyUnicode_FromFormat(")");
1028 1
    PyUString_ConcatAndDel(&ret, tmp);
1029 1
    return ret;
1030
}
1031

1032
/*NUMPY_API
1033
 *
1034
 * Removes the axes flagged as True from the array,
1035
 * modifying it in place. If an axis flagged for removal
1036
 * has a shape entry bigger than one, this effectively selects
1037
 * index zero for that axis.
1038
 *
1039
 * WARNING: If an axis flagged for removal has a shape equal to zero,
1040
 *          the array will point to invalid memory. The caller must
1041
 *          validate this!
1042
 *          If an axis flagged for removal has a shape larger than one,
1043
 *          the aligned flag (and in the future the contiguous flags),
1044
 *          may need explicit update.
1045
 *          (check also NPY_RELAXED_STRIDES_CHECKING)
1046
 *
1047
 * For example, this can be used to remove the reduction axes
1048
 * from a reduction result once its computation is complete.
1049
 */
1050
NPY_NO_EXPORT void
1051 1
PyArray_RemoveAxesInPlace(PyArrayObject *arr, const npy_bool *flags)
1052
{
1053 1
    PyArrayObject_fields *fa = (PyArrayObject_fields *)arr;
1054 1
    npy_intp *shape = fa->dimensions, *strides = fa->strides;
1055 1
    int idim, ndim = fa->nd, idim_out = 0;
1056

1057
    /* Compress the dimensions and strides */
1058 1
    for (idim = 0; idim < ndim; ++idim) {
1059 1
        if (!flags[idim]) {
1060 1
            shape[idim_out] = shape[idim];
1061 1
            strides[idim_out] = strides[idim];
1062 1
            ++idim_out;
1063
        }
1064
    }
1065

1066
    /* The final number of dimensions */
1067 1
    fa->nd = idim_out;
1068

1069
    /* May not be necessary for NPY_RELAXED_STRIDES_CHECKING (see comment) */
1070 1
    PyArray_UpdateFlags(arr, NPY_ARRAY_C_CONTIGUOUS | NPY_ARRAY_F_CONTIGUOUS);
1071
}

Read our documentation on viewing source code .

Loading