We are currently solving the following issue within our repository. Here is the issue text:
--- BEGIN ISSUE ---
matrix power is missing from cupy.linalg.
There is a corresponding `numpy.linalg.matrix_power` in NumPy. Is there a reason why this has not made it into CuPy?

If it helps, I am willing to take a shot at this one, though I might need a bit of hand-holding as I've not much experience with Cython code.

--- END ISSUE ---

Below are some code segments, each from a relevant file. One or more of these files may contain bugs.
--- BEGIN FILES ---
Path: `cupy/linalg/__init__.py`
Content:
```
1 # Functions from the following NumPy document
2 # https://docs.scipy.org/doc/numpy/reference/routines.linalg.html
3 
4 # "NOQA" to suppress flake8 warning
5 from cupy.linalg import decomposition  # NOQA
6 from cupy.linalg import eigenvalue  # NOQA
7 from cupy.linalg import einsum  # NOQA
8 from cupy.linalg import norms  # NOQA
9 from cupy.linalg.norms import det  # NOQA
10 from cupy.linalg.norms import matrix_rank  # NOQA
11 from cupy.linalg.norms import norm  # NOQA
12 from cupy.linalg.norms import slogdet  # NOQA
13 from cupy.linalg import product  # NOQA
14 from cupy.linalg import solve  # NOQA
15 
16 from cupy.linalg.decomposition import cholesky  # NOQA
17 from cupy.linalg.decomposition import qr  # NOQA
18 from cupy.linalg.decomposition import svd  # NOQA
19 
20 from cupy.linalg.eigenvalue import eigh  # NOQA
21 from cupy.linalg.eigenvalue import eigvalsh  # NOQA
22 
23 from cupy.linalg.solve import inv  # NOQA
24 from cupy.linalg.solve import pinv  # NOQA
25 from cupy.linalg.solve import solve  # NOQA
26 from cupy.linalg.solve import tensorinv  # NOQA
27 from cupy.linalg.solve import tensorsolve  # NOQA
28
```

Path: `cupy/linalg/product.py`
Content:
```
1 import collections
2 
3 import numpy
4 import six
5 
6 import cupy
7 from cupy import core
8 from cupy import internal
9 
10 
11 matmul = core.matmul
12 
13 
14 def dot(a, b, out=None):
15     """Returns a dot product of two arrays.
16 
17     For arrays with more than one axis, it computes the dot product along the
18     last axis of ``a`` and the second-to-last axis of ``b``. This is just a
19     matrix product if the both arrays are 2-D. For 1-D arrays, it uses their
20     unique axis as an axis to take dot product over.
21 
22     Args:
23         a (cupy.ndarray): The left argument.
24         b (cupy.ndarray): The right argument.
25         out (cupy.ndarray): Output array.
26 
27     Returns:
28         cupy.ndarray: The dot product of ``a`` and ``b``.
29 
30     .. seealso:: :func:`numpy.dot`
31 
32     """
33     # TODO(okuta): check type
34     return a.dot(b, out)
35 
36 
37 def vdot(a, b):
38     """Returns the dot product of two vectors.
39 
40     The input arrays are flattened into 1-D vectors and then it performs inner
41     product of these vectors.
42 
43     Args:
44         a (cupy.ndarray): The first argument.
45         b (cupy.ndarray): The second argument.
46 
47     Returns:
48         cupy.ndarray: Zero-dimensional array of the dot product result.
49 
50     .. seealso:: :func:`numpy.vdot`
51 
52     """
53     if a.size != b.size:
54         raise ValueError('Axis dimension mismatch')
55     if a.dtype.kind == 'c':
56         a = a.conj()
57 
58     return core.tensordot_core(a, b, None, 1, 1, a.size, ())
59 
60 
61 def inner(a, b):
62     """Returns the inner product of two arrays.
63 
64     It uses the last axis of each argument to take sum product.
65 
66     Args:
67         a (cupy.ndarray): The first argument.
68         b (cupy.ndarray): The second argument.
69 
70     Returns:
71         cupy.ndarray: The inner product of ``a`` and ``b``.
72 
73     .. seealso:: :func:`numpy.inner`
74 
75     """
76     a_ndim = a.ndim
77     b_ndim = b.ndim
78     if a_ndim == 0 or b_ndim == 0:
79         return cupy.multiply(a, b)
80 
81     a_axis = a_ndim - 1
82     b_axis = b_ndim - 1
83 
84     if a.shape[-1] != b.shape[-1]:
85         raise ValueError('Axis dimension mismatch')
86 
87     if a_axis:
88         a = cupy.rollaxis(a, a_axis, 0)
89     if b_axis:
90         b = cupy.rollaxis(b, b_axis, 0)
91 
92     ret_shape = a.shape[1:] + b.shape[1:]
93 
94     k = a.shape[0]
95     n = a.size // k
96     m = b.size // k
97 
98     return core.tensordot_core(a, b, None, n, m, k, ret_shape)
99 
100 
101 def outer(a, b, out=None):
102     """Returns the outer product of two vectors.
103 
104     The input arrays are flattened into 1-D vectors and then it performs outer
105     product of these vectors.
106 
107     Args:
108         a (cupy.ndarray): The first argument.
109         b (cupy.ndarray): The second argument.
110         out (cupy.ndarray): Output array.
111 
112     Returns:
113         cupy.ndarray: 2-D array of the outer product of ``a`` and ``b``.
114 
115     .. seealso:: :func:`numpy.outer`
116 
117     """
118     n = a.size
119     m = b.size
120     ret_shape = (n, m)
121 
122     if out is None:
123         return core.tensordot_core(a, b, None, n, m, 1, ret_shape)
124 
125     if out.size != n * m:
126         raise ValueError('Output array has an invalid size')
127     if out.flags.c_contiguous:
128         return core.tensordot_core(a, b, out, n, m, 1, ret_shape)
129     else:
130         out[:] = core.tensordot_core(a, b, None, n, m, 1, ret_shape)
131         return out
132 
133 
134 def tensordot(a, b, axes=2):
135     """Returns the tensor dot product of two arrays along specified axes.
136 
137     This is equivalent to compute dot product along the specified axes which
138     are treated as one axis by reshaping.
139 
140     Args:
141         a (cupy.ndarray): The first argument.
142         b (cupy.ndarray): The second argument.
143         axes:
144             - If it is an integer, then ``axes`` axes at the last of ``a`` and
145               the first of ``b`` are used.
146             - If it is a pair of sequences of integers, then these two
147               sequences specify the list of axes for ``a`` and ``b``. The
148               corresponding axes are paired for sum-product.
149 
150     Returns:
151         cupy.ndarray: The tensor dot product of ``a`` and ``b`` along the
152         axes specified by ``axes``.
153 
154     .. seealso:: :func:`numpy.tensordot`
155 
156     """
157     a_ndim = a.ndim
158     b_ndim = b.ndim
159     if a_ndim == 0 or b_ndim == 0:
160         if axes != 0 and axes != ((), ()):
161             raise ValueError('An input is zero-dim while axes has dimensions')
162         return cupy.multiply(a, b)
163 
164     if isinstance(axes, collections.Sequence):
165         if len(axes) != 2:
166             raise ValueError('Axes must consist of two arrays.')
167         a_axes, b_axes = axes
168         if numpy.isscalar(a_axes):
169             a_axes = a_axes,
170         if numpy.isscalar(b_axes):
171             b_axes = b_axes,
172     else:
173         a_axes = tuple(six.moves.range(a_ndim - axes, a_ndim))
174         b_axes = tuple(six.moves.range(axes))
175 
176     sum_ndim = len(a_axes)
177     if sum_ndim != len(b_axes):
178         raise ValueError('Axes length mismatch')
179 
180     for a_axis, b_axis in zip(a_axes, b_axes):
181         if a.shape[a_axis] != b.shape[b_axis]:
182             raise ValueError('Axis dimension mismatch')
183 
184     # Make the axes non-negative
185     a = _move_axes_to_head(a, [axis % a_ndim for axis in a_axes])
186     b = _move_axes_to_head(b, [axis % b_ndim for axis in b_axes])
187 
188     ret_shape = a.shape[sum_ndim:] + b.shape[sum_ndim:]
189 
190     k = internal.prod(a.shape[:sum_ndim])
191     n = a.size // k
192     m = b.size // k
193 
194     return core.tensordot_core(a, b, None, n, m, k, ret_shape)
195 
196 
197 # TODO(okuta): Implement matrix_power
198 
199 
200 def kron(a, b):
201     """Returns the kronecker product of two arrays.
202 
203     Args:
204         a (~cupy.ndarray): The first argument.
205         b (~cupy.ndarray): The second argument.
206 
207     Returns:
208         ~cupy.ndarray: Output array.
209 
210     .. seealso:: :func:`numpy.kron`
211 
212     """
213     a_ndim = a.ndim
214     b_ndim = b.ndim
215     if a_ndim == 0 or b_ndim == 0:
216         return cupy.multiply(a, b)
217 
218     ndim = b_ndim
219     a_shape = a.shape
220     b_shape = b.shape
221     if a_ndim != b_ndim:
222         if b_ndim > a_ndim:
223             a_shape = (1,) * (b_ndim - a_ndim) + a_shape
224         else:
225             b_shape = (1,) * (a_ndim - b_ndim) + b_shape
226             ndim = a_ndim
227 
228     axis = ndim - 1
229     out = core.tensordot_core(a, b, None, a.size, b.size, 1, a_shape + b_shape)
230     for _ in six.moves.range(ndim):
231         out = core.concatenate_method(out, axis=axis)
232 
233     return out
234 
235 
236 def _move_axes_to_head(a, axes):
237     # This function moves the axes of ``s`` to the head of the shape.
238     for idx, axis in enumerate(axes):
239         if idx != axis:
240             break
241     else:
242         return a
243 
244     return a.transpose(
245         axes + [i for i in six.moves.range(a.ndim) if i not in axes])
246
```
--- END FILES ---

Please first localize the bug based on the issue statement, and then generate a patch according to the `git diff` format fenced by three backticks.

Here is an example:

```diff
diff --git a/examples/server_async.py b/examples/server_async.py
--- a/examples/server_async.py
+++ b/examples/server_async.py
@@ -313,4 +313,4 @@
 
 
 if __name__ == "__main__":
-    asyncio.run(run_async_server("."), debug=True)
+    asyncio.run(run_async_server(), debug=True)
diff --git a/examples/server_sync.py b/examples/server_sync.py
--- a/examples/server_sync.py
+++ b/examples/server_sync.py
@@ -313,5 +313,5 @@
 
 
 if __name__ == "__main__":
-    server = run_sync_server(".")
+    server = run_sync_server()
     server.shutdown()

```

diff --git a/cupy/linalg/__init__.py b/cupy/linalg/__init__.py
--- a/cupy/linalg/__init__.py
+++ b/cupy/linalg/__init__.py
@@ -25,3 +25,5 @@
 from cupy.linalg.solve import solve  # NOQA
 from cupy.linalg.solve import tensorinv  # NOQA
 from cupy.linalg.solve import tensorsolve  # NOQA
+
+from cupy.linalg.product import matrix_power  # NOQA
diff --git a/cupy/linalg/product.py b/cupy/linalg/product.py
--- a/cupy/linalg/product.py
+++ b/cupy/linalg/product.py
@@ -7,6 +7,7 @@
 from cupy import core
 from cupy import internal
 
+from cupy.linalg.solve import inv
 
 matmul = core.matmul
 
@@ -194,7 +195,49 @@
     return core.tensordot_core(a, b, None, n, m, k, ret_shape)
 
 
-# TODO(okuta): Implement matrix_power
+def matrix_power(M, n):
+    """Raise a square matrix to the (integer) power `n`.
+
+    Args:
+        M (~cupy.ndarray): Matrix to raise by power n.
+        n (~int): Power to raise matrix to.
+
+    Returns:
+        ~cupy.ndarray: Output array.
+
+    .. note:: M must be of dtype `float32` or `float64`.
+
+    ..seealso:: :func:`numpy.linalg.matrix_power`
+    """
+    if M.ndim != 2 or M.shape[0] != M.shape[1]:
+        raise ValueError("input must be a square array")
+    if not isinstance(n, six.integer_types):
+        raise TypeError("exponent must be an integer")
+
+    if n == 0:
+        return cupy.identity(M.shape[0], dtype=M.dtype)
+    elif n < 0:
+        M = inv(M)
+        n *= -1
+
+    # short-cuts
+    if n <= 3:
+        if n == 1:
+            return M
+        elif n == 2:
+            return cupy.matmul(M, M)
+        else:
+            return cupy.matmul(cupy.matmul(M, M), M)
+
+    # binary decomposition to reduce the number of Matrix
+    # multiplications for n > 3.
+    result, Z = None, None
+    for b in cupy.binary_repr(n)[::-1]:
+        Z = M if Z is None else cupy.matmul(Z, Z)
+        if b == '1':
+            result = Z if result is None else cupy.matmul(result, Z)
+
+    return result
 
 
 def kron(a, b):

{"golden_diff": "diff --git a/cupy/linalg/__init__.py b/cupy/linalg/__init__.py\n--- a/cupy/linalg/__init__.py\n+++ b/cupy/linalg/__init__.py\n@@ -25,3 +25,5 @@\n from cupy.linalg.solve import solve  # NOQA\n from cupy.linalg.solve import tensorinv  # NOQA\n from cupy.linalg.solve import tensorsolve  # NOQA\n+\n+from cupy.linalg.product import matrix_power  # NOQA\ndiff --git a/cupy/linalg/product.py b/cupy/linalg/product.py\n--- a/cupy/linalg/product.py\n+++ b/cupy/linalg/product.py\n@@ -7,6 +7,7 @@\n from cupy import core\n from cupy import internal\n \n+from cupy.linalg.solve import inv\n \n matmul = core.matmul\n \n@@ -194,7 +195,49 @@\n     return core.tensordot_core(a, b, None, n, m, k, ret_shape)\n \n \n-# TODO(okuta): Implement matrix_power\n+def matrix_power(M, n):\n+    \"\"\"Raise a square matrix to the (integer) power `n`.\n+\n+    Args:\n+        M (~cupy.ndarray): Matrix to raise by power n.\n+        n (~int): Power to raise matrix to.\n+\n+    Returns:\n+        ~cupy.ndarray: Output array.\n+\n+    .. note:: M must be of dtype `float32` or `float64`.\n+\n+    ..seealso:: :func:`numpy.linalg.matrix_power`\n+    \"\"\"\n+    if M.ndim != 2 or M.shape[0] != M.shape[1]:\n+        raise ValueError(\"input must be a square array\")\n+    if not isinstance(n, six.integer_types):\n+        raise TypeError(\"exponent must be an integer\")\n+\n+    if n == 0:\n+        return cupy.identity(M.shape[0], dtype=M.dtype)\n+    elif n < 0:\n+        M = inv(M)\n+        n *= -1\n+\n+    # short-cuts\n+    if n <= 3:\n+        if n == 1:\n+            return M\n+        elif n == 2:\n+            return cupy.matmul(M, M)\n+        else:\n+            return cupy.matmul(cupy.matmul(M, M), M)\n+\n+    # binary decomposition to reduce the number of Matrix\n+    # multiplications for n > 3.\n+    result, Z = None, None\n+    for b in cupy.binary_repr(n)[::-1]:\n+        Z = M if Z is None else cupy.matmul(Z, Z)\n+        if b == '1':\n+            result = Z if result is None else cupy.matmul(result, Z)\n+\n+    return result\n \n \n def kron(a, b):\n", "issue": "matrix power is missing from cupy.linalg.\nThere is a corresponding `numpy.linalg.matrix_power` in NumPy. Is there a reason why this has not made it into CuPy?\r\n\r\nIf it helps, I am willing to take a shot at this one, though I might need a bit of hand-holding as I've not much experience with Cython code.\n", "before_files": [{"content": "# Functions from the following NumPy document\n# https://docs.scipy.org/doc/numpy/reference/routines.linalg.html\n\n# \"NOQA\" to suppress flake8 warning\nfrom cupy.linalg import decomposition  # NOQA\nfrom cupy.linalg import eigenvalue  # NOQA\nfrom cupy.linalg import einsum  # NOQA\nfrom cupy.linalg import norms  # NOQA\nfrom cupy.linalg.norms import det  # NOQA\nfrom cupy.linalg.norms import matrix_rank  # NOQA\nfrom cupy.linalg.norms import norm  # NOQA\nfrom cupy.linalg.norms import slogdet  # NOQA\nfrom cupy.linalg import product  # NOQA\nfrom cupy.linalg import solve  # NOQA\n\nfrom cupy.linalg.decomposition import cholesky  # NOQA\nfrom cupy.linalg.decomposition import qr  # NOQA\nfrom cupy.linalg.decomposition import svd  # NOQA\n\nfrom cupy.linalg.eigenvalue import eigh  # NOQA\nfrom cupy.linalg.eigenvalue import eigvalsh  # NOQA\n\nfrom cupy.linalg.solve import inv  # NOQA\nfrom cupy.linalg.solve import pinv  # NOQA\nfrom cupy.linalg.solve import solve  # NOQA\nfrom cupy.linalg.solve import tensorinv  # NOQA\nfrom cupy.linalg.solve import tensorsolve  # NOQA\n", "path": "cupy/linalg/__init__.py"}, {"content": "import collections\n\nimport numpy\nimport six\n\nimport cupy\nfrom cupy import core\nfrom cupy import internal\n\n\nmatmul = core.matmul\n\n\ndef dot(a, b, out=None):\n    \"\"\"Returns a dot product of two arrays.\n\n    For arrays with more than one axis, it computes the dot product along the\n    last axis of ``a`` and the second-to-last axis of ``b``. This is just a\n    matrix product if the both arrays are 2-D. For 1-D arrays, it uses their\n    unique axis as an axis to take dot product over.\n\n    Args:\n        a (cupy.ndarray): The left argument.\n        b (cupy.ndarray): The right argument.\n        out (cupy.ndarray): Output array.\n\n    Returns:\n        cupy.ndarray: The dot product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.dot`\n\n    \"\"\"\n    # TODO(okuta): check type\n    return a.dot(b, out)\n\n\ndef vdot(a, b):\n    \"\"\"Returns the dot product of two vectors.\n\n    The input arrays are flattened into 1-D vectors and then it performs inner\n    product of these vectors.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n\n    Returns:\n        cupy.ndarray: Zero-dimensional array of the dot product result.\n\n    .. seealso:: :func:`numpy.vdot`\n\n    \"\"\"\n    if a.size != b.size:\n        raise ValueError('Axis dimension mismatch')\n    if a.dtype.kind == 'c':\n        a = a.conj()\n\n    return core.tensordot_core(a, b, None, 1, 1, a.size, ())\n\n\ndef inner(a, b):\n    \"\"\"Returns the inner product of two arrays.\n\n    It uses the last axis of each argument to take sum product.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n\n    Returns:\n        cupy.ndarray: The inner product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.inner`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        return cupy.multiply(a, b)\n\n    a_axis = a_ndim - 1\n    b_axis = b_ndim - 1\n\n    if a.shape[-1] != b.shape[-1]:\n        raise ValueError('Axis dimension mismatch')\n\n    if a_axis:\n        a = cupy.rollaxis(a, a_axis, 0)\n    if b_axis:\n        b = cupy.rollaxis(b, b_axis, 0)\n\n    ret_shape = a.shape[1:] + b.shape[1:]\n\n    k = a.shape[0]\n    n = a.size // k\n    m = b.size // k\n\n    return core.tensordot_core(a, b, None, n, m, k, ret_shape)\n\n\ndef outer(a, b, out=None):\n    \"\"\"Returns the outer product of two vectors.\n\n    The input arrays are flattened into 1-D vectors and then it performs outer\n    product of these vectors.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n        out (cupy.ndarray): Output array.\n\n    Returns:\n        cupy.ndarray: 2-D array of the outer product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.outer`\n\n    \"\"\"\n    n = a.size\n    m = b.size\n    ret_shape = (n, m)\n\n    if out is None:\n        return core.tensordot_core(a, b, None, n, m, 1, ret_shape)\n\n    if out.size != n * m:\n        raise ValueError('Output array has an invalid size')\n    if out.flags.c_contiguous:\n        return core.tensordot_core(a, b, out, n, m, 1, ret_shape)\n    else:\n        out[:] = core.tensordot_core(a, b, None, n, m, 1, ret_shape)\n        return out\n\n\ndef tensordot(a, b, axes=2):\n    \"\"\"Returns the tensor dot product of two arrays along specified axes.\n\n    This is equivalent to compute dot product along the specified axes which\n    are treated as one axis by reshaping.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n        axes:\n            - If it is an integer, then ``axes`` axes at the last of ``a`` and\n              the first of ``b`` are used.\n            - If it is a pair of sequences of integers, then these two\n              sequences specify the list of axes for ``a`` and ``b``. The\n              corresponding axes are paired for sum-product.\n\n    Returns:\n        cupy.ndarray: The tensor dot product of ``a`` and ``b`` along the\n        axes specified by ``axes``.\n\n    .. seealso:: :func:`numpy.tensordot`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        if axes != 0 and axes != ((), ()):\n            raise ValueError('An input is zero-dim while axes has dimensions')\n        return cupy.multiply(a, b)\n\n    if isinstance(axes, collections.Sequence):\n        if len(axes) != 2:\n            raise ValueError('Axes must consist of two arrays.')\n        a_axes, b_axes = axes\n        if numpy.isscalar(a_axes):\n            a_axes = a_axes,\n        if numpy.isscalar(b_axes):\n            b_axes = b_axes,\n    else:\n        a_axes = tuple(six.moves.range(a_ndim - axes, a_ndim))\n        b_axes = tuple(six.moves.range(axes))\n\n    sum_ndim = len(a_axes)\n    if sum_ndim != len(b_axes):\n        raise ValueError('Axes length mismatch')\n\n    for a_axis, b_axis in zip(a_axes, b_axes):\n        if a.shape[a_axis] != b.shape[b_axis]:\n            raise ValueError('Axis dimension mismatch')\n\n    # Make the axes non-negative\n    a = _move_axes_to_head(a, [axis % a_ndim for axis in a_axes])\n    b = _move_axes_to_head(b, [axis % b_ndim for axis in b_axes])\n\n    ret_shape = a.shape[sum_ndim:] + b.shape[sum_ndim:]\n\n    k = internal.prod(a.shape[:sum_ndim])\n    n = a.size // k\n    m = b.size // k\n\n    return core.tensordot_core(a, b, None, n, m, k, ret_shape)\n\n\n# TODO(okuta): Implement matrix_power\n\n\ndef kron(a, b):\n    \"\"\"Returns the kronecker product of two arrays.\n\n    Args:\n        a (~cupy.ndarray): The first argument.\n        b (~cupy.ndarray): The second argument.\n\n    Returns:\n        ~cupy.ndarray: Output array.\n\n    .. seealso:: :func:`numpy.kron`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        return cupy.multiply(a, b)\n\n    ndim = b_ndim\n    a_shape = a.shape\n    b_shape = b.shape\n    if a_ndim != b_ndim:\n        if b_ndim > a_ndim:\n            a_shape = (1,) * (b_ndim - a_ndim) + a_shape\n        else:\n            b_shape = (1,) * (a_ndim - b_ndim) + b_shape\n            ndim = a_ndim\n\n    axis = ndim - 1\n    out = core.tensordot_core(a, b, None, a.size, b.size, 1, a_shape + b_shape)\n    for _ in six.moves.range(ndim):\n        out = core.concatenate_method(out, axis=axis)\n\n    return out\n\n\ndef _move_axes_to_head(a, axes):\n    # This function moves the axes of ``s`` to the head of the shape.\n    for idx, axis in enumerate(axes):\n        if idx != axis:\n            break\n    else:\n        return a\n\n    return a.transpose(\n        axes + [i for i in six.moves.range(a.ndim) if i not in axes])\n", "path": "cupy/linalg/product.py"}], "after_files": [{"content": "# Functions from the following NumPy document\n# https://docs.scipy.org/doc/numpy/reference/routines.linalg.html\n\n# \"NOQA\" to suppress flake8 warning\nfrom cupy.linalg import decomposition  # NOQA\nfrom cupy.linalg import eigenvalue  # NOQA\nfrom cupy.linalg import einsum  # NOQA\nfrom cupy.linalg import norms  # NOQA\nfrom cupy.linalg.norms import det  # NOQA\nfrom cupy.linalg.norms import matrix_rank  # NOQA\nfrom cupy.linalg.norms import norm  # NOQA\nfrom cupy.linalg.norms import slogdet  # NOQA\nfrom cupy.linalg import product  # NOQA\nfrom cupy.linalg import solve  # NOQA\n\nfrom cupy.linalg.decomposition import cholesky  # NOQA\nfrom cupy.linalg.decomposition import qr  # NOQA\nfrom cupy.linalg.decomposition import svd  # NOQA\n\nfrom cupy.linalg.eigenvalue import eigh  # NOQA\nfrom cupy.linalg.eigenvalue import eigvalsh  # NOQA\n\nfrom cupy.linalg.solve import inv  # NOQA\nfrom cupy.linalg.solve import pinv  # NOQA\nfrom cupy.linalg.solve import solve  # NOQA\nfrom cupy.linalg.solve import tensorinv  # NOQA\nfrom cupy.linalg.solve import tensorsolve  # NOQA\n\nfrom cupy.linalg.product import matrix_power  # NOQA\n", "path": "cupy/linalg/__init__.py"}, {"content": "import collections\n\nimport numpy\nimport six\n\nimport cupy\nfrom cupy import core\nfrom cupy import internal\n\nfrom cupy.linalg.solve import inv\n\nmatmul = core.matmul\n\n\ndef dot(a, b, out=None):\n    \"\"\"Returns a dot product of two arrays.\n\n    For arrays with more than one axis, it computes the dot product along the\n    last axis of ``a`` and the second-to-last axis of ``b``. This is just a\n    matrix product if the both arrays are 2-D. For 1-D arrays, it uses their\n    unique axis as an axis to take dot product over.\n\n    Args:\n        a (cupy.ndarray): The left argument.\n        b (cupy.ndarray): The right argument.\n        out (cupy.ndarray): Output array.\n\n    Returns:\n        cupy.ndarray: The dot product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.dot`\n\n    \"\"\"\n    # TODO(okuta): check type\n    return a.dot(b, out)\n\n\ndef vdot(a, b):\n    \"\"\"Returns the dot product of two vectors.\n\n    The input arrays are flattened into 1-D vectors and then it performs inner\n    product of these vectors.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n\n    Returns:\n        cupy.ndarray: Zero-dimensional array of the dot product result.\n\n    .. seealso:: :func:`numpy.vdot`\n\n    \"\"\"\n    if a.size != b.size:\n        raise ValueError('Axis dimension mismatch')\n    if a.dtype.kind == 'c':\n        a = a.conj()\n\n    return core.tensordot_core(a, b, None, 1, 1, a.size, ())\n\n\ndef inner(a, b):\n    \"\"\"Returns the inner product of two arrays.\n\n    It uses the last axis of each argument to take sum product.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n\n    Returns:\n        cupy.ndarray: The inner product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.inner`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        return cupy.multiply(a, b)\n\n    a_axis = a_ndim - 1\n    b_axis = b_ndim - 1\n\n    if a.shape[-1] != b.shape[-1]:\n        raise ValueError('Axis dimension mismatch')\n\n    if a_axis:\n        a = cupy.rollaxis(a, a_axis, 0)\n    if b_axis:\n        b = cupy.rollaxis(b, b_axis, 0)\n\n    ret_shape = a.shape[1:] + b.shape[1:]\n\n    k = a.shape[0]\n    n = a.size // k\n    m = b.size // k\n\n    return core.tensordot_core(a, b, None, n, m, k, ret_shape)\n\n\ndef outer(a, b, out=None):\n    \"\"\"Returns the outer product of two vectors.\n\n    The input arrays are flattened into 1-D vectors and then it performs outer\n    product of these vectors.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n        out (cupy.ndarray): Output array.\n\n    Returns:\n        cupy.ndarray: 2-D array of the outer product of ``a`` and ``b``.\n\n    .. seealso:: :func:`numpy.outer`\n\n    \"\"\"\n    n = a.size\n    m = b.size\n    ret_shape = (n, m)\n\n    if out is None:\n        return core.tensordot_core(a, b, None, n, m, 1, ret_shape)\n\n    if out.size != n * m:\n        raise ValueError('Output array has an invalid size')\n    if out.flags.c_contiguous:\n        return core.tensordot_core(a, b, out, n, m, 1, ret_shape)\n    else:\n        out[:] = core.tensordot_core(a, b, None, n, m, 1, ret_shape)\n        return out\n\n\ndef tensordot(a, b, axes=2):\n    \"\"\"Returns the tensor dot product of two arrays along specified axes.\n\n    This is equivalent to compute dot product along the specified axes which\n    are treated as one axis by reshaping.\n\n    Args:\n        a (cupy.ndarray): The first argument.\n        b (cupy.ndarray): The second argument.\n        axes:\n            - If it is an integer, then ``axes`` axes at the last of ``a`` and\n              the first of ``b`` are used.\n            - If it is a pair of sequences of integers, then these two\n              sequences specify the list of axes for ``a`` and ``b``. The\n              corresponding axes are paired for sum-product.\n\n    Returns:\n        cupy.ndarray: The tensor dot product of ``a`` and ``b`` along the\n        axes specified by ``axes``.\n\n    .. seealso:: :func:`numpy.tensordot`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        if axes != 0 and axes != ((), ()):\n            raise ValueError('An input is zero-dim while axes has dimensions')\n        return cupy.multiply(a, b)\n\n    if isinstance(axes, collections.Sequence):\n        if len(axes) != 2:\n            raise ValueError('Axes must consist of two arrays.')\n        a_axes, b_axes = axes\n        if numpy.isscalar(a_axes):\n            a_axes = a_axes,\n        if numpy.isscalar(b_axes):\n            b_axes = b_axes,\n    else:\n        a_axes = tuple(six.moves.range(a_ndim - axes, a_ndim))\n        b_axes = tuple(six.moves.range(axes))\n\n    sum_ndim = len(a_axes)\n    if sum_ndim != len(b_axes):\n        raise ValueError('Axes length mismatch')\n\n    for a_axis, b_axis in zip(a_axes, b_axes):\n        if a.shape[a_axis] != b.shape[b_axis]:\n            raise ValueError('Axis dimension mismatch')\n\n    # Make the axes non-negative\n    a = _move_axes_to_head(a, [axis % a_ndim for axis in a_axes])\n    b = _move_axes_to_head(b, [axis % b_ndim for axis in b_axes])\n\n    ret_shape = a.shape[sum_ndim:] + b.shape[sum_ndim:]\n\n    k = internal.prod(a.shape[:sum_ndim])\n    n = a.size // k\n    m = b.size // k\n\n    return core.tensordot_core(a, b, None, n, m, k, ret_shape)\n\n\ndef matrix_power(M, n):\n    \"\"\"Raise a square matrix to the (integer) power `n`.\n\n    Args:\n        M (~cupy.ndarray): Matrix to raise by power n.\n        n (~int): Power to raise matrix to.\n\n    Returns:\n        ~cupy.ndarray: Output array.\n\n    .. note:: M must be of dtype `float32` or `float64`.\n\n    ..seealso:: :func:`numpy.linalg.matrix_power`\n    \"\"\"\n    if M.ndim != 2 or M.shape[0] != M.shape[1]:\n        raise ValueError(\"input must be a square array\")\n    if not isinstance(n, six.integer_types):\n        raise TypeError(\"exponent must be an integer\")\n\n    if n == 0:\n        return cupy.identity(M.shape[0], dtype=M.dtype)\n    elif n < 0:\n        M = inv(M)\n        n *= -1\n\n    # short-cuts\n    if n <= 3:\n        if n == 1:\n            return M\n        elif n == 2:\n            return cupy.matmul(M, M)\n        else:\n            return cupy.matmul(cupy.matmul(M, M), M)\n\n    # binary decomposition to reduce the number of Matrix\n    # multiplications for n > 3.\n    result, Z = None, None\n    for b in cupy.binary_repr(n)[::-1]:\n        Z = M if Z is None else cupy.matmul(Z, Z)\n        if b == '1':\n            result = Z if result is None else cupy.matmul(result, Z)\n\n    return result\n\n\ndef kron(a, b):\n    \"\"\"Returns the kronecker product of two arrays.\n\n    Args:\n        a (~cupy.ndarray): The first argument.\n        b (~cupy.ndarray): The second argument.\n\n    Returns:\n        ~cupy.ndarray: Output array.\n\n    .. seealso:: :func:`numpy.kron`\n\n    \"\"\"\n    a_ndim = a.ndim\n    b_ndim = b.ndim\n    if a_ndim == 0 or b_ndim == 0:\n        return cupy.multiply(a, b)\n\n    ndim = b_ndim\n    a_shape = a.shape\n    b_shape = b.shape\n    if a_ndim != b_ndim:\n        if b_ndim > a_ndim:\n            a_shape = (1,) * (b_ndim - a_ndim) + a_shape\n        else:\n            b_shape = (1,) * (a_ndim - b_ndim) + b_shape\n            ndim = a_ndim\n\n    axis = ndim - 1\n    out = core.tensordot_core(a, b, None, a.size, b.size, 1, a_shape + b_shape)\n    for _ in six.moves.range(ndim):\n        out = core.concatenate_method(out, axis=axis)\n\n    return out\n\n\ndef _move_axes_to_head(a, axes):\n    # This function moves the axes of ``s`` to the head of the shape.\n    for idx, axis in enumerate(axes):\n        if idx != axis:\n            break\n    else:\n        return a\n\n    return a.transpose(\n        axes + [i for i in six.moves.range(a.ndim) if i not in axes])\n", "path": "cupy/linalg/product.py"}]}