123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868 |
- from typing import Tuple as tTuple
- from functools import wraps
- from sympy.core import S, Integer, Basic, Mul, Add
- from sympy.core.assumptions import check_assumptions
- from sympy.core.decorators import call_highest_priority
- from sympy.core.expr import Expr, ExprBuilder
- from sympy.core.logic import FuzzyBool
- from sympy.core.symbol import Str, Dummy, symbols, Symbol
- from sympy.core.sympify import SympifyError, _sympify
- from sympy.external.gmpy import SYMPY_INTS
- from sympy.functions import conjugate, adjoint
- from sympy.functions.special.tensor_functions import KroneckerDelta
- from sympy.matrices.common import NonSquareMatrixError
- from sympy.matrices.matrices import MatrixKind, MatrixBase
- from sympy.multipledispatch import dispatch
- 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)
- return func(a, b)
- except SympifyError:
- return retval
- return __sympifyit_wrapper
- return deco
- class MatrixExpr(Expr):
- """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
- See Also
- ========
- MatrixSymbol, MatAdd, MatMul, Transpose, Inverse
- """
- # Should not be considered iterable by the
- # sympy.utilities.iterables.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 # type: bool
- is_MatrixExpr = True # type: bool
- is_Identity = None # type: FuzzyBool
- is_Inverse = False
- is_Transpose = False
- is_ZeroMatrix = False
- is_MatAdd = False
- is_MatMul = False
- is_commutative = False
- is_number = False
- is_symbol = False
- is_scalar = False
- kind: MatrixKind = MatrixKind()
- def __new__(cls, *args, **kwargs):
- args = map(_sympify, args)
- return Basic.__new__(cls, *args, **kwargs)
- # The following is adapted from the core Expr object
- @property
- def shape(self) -> tTuple[Expr, Expr]:
- raise NotImplementedError
- @property
- def _add_handler(self):
- return MatAdd
- @property
- def _mul_handler(self):
- return MatMul
- def __neg__(self):
- return MatMul(S.NegativeOne, self).doit()
- def __abs__(self):
- raise NotImplementedError
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__radd__')
- def __add__(self, other):
- return MatAdd(self, other, check=True).doit()
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__add__')
- def __radd__(self, other):
- return MatAdd(other, self, check=True).doit()
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__rsub__')
- def __sub__(self, other):
- return MatAdd(self, -other, check=True).doit()
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__sub__')
- def __rsub__(self, other):
- return MatAdd(other, -self, check=True).doit()
- @_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):
- return MatPow(self, other).doit()
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__pow__')
- def __rpow__(self, other):
- raise NotImplementedError("Matrix Power not defined")
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__rtruediv__')
- def __truediv__(self, other):
- return self * other**S.NegativeOne
- @_sympifyit('other', NotImplemented)
- @call_highest_priority('__truediv__')
- def __rtruediv__(self, other):
- raise NotImplementedError()
- #return MatMul(other, Pow(self, S.NegativeOne))
- @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.adjoint import Adjoint
- return Adjoint(Transpose(self))
- def as_real_imag(self, deep=True, **hints):
- real = S.Half * (self + self._eval_conjugate())
- im = (self - self._eval_conjugate())/(2*S.ImaginaryUnit)
- return (real, im)
- def _eval_inverse(self):
- return Inverse(self)
- def _eval_determinant(self):
- return Determinant(self)
- def _eval_transpose(self):
- return Transpose(self)
- def _eval_power(self, exp):
- """
- Override this in sub-classes to implement simplification of powers. The cases where the exponent
- is -1, 0, 1 are already covered in MatPow.doit(), so implementations can exclude these cases.
- """
- return MatPow(self, exp)
- def _eval_simplify(self, **kwargs):
- if self.is_Atom:
- return self
- else:
- return self.func(*[simplify(x, **kwargs) for x in self.args])
- def _eval_adjoint(self):
- from sympy.matrices.expressions.adjoint import Adjoint
- return Adjoint(self)
- def _eval_derivative_n_times(self, x, n):
- return Basic._eval_derivative_n_times(self, x, n)
- def _eval_derivative(self, x):
- # `x` is a scalar:
- if self.has(x):
- # See if there are other methods using it:
- return super()._eval_derivative(x)
- else:
- return ZeroMatrix(*self.shape)
- @classmethod
- def _check_dim(cls, dim):
- """Helper function to check invalid matrix dimensions"""
- ok = check_assumptions(dim, integer=True, nonnegative=True)
- if ok is False:
- raise ValueError(
- "The dimension specification {} should be "
- "a nonnegative integer.".format(dim))
- def _entry(self, i, j, **kwargs):
- raise NotImplementedError(
- "Indexing not implemented for %s" % self.__class__.__name__)
- def adjoint(self):
- return adjoint(self)
- 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)
- @property
- def T(self):
- '''Matrix transposition'''
- return self.transpose()
- def inverse(self):
- if not self.is_square:
- raise NonSquareMatrixError('Inverse of non-square matrix')
- return self._eval_inverse()
- def inv(self):
- return self.inverse()
- def det(self):
- from sympy.matrices.expressions.determinant import det
- return det(self)
- @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, (SYMPY_INTS, 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)
- def _is_shape_symbolic(self) -> bool:
- return (not isinstance(self.rows, (SYMPY_INTS, Integer))
- or not isinstance(self.cols, (SYMPY_INTS, Integer)))
- 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]])
- See Also
- ========
- as_mutable: returns mutable Matrix type
- """
- if self._is_shape_symbolic():
- raise ValueError(
- 'Matrix with symbolic shape '
- 'cannot be represented explicitly.')
- from sympy.matrices.immutable import ImmutableDenseMatrix
- return ImmutableDenseMatrix([[self[i, j]
- for j in range(self.cols)]
- for i in range(self.rows)])
- 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]])
- See Also
- ========
- 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
- 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)
- @staticmethod
- def from_index_summation(expr, first_index=None, last_index=None, dimensions=None):
- r"""
- Parse expression of matrices with explicitly summed indices into a
- matrix expression without indices, if possible.
- This transformation expressed in mathematical notation:
- `\sum_{j=0}^{N-1} A_{i,j} B_{j,k} \Longrightarrow \mathbf{A}\cdot \mathbf{B}`
- Optional parameter ``first_index``: specify which free index to use as
- the index starting the expression.
- Examples
- ========
- >>> from sympy import MatrixSymbol, MatrixExpr, Sum
- >>> from sympy.abc import i, j, k, l, N
- >>> A = MatrixSymbol("A", N, N)
- >>> B = MatrixSymbol("B", N, N)
- >>> expr = Sum(A[i, j]*B[j, k], (j, 0, N-1))
- >>> MatrixExpr.from_index_summation(expr)
- A*B
- Transposition is detected:
- >>> expr = Sum(A[j, i]*B[j, k], (j, 0, N-1))
- >>> MatrixExpr.from_index_summation(expr)
- A.T*B
- Detect the trace:
- >>> expr = Sum(A[i, i], (i, 0, N-1))
- >>> MatrixExpr.from_index_summation(expr)
- Trace(A)
- More complicated expressions:
- >>> expr = Sum(A[i, j]*B[k, j]*A[l, k], (j, 0, N-1), (k, 0, N-1))
- >>> MatrixExpr.from_index_summation(expr)
- A*B.T*A.T
- """
- from sympy.tensor.array.expressions.conv_indexed_to_array import convert_indexed_to_array
- from sympy.tensor.array.expressions.conv_array_to_matrix import convert_array_to_matrix
- first_indices = []
- if first_index is not None:
- first_indices.append(first_index)
- if last_index is not None:
- first_indices.append(last_index)
- arr = convert_indexed_to_array(expr, first_indices=first_indices)
- return convert_array_to_matrix(arr)
- def applyfunc(self, func):
- from .applyfunc import ElementwiseApplyFunction
- return ElementwiseApplyFunction(func, self)
- @dispatch(MatrixExpr, Expr)
- def _eval_is_eq(lhs, rhs): # noqa:F811
- return False
- @dispatch(MatrixExpr, MatrixExpr) # type: ignore
- def _eval_is_eq(lhs, rhs): # noqa:F811
- if lhs.shape != rhs.shape:
- return False
- if (lhs - rhs).is_ZeroMatrix:
- return True
- def get_postprocessor(cls):
- def _postprocessor(expr):
- # To avoid circular imports, we can't have MatMul/MatAdd on the top level
- mat_class = {Mul: MatMul, Add: MatAdd}[cls]
- nonmatrices = []
- matrices = []
- for term in expr.args:
- if isinstance(term, MatrixExpr):
- matrices.append(term)
- else:
- nonmatrices.append(term)
- if not matrices:
- return cls._from_args(nonmatrices)
- if nonmatrices:
- if cls == Mul:
- for i in range(len(matrices)):
- if not matrices[i].is_MatrixExpr:
- # If one of the matrices explicit, absorb the scalar into it
- # (doit will combine all explicit matrices into one, so it
- # doesn't matter which)
- matrices[i] = matrices[i].__mul__(cls._from_args(nonmatrices))
- nonmatrices = []
- break
- else:
- # Maintain the ability to create Add(scalar, matrix) without
- # raising an exception. That way different algorithms can
- # replace matrix expressions with non-commutative symbols to
- # manipulate them like non-commutative scalars.
- return cls._from_args(nonmatrices + [mat_class(*matrices).doit(deep=False)])
- if mat_class == MatAdd:
- return mat_class(*matrices).doit(deep=False)
- return mat_class(cls._from_args(nonmatrices), *matrices).doit(deep=False)
- return _postprocessor
- Basic._constructor_postprocessor_mapping[MatrixExpr] = {
- "Mul": [get_postprocessor(Mul)],
- "Add": [get_postprocessor(Add)],
- }
- def _matrix_derivative(expr, x, old_algorithm=False):
- if isinstance(expr, MatrixBase) or isinstance(x, MatrixBase):
- # Do not use array expressions for explicit matrices:
- old_algorithm = True
- if old_algorithm:
- return _matrix_derivative_old_algorithm(expr, x)
- from sympy.tensor.array.expressions.conv_matrix_to_array import convert_matrix_to_array
- from sympy.tensor.array.expressions.arrayexpr_derivatives import array_derive
- from sympy.tensor.array.expressions.conv_array_to_matrix import convert_array_to_matrix
- array_expr = convert_matrix_to_array(expr)
- diff_array_expr = array_derive(array_expr, x)
- diff_matrix_expr = convert_array_to_matrix(diff_array_expr)
- return diff_matrix_expr
- def _matrix_derivative_old_algorithm(expr, x):
- from sympy.tensor.array.array_derivatives import ArrayDerivative
- lines = expr._eval_derivative_matrix_lines(x)
- parts = [i.build() for i in lines]
- from sympy.tensor.array.expressions.conv_array_to_matrix import convert_array_to_matrix
- parts = [[convert_array_to_matrix(j) for j in i] for i in parts]
- def _get_shape(elem):
- if isinstance(elem, MatrixExpr):
- return elem.shape
- return 1, 1
- def get_rank(parts):
- return sum([j not in (1, None) for i in parts for j in _get_shape(i)])
- ranks = [get_rank(i) for i in parts]
- rank = ranks[0]
- def contract_one_dims(parts):
- if len(parts) == 1:
- return parts[0]
- else:
- p1, p2 = parts[:2]
- if p2.is_Matrix:
- p2 = p2.T
- if p1 == Identity(1):
- pbase = p2
- elif p2 == Identity(1):
- pbase = p1
- else:
- pbase = p1*p2
- if len(parts) == 2:
- return pbase
- else: # len(parts) > 2
- if pbase.is_Matrix:
- raise ValueError("")
- return pbase*Mul.fromiter(parts[2:])
- if rank <= 2:
- return Add.fromiter([contract_one_dims(i) for i in parts])
- return ArrayDerivative(expr, x)
- 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.matrices.matrices import MatrixBase
- if isinstance(name, (MatrixBase,)):
- if n.is_Integer and m.is_Integer:
- return name[n, m]
- if isinstance(name, str):
- name = Symbol(name)
- else:
- name = _sympify(name)
- if not isinstance(name.kind, MatrixKind):
- raise TypeError("First argument of MatrixElement should be a matrix")
- 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]]
- @property
- def indices(self):
- return self.args[1:]
- def _eval_derivative(self, v):
- if not isinstance(v, MatrixElement):
- from sympy.matrices.matrices import MatrixBase
- if isinstance(self.parent, MatrixBase):
- return self.parent.diff(v)[self.i, self.j]
- return S.Zero
- M = self.args[0]
- m, n = self.parent.shape
- if M == v.args[0]:
- return KroneckerDelta(self.args[1], v.args[1], (0, m-1)) * \
- KroneckerDelta(self.args[2], v.args[2], (0, n-1))
- if isinstance(M, Inverse):
- from sympy.concrete.summations import Sum
- i, j = self.args[1:]
- i1, i2 = symbols("z1, z2", cls=Dummy)
- Y = M.args[0]
- r1, r2 = Y.shape
- return -Sum(M[i, i1]*Y[i1, i2].diff(v)*M[i2, j], (i1, 0, r1-1), (i2, 0, r2-1))
- if self.has(v.args[0]):
- return None
- return S.Zero
- 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
- Examples
- ========
- >>> 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
- is_symbol = True
- _diff_wrt = True
- def __new__(cls, name, n, m):
- n, m = _sympify(n), _sympify(m)
- cls._check_dim(m)
- cls._check_dim(n)
- if isinstance(name, str):
- name = Str(name)
- obj = Basic.__new__(cls, name, n, m)
- return obj
- @property
- def shape(self):
- return self.args[1], self.args[2]
- @property
- def name(self):
- return self.args[0].name
- def _entry(self, i, j, **kwargs):
- return MatrixElement(self, i, j)
- @property
- def free_symbols(self):
- return {self}
- def _eval_simplify(self, **kwargs):
- return self
- def _eval_derivative(self, x):
- # x is a scalar:
- return ZeroMatrix(self.shape[0], self.shape[1])
- def _eval_derivative_matrix_lines(self, x):
- if self != x:
- first = ZeroMatrix(x.shape[0], self.shape[0]) if self.shape[0] != 1 else S.Zero
- second = ZeroMatrix(x.shape[1], self.shape[1]) if self.shape[1] != 1 else S.Zero
- return [_LeftRightArgs(
- [first, second],
- )]
- else:
- first = Identity(self.shape[0]) if self.shape[0] != 1 else S.One
- second = Identity(self.shape[1]) if self.shape[1] != 1 else S.One
- return [_LeftRightArgs(
- [first, second],
- )]
- def matrix_symbols(expr):
- return [sym for sym in expr.free_symbols if sym.is_Matrix]
- class _LeftRightArgs:
- r"""
- Helper class to compute matrix derivatives.
- The logic: when an expression is derived by a matrix `X_{mn}`, two lines of
- matrix multiplications are created: the one contracted to `m` (first line),
- and the one contracted to `n` (second line).
- Transposition flips the side by which new matrices are connected to the
- lines.
- The trace connects the end of the two lines.
- """
- def __init__(self, lines, higher=S.One):
- self._lines = [i for i in lines]
- self._first_pointer_parent = self._lines
- self._first_pointer_index = 0
- self._first_line_index = 0
- self._second_pointer_parent = self._lines
- self._second_pointer_index = 1
- self._second_line_index = 1
- self.higher = higher
- @property
- def first_pointer(self):
- return self._first_pointer_parent[self._first_pointer_index]
- @first_pointer.setter
- def first_pointer(self, value):
- self._first_pointer_parent[self._first_pointer_index] = value
- @property
- def second_pointer(self):
- return self._second_pointer_parent[self._second_pointer_index]
- @second_pointer.setter
- def second_pointer(self, value):
- self._second_pointer_parent[self._second_pointer_index] = value
- def __repr__(self):
- built = [self._build(i) for i in self._lines]
- return "_LeftRightArgs(lines=%s, higher=%s)" % (
- built,
- self.higher,
- )
- def transpose(self):
- self._first_pointer_parent, self._second_pointer_parent = self._second_pointer_parent, self._first_pointer_parent
- self._first_pointer_index, self._second_pointer_index = self._second_pointer_index, self._first_pointer_index
- self._first_line_index, self._second_line_index = self._second_line_index, self._first_line_index
- return self
- @staticmethod
- def _build(expr):
- if isinstance(expr, ExprBuilder):
- return expr.build()
- if isinstance(expr, list):
- if len(expr) == 1:
- return expr[0]
- else:
- return expr[0](*[_LeftRightArgs._build(i) for i in expr[1]])
- else:
- return expr
- def build(self):
- data = [self._build(i) for i in self._lines]
- if self.higher != 1:
- data += [self._build(self.higher)]
- data = [i for i in data]
- return data
- def matrix_form(self):
- if self.first != 1 and self.higher != 1:
- raise ValueError("higher dimensional array cannot be represented")
- def _get_shape(elem):
- if isinstance(elem, MatrixExpr):
- return elem.shape
- return (None, None)
- if _get_shape(self.first)[1] != _get_shape(self.second)[1]:
- # Remove one-dimensional identity matrices:
- # (this is needed by `a.diff(a)` where `a` is a vector)
- if _get_shape(self.second) == (1, 1):
- return self.first*self.second[0, 0]
- if _get_shape(self.first) == (1, 1):
- return self.first[1, 1]*self.second.T
- raise ValueError("incompatible shapes")
- if self.first != 1:
- return self.first*self.second.T
- else:
- return self.higher
- def rank(self):
- """
- Number of dimensions different from trivial (warning: not related to
- matrix rank).
- """
- rank = 0
- if self.first != 1:
- rank += sum([i != 1 for i in self.first.shape])
- if self.second != 1:
- rank += sum([i != 1 for i in self.second.shape])
- if self.higher != 1:
- rank += 2
- return rank
- def _multiply_pointer(self, pointer, other):
- from ...tensor.array.expressions.array_expressions import ArrayTensorProduct
- from ...tensor.array.expressions.array_expressions import ArrayContraction
- subexpr = ExprBuilder(
- ArrayContraction,
- [
- ExprBuilder(
- ArrayTensorProduct,
- [
- pointer,
- other
- ]
- ),
- (1, 2)
- ],
- validator=ArrayContraction._validate
- )
- return subexpr
- def append_first(self, other):
- self.first_pointer *= other
- def append_second(self, other):
- self.second_pointer *= other
- def _make_matrix(x):
- from sympy.matrices.immutable import ImmutableDenseMatrix
- if isinstance(x, MatrixExpr):
- return x
- return ImmutableDenseMatrix([[x]])
- from .matmul import MatMul
- from .matadd import MatAdd
- from .matpow import MatPow
- from .transpose import Transpose
- from .inverse import Inverse
- from .special import ZeroMatrix, Identity
- from .determinant import Determinant
|