# Source code for sympy.matrices.expressions.matexpr

from __future__ import print_function, division

from functools import wraps

from sympy.core import S, Symbol, Tuple, Integer, Basic, Expr, Eq
from sympy.core.decorators import call_highest_priority
from sympy.core.compatibility import range
from sympy.core.sympify import SympifyError, sympify
from sympy.functions.special.tensor_functions import KroneckerDelta
from sympy.matrices import ShapeError
from sympy.simplify import simplify
from sympy.utilities.misc import filldedent

def _sympifyit(arg, retval=None):
# This version of _sympifyit sympifies MutableMatrix objects
def deco(func):
@wraps(func)
def __sympifyit_wrapper(a, b):
try:
b = sympify(b, strict=True)
return func(a, b)
except SympifyError:
return retval

return __sympifyit_wrapper

return deco

[docs]class MatrixExpr(Basic):
""" Superclass for Matrix Expressions

MatrixExprs represent abstract matrices, linear transformations represented
within a particular basis.

Examples
========

>>> from sympy import MatrixSymbol
>>> A = MatrixSymbol('A', 3, 3)
>>> y = MatrixSymbol('y', 3, 1)
>>> x = (A.T*A).I * A * y

========
MatrixSymbol
MatMul
Transpose
Inverse
"""

# Should not be considered iterable by the
# sympy.core.compatibility.iterable function. Subclass that actually are
# iterable (i.e., explicit matrices) should set this to True.
_iterable = False

_op_priority = 11.0

is_Matrix = True
is_MatrixExpr = True
is_Identity = None
is_Inverse = False
is_Transpose = False
is_ZeroMatrix = False
is_MatMul = False

is_commutative = False

def __new__(cls, *args, **kwargs):
args = map(sympify, args)
return Basic.__new__(cls, *args, **kwargs)

# The following is adapted from the core Expr object
def __neg__(self):
return MatMul(S.NegativeOne, self).doit()

def __abs__(self):
raise NotImplementedError

@_sympifyit('other', NotImplemented)

@_sympifyit('other', NotImplemented)

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rsub__')
def __sub__(self, other):

@_sympifyit('other', NotImplemented)
@call_highest_priority('__sub__')
def __rsub__(self, other):

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rmul__')
def __mul__(self, other):
return MatMul(self, other).doit()

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rmul__')
def __matmul__(self, other):
return MatMul(self, other).doit()

@_sympifyit('other', NotImplemented)
@call_highest_priority('__mul__')
def __rmul__(self, other):
return MatMul(other, self).doit()

@_sympifyit('other', NotImplemented)
@call_highest_priority('__mul__')
def __rmatmul__(self, other):
return MatMul(other, self).doit()

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rpow__')
def __pow__(self, other):
if not self.is_square:
raise ShapeError("Power of non-square matrix %s" % self)
elif self.is_Identity:
return self
elif other is S.NegativeOne:
return Inverse(self)
elif other is S.Zero:
return Identity(self.rows)
elif other is S.One:
return self
return MatPow(self, other)

@_sympifyit('other', NotImplemented)
@call_highest_priority('__pow__')
def __rpow__(self, other):
raise NotImplementedError("Matrix Power not defined")

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rdiv__')
def __div__(self, other):
return self * other**S.NegativeOne

@_sympifyit('other', NotImplemented)
@call_highest_priority('__div__')
def __rdiv__(self, other):
raise NotImplementedError()
#return MatMul(other, Pow(self, S.NegativeOne))

__truediv__ = __div__
__rtruediv__ = __rdiv__

@property
def rows(self):
return self.shape[0]

@property
def cols(self):
return self.shape[1]

@property
def is_square(self):
return self.rows == self.cols

def _eval_conjugate(self):
from sympy.matrices.expressions.transpose import Transpose

def as_real_imag(self):
from sympy import I
real = (S(1)/2) * (self + self._eval_conjugate())
im = (self - self._eval_conjugate())/(2*I)
return (real, im)

def _eval_inverse(self):
from sympy.matrices.expressions.inverse import Inverse
return Inverse(self)

def _eval_transpose(self):
return Transpose(self)

def _eval_power(self, exp):
return MatPow(self, exp)

def _eval_simplify(self, **kwargs):
if self.is_Atom:
return self
else:
return self.__class__(*[simplify(x, **kwargs) for x in self.args])

def _entry(self, i, j):
raise NotImplementedError(
"Indexing not implemented for %s" % self.__class__.__name__)

[docs]    def as_coeff_Mul(self, rational=False):
"""Efficiently extract the coefficient of a product. """
return S.One, self

def conjugate(self):
return conjugate(self)

def transpose(self):
from sympy.matrices.expressions.transpose import transpose
return transpose(self)

T = property(transpose, None, None, 'Matrix transposition.')

def inverse(self):
return self._eval_inverse()

@property
def I(self):
return self.inverse()

def valid_index(self, i, j):
def is_valid(idx):
return isinstance(idx, (int, Integer, Symbol, Expr))
return (is_valid(i) and is_valid(j) and
(self.rows is None or
(0 <= i) != False and (i < self.rows) != False) and
(0 <= j) != False and (j < self.cols) != False)

def __getitem__(self, key):
if not isinstance(key, tuple) and isinstance(key, slice):
from sympy.matrices.expressions.slice import MatrixSlice
return MatrixSlice(self, key, (0, None, 1))
if isinstance(key, tuple) and len(key) == 2:
i, j = key
if isinstance(i, slice) or isinstance(j, slice):
from sympy.matrices.expressions.slice import MatrixSlice
return MatrixSlice(self, i, j)
i, j = sympify(i), sympify(j)
if self.valid_index(i, j) != False:
return self._entry(i, j)
else:
raise IndexError("Invalid indices (%s, %s)" % (i, j))
elif isinstance(key, (int, Integer)):
# row-wise decomposition of matrix
rows, cols = self.shape
# allow single indexing if number of columns is known
if not isinstance(cols, Integer):
raise IndexError(filldedent('''
Single indexing is only supported when the number
of columns is known.'''))
key = sympify(key)
i = key // cols
j = key % cols
if self.valid_index(i, j) != False:
return self._entry(i, j)
else:
raise IndexError("Invalid index %s" % key)
elif isinstance(key, (Symbol, Expr)):
raise IndexError(filldedent('''
Only integers may be used when addressing the matrix
with a single index.'''))
raise IndexError("Invalid index, wanted %s[i,j]" % self)

[docs]    def as_explicit(self):
"""
Returns a dense Matrix with elements represented explicitly

Returns an object of type ImmutableDenseMatrix.

Examples
========

>>> from sympy import Identity
>>> I = Identity(3)
>>> I
I
>>> I.as_explicit()
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 1]])

========
as_mutable: returns mutable Matrix type

"""
from sympy.matrices.immutable import ImmutableDenseMatrix
return ImmutableDenseMatrix([[    self[i, j]
for j in range(self.cols)]
for i in range(self.rows)])

[docs]    def as_mutable(self):
"""
Returns a dense, mutable matrix with elements represented explicitly

Examples
========

>>> from sympy import Identity
>>> I = Identity(3)
>>> I
I
>>> I.shape
(3, 3)
>>> I.as_mutable()
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 1]])

========
as_explicit: returns ImmutableDenseMatrix
"""
return self.as_explicit().as_mutable()

def __array__(self):
from numpy import empty
a = empty(self.shape, dtype=object)
for i in range(self.rows):
for j in range(self.cols):
a[i, j] = self[i, j]
return a

[docs]    def equals(self, other):
"""
Test elementwise equality between matrices, potentially of different
types

>>> from sympy import Identity, eye
>>> Identity(3).equals(eye(3))
True
"""
return self.as_explicit().equals(other)

def canonicalize(self):
return self

def as_coeff_mmul(self):
return 1, MatMul(self)

class MatrixElement(Expr):
parent = property(lambda self: self.args[0])
i = property(lambda self: self.args[1])
j = property(lambda self: self.args[2])
_diff_wrt = True
is_symbol = True
is_commutative = True

def __new__(cls, name, n, m):
n, m = map(sympify, (n, m))
from sympy import MatrixBase
if isinstance(name, (MatrixBase,)):
if n.is_Integer and m.is_Integer:
return name[n, m]
name = sympify(name)
obj = Expr.__new__(cls, name, n, m)
return obj

def doit(self, **kwargs):
deep = kwargs.get('deep', True)
if deep:
args = [arg.doit(**kwargs) for arg in self.args]
else:
args = self.args
return args[0][args[1], args[2]]

def _eval_derivative(self, v):
if not isinstance(v, MatrixElement):
from sympy import MatrixBase
if isinstance(self.parent, MatrixBase):
return self.parent.diff(v)[self.i, self.j]
return S.Zero

if self.args[0] != v.args[0]:
return S.Zero

return KroneckerDelta(self.args[1], v.args[1])*KroneckerDelta(self.args[2], v.args[2])

[docs]class MatrixSymbol(MatrixExpr):
"""Symbolic representation of a Matrix object

Creates a SymPy Symbol to represent a Matrix. This matrix has a shape and
can be included in Matrix Expressions

>>> from sympy import MatrixSymbol, Identity
>>> A = MatrixSymbol('A', 3, 4) # A 3 by 4 Matrix
>>> B = MatrixSymbol('B', 4, 3) # A 4 by 3 Matrix
>>> A.shape
(3, 4)
>>> 2*A*B + Identity(3)
I + 2*A*B
"""
is_commutative = False

def __new__(cls, name, n, m):
n, m = sympify(n), sympify(m)
obj = Basic.__new__(cls, name, n, m)
return obj

def _hashable_content(self):
return(self.name, self.shape)

@property
def shape(self):
return self.args[1:3]

@property
def name(self):
return self.args[0]

def _eval_subs(self, old, new):
# only do substitutions in shape
shape = Tuple(*self.shape)._subs(old, new)
return MatrixSymbol(self.name, *shape)

def __call__(self, *args):
raise TypeError( "%s object is not callable" % self.__class__ )

def _entry(self, i, j):
return MatrixElement(self, i, j)

@property
def free_symbols(self):
return set((self,))

def doit(self, **hints):
if hints.get('deep', True):
return type(self)(self.name, self.args[1].doit(**hints),
self.args[2].doit(**hints))
else:
return self

def _eval_simplify(self, **kwargs):
return self

[docs]class Identity(MatrixExpr):
"""The Matrix Identity I - multiplicative identity

>>> from sympy.matrices import Identity, MatrixSymbol
>>> A = MatrixSymbol('A', 3, 5)
>>> I = Identity(3)
>>> I*A
A
"""

is_Identity = True

def __new__(cls, n):
return super(Identity, cls).__new__(cls, sympify(n))

@property
def rows(self):
return self.args[0]

@property
def cols(self):
return self.args[0]

@property
def shape(self):
return (self.args[0], self.args[0])

def _eval_transpose(self):
return self

def _eval_trace(self):
return self.rows

def _eval_inverse(self):
return self

def conjugate(self):
return self

def _entry(self, i, j):
eq = Eq(i, j)
if eq is S.true:
return S.One
elif eq is S.false:
return S.Zero
return KroneckerDelta(i, j)

def _eval_determinant(self):
return S.One

[docs]class ZeroMatrix(MatrixExpr):
"""The Matrix Zero 0 - additive identity

>>> from sympy import MatrixSymbol, ZeroMatrix
>>> A = MatrixSymbol('A', 3, 5)
>>> Z = ZeroMatrix(3, 5)
>>> A+Z
A
>>> Z*A.T
0
"""
is_ZeroMatrix = True

def __new__(cls, m, n):
return super(ZeroMatrix, cls).__new__(cls, m, n)

@property
def shape(self):
return (self.args[0], self.args[1])

@_sympifyit('other', NotImplemented)
@call_highest_priority('__rpow__')
def __pow__(self, other):
if other != 1 and not self.is_square:
raise ShapeError("Power of non-square matrix %s" % self)
if other == 0:
return Identity(self.rows)
if other < 1:
raise ValueError("Matrix det == 0; not invertible.")
return self

def _eval_transpose(self):
return ZeroMatrix(self.cols, self.rows)

def _eval_trace(self):
return S.Zero

def _eval_determinant(self):
return S.Zero

def conjugate(self):
return self

def _entry(self, i, j):
return S.Zero

def __nonzero__(self):
return False

__bool__ = __nonzero__

def matrix_symbols(expr):
return [sym for sym in expr.free_symbols if sym.is_Matrix]

from .matmul import MatMul