Opened 15 months ago

# Use _matmul_ operator (@) for matrix multiplication, map composition, tensor contraction

Reported by: Owned by: mkoeppe major sage-9.5 linear algebra tscrim, egourgoulhon Matthias Koeppe N/A u/mkoeppe/use__matmul__operator____ da5104cc5c3b00d2ccc68b76b65c2d1c2c56fa98

#22760 added support for __matmul__ in the coercion model.

We should start using it.

First step: Review the semantics of this operator in major Python software for matrix and tensor computation (NumPy, Numba, TensorFlow, PyTorch, ...) so that we do not paint ourselves into a corner.

Follow-up tickets:

• #32212 sage.geometry.hyperbolic_space

### comment:1 follow-up: ↓ 2 Changed 15 months ago by nbruin

I think in most cases we *don't* need it. In a matrix algebra, it's pretty clear that multiplication means matrix multiplication; not Hadamard product.

The first place I'd see a possible application is for self-maps of a ring, where one may want to consider both composition and point-wise products.

### comment:2 in reply to: ↑ 1 ; follow-up: ↓ 8 Changed 15 months ago by mkoeppe

The first place I'd see a possible application is for self-maps of a ring, where one may want to consider both composition and point-wise products.

+1. Wondering now if we need a category of @-itive semigroups.

### comment:3 follow-ups: ↓ 4 ↓ 20 Changed 15 months ago by jhpalmieri

I think using @ for tensor products would be nicer than tensor([M1,M2,M3]). Aside from the fact that it's shorter, Sage's current implementation is a bit of a mess. For example, I can never remember that it should be tensor([M1, M2, M3]), not tensor(M1, M2, M3), or you could use M1.tensor(M2, M3) but not M1.tensor([M2, M3]).

### comment:4 in reply to: ↑ 3 Changed 15 months ago by tscrim

I think using @ for tensor products would be nicer than tensor([M1,M2,M3]). Aside from the fact that it's shorter, Sage's current implementation is a bit of a mess. For example, I can never remember that it should be tensor([M1, M2, M3]), not tensor(M1, M2, M3), or you could use M1.tensor(M2, M3) but not M1.tensor([M2, M3]).

+1 for this. This would be a good place to start providing an implementation.

### comment:6 Changed 15 months ago by nbruin

I think '@' for tensor products would be very confusing for other uses in python. While it was introduced non-prescriptively as another binary operator, the implementation name __matmul__ suggests otherwise, and the numpy reason for requesting it was very solidly for matrix multiplication, i.e., map composition.

It would seem to me that the place where having operator notation for tensor products might be interesting is for Kronecker products of matrices. That's exactly the place where it directly clashes with numpy notation!! I also think it's very inefficient for multivariate Kronecker products, because the infix notation necessitates the construction of intermediate results.

For modules, I would expect tensor products and homs to be equally prevalent, so keeping a symmetry in their notation seems like a desirable thing to have.

So, I think "@" is a poor fit for tensor product.

### comment:7 follow-up: ↓ 16 Changed 15 months ago by mkoeppe

+1 on giving the "map composition" priority for this operator.

In the realm of tensors, matrix multiplication would generalize to tensor contraction, not tensor product.

### comment:8 in reply to: ↑ 2 Changed 15 months ago by mkoeppe

Wondering now if we need a category of @-itive semigroups.

Any thoughts on this? Do we need a category SemigroupsWithRespectToMatmul?

### comment:9 Changed 15 months ago by nbruin

I'd say no, for the same reason why we don't have this for "|", "&", "<<", ">>", and other binary operators that python provides. I think we first need a convincing use-case before we start building infrastructure. Code is a burden, not an asset :-).

### comment:10 Changed 15 months ago by mkoeppe

Note that these other binary operators do not actually participate in the coercion framework...

But great point, of course, that we should start with something concrete first.

### comment:11 Changed 15 months ago by mkoeppe

• Branch set to u/mkoeppe/use__matmul__operator____

### comment:12 Changed 15 months ago by nbruin

• Commit set to a982dd81253277c909d7f29b2e0615d5e90e9d86

I think you may want to reconsider #22760 in the light of possible usage scenarios. The coercion framework is particularly designed to figure out common parents into which the operands can be coerced so that the operation can be applied. This does not apply to all cases; for instance, for actions there is no appropriate *common* parent. For actions, other procedures are followed (and generally, less powerful ones; which is appropriate).

If @ is going to be composition, it's going to be mostly a partial operation if regarded as, say, an operation on homomorphisms between modules over a field. Alternatively, it's an operation that combines objects in DIFFERENT parents (e.g., a pairing $\mathrm{Hom}(A,B) \times \mathrm{Hom}(B,C) \to \mathrm{Hom(}A,C)$) in which case the coercion framework probably doesn't have an appropriate setting yet. I think you really want to know what the coercion framework is supposed to accomplish for you before you try and hook @ into it. If the only scenarios where things work are the cases where _matmul_ already knows what to do, then there's no benefit from an extra indirection layer: you could just put the login into __matmul__ directly.

New commits:

 ​ea74b6f Add support for __matmul__ in the coercion model ​d939b2c Update doctests for py3 ​e8d7924 Merge branch 't/22760/add_support_for___matmul___in_the_coercion_model' into t/30244/use__matmul__operator____ ​a982dd8 sage.categories.map.Map: Add __matmul__

### comment:13 Changed 15 months ago by mkoeppe

This is a great point.

On this branch, as you perhaps saw, I am already overriding the (double-underscore) __matmul__ operator, so the coercion framework is not actually involved at all when it comes to the @ operation between two Maps. So in this case, I don't think there is actually any overhead/indirection.

### comment:14 Changed 15 months ago by mkoeppe

But a concern is that by overriding it, it is disabling coerce actions...

### comment:15 follow-up: ↓ 17 Changed 15 months ago by tscrim

I don't think there is any harm in #22760 by adding the hook for it. However, I thought the coercion framework was designed to also handle actions. For example, you need coercion for z * A when z is an integer but A only has a QQ-action. I also thought the coercion framework was what called the different actions to see what was appropriate too.

That being said, I see your point about _matmul_ potentially not being so useful by itself. However, we might as well include it in case someone does have a use case for it.

### comment:16 in reply to: ↑ 7 Changed 15 months ago by egourgoulhon

In the realm of tensors, matrix multiplication would generalize to tensor contraction, not tensor product.

+1

### comment:17 in reply to: ↑ 15 Changed 15 months ago by nbruin

I don't think there is any harm in #22760 by adding the hook for it. However, I thought the coercion framework was designed to also handle actions. For example, you need coercion for z * A when z is an integer but A only has a QQ-action. I also thought the coercion framework was what called the different actions to see what was appropriate too.

Indeed, some coercion steps are possible for actions, but not nearly as much as for operations internal to structures. For instance, for addition between ZZ[x,y] and QQ[z], the coercion system will construct a common covering structure QQ[x,y,z] by combining a sequence of "construction functors" according to certain (heuristic!) rules. I don't think the rules for actions are nearly as advanced -- probably no coercions on the acted-upon set are tried at all; and probably shouldn't.

I'd hope there are better tools available for exploring what coercion to take than to query and see "what works". If that's actually what happens, then applying it to partial operators such as composition is definitely inappropriate.

In general, I'm not so sure coercion will help for @ if it goes the composition route, and then having the hook in the system is going to be counterproductive, because people will stumble on it and do unhelpful things with it. So I disagree with the idea that putting hooks just in case someone finds a use for it is harmless. Not making a design decision on it now can also mean not making a design mistake now. Avoiding mistakes has benefits.

### comment:18 Changed 15 months ago by mkoeppe

I think I will try out an alternative implementation of composition as a coerce-action (keeping the inherited __mul__ operator). Then we can experiment with this a little to gather more insights.

### comment:19 Changed 15 months ago by mkoeppe

• Summary changed from Use _matmul_ operator (@) to Use _matmul_ operator (@) for matrix multiplication, map composition, tensor contraction

### comment:20 in reply to: ↑ 3 Changed 15 months ago by mkoeppe

I think using @ for tensor products would be nicer than tensor([M1,M2,M3]). Aside from the fact that it's shorter, Sage's current implementation is a bit of a mess. For example, I can never remember that it should be tensor([M1, M2, M3]), not tensor(M1, M2, M3), or you could use M1.tensor(M2, M3) but not M1.tensor([M2, M3]).

As notation for tensor products, perhaps we can use a different operator ... how about the bitwise-and operator &?

### comment:21 Changed 15 months ago by mkoeppe

• Milestone changed from sage-9.2 to sage-9.3

### comment:22 Changed 8 months ago by mkoeppe

• Milestone changed from sage-9.3 to sage-9.4

Setting new milestone based on a cursory review of ticket status, priority, and last modification date.

### comment:23 Changed 3 months ago by git

• Commit changed from a982dd81253277c909d7f29b2e0615d5e90e9d86 to c08e5b5718a1f7a3e3ae2da578a92ab2c4848cf8

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

 ​c08e5b5 sage.categories.map.Map: Add __matmul__

### comment:24 Changed 3 months ago by git

• Commit changed from c08e5b5718a1f7a3e3ae2da578a92ab2c4848cf8 to 63b8892628388ebb0dd5e6b313b211b00c223ff0

Branch pushed to git repo; I updated commit sha1. Last 10 new commits:

 ​d263307 IdentityMorphism: Override __matmul__, not __mul__ ​b9808f8 PoorManMap: Add __matmul__, let __mul__ just call __matmul__ ​3efa472 ActionEndomorphism: Override __matmul__, not __mul__ ​430dad5 ConstructionFunctor: Switch to __matmul__ for functor composition; delegate from __mul__ to __matmul__ ​14ad8cd WordMorphism: Also define __matmul__ ​9c85802 MatrixSpace._get_action_: Handle matmul for matrix-matrix actions ​9d1dccd Matrix.__matmul__: New ​e84a7dc MatrixMorphism_abstract: Define __matmul__, not __mul__; use @ in examples ​66206f9 FanMorphism: Define __matmul__, not __mul__ ​63b8892 HyperbolicIsometry: Split out _composition from __mul__ so that also @ between isometries works

### comment:25 Changed 3 months ago by git

• Commit changed from 63b8892628388ebb0dd5e6b313b211b00c223ff0 to 9588f52a880461173a4c48475ba933cb74cdac90

Branch pushed to git repo; I updated commit sha1. New commits:

 ​9588f52 FiniteSetEndoMap*: Define __matmul__, delegate to it from __mul__

### comment:26 Changed 3 months ago by mkoeppe

• Authors set to Matthias Koeppe
• Status changed from new to needs_review

### comment:27 Changed 3 months ago by git

• Commit changed from 9588f52a880461173a4c48475ba933cb74cdac90 to f62639491c5605ba478791727f44bb3bca14ecdf

Branch pushed to git repo; I updated commit sha1. New commits:

 ​f626394 TensorWithIndices: Make __matmul__ an alias of __mul__

### comment:28 Changed 3 months ago by git

• Commit changed from f62639491c5605ba478791727f44bb3bca14ecdf to 451eb95ea1467a4f0ceb468e2c081686efccbcc5

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

 ​3006833 FiniteSetEndoMap*: Define __matmul__, delegate to it from __mul__ ​de14c34 TensorWithIndices: Make __matmul__ an alias of __mul__ ​451eb95 Map.__mul__: Add doctest output

### comment:29 Changed 3 months ago by git

• Commit changed from 451eb95ea1467a4f0ceb468e2c081686efccbcc5 to 348f6805894373fa8036afbd0db7423af127bb2d

Branch pushed to git repo; I updated commit sha1. New commits:

 ​348f680 TensorWithIndices: Update doctests

### comment:30 Changed 3 months ago by git

• Commit changed from 348f6805894373fa8036afbd0db7423af127bb2d to da5104cc5c3b00d2ccc68b76b65c2d1c2c56fa98

Branch pushed to git repo; I updated commit sha1. New commits:

 ​878e923 CoercionModel.verify_coercion_maps: Use @ instead of * for composition ​da5104c PrecomposedAction.__init__: Use @ instead of * for composition

### comment:31 Changed 3 months ago by mkoeppe

• Description modified (diff)

### comment:32 Changed 3 months ago by mkoeppe

Not sure if @ should also be used for matrix-vector multiplication. The current code on the branch does not do this; but scipy seems to think so (see for example https://docs.scipy.org/doc/scipy/reference/optimize.linprog-highs.html)

### comment:33 Changed 3 months ago by mkoeppe

• Milestone changed from sage-9.4 to sage-9.5

### comment:34 Changed 6 weeks ago by mkoeppe

• Status changed from needs_review to needs_work
Note: See TracTickets for help on using tickets.