Opened 15 months ago

Last modified 6 weeks ago

#30244 needs_work enhancement

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

Reported by: mkoeppe Owned by:
Priority: major Milestone: sage-9.5
Component: linear algebra Keywords:
Cc: tscrim, egourgoulhon Merged in:
Authors: Matthias Koeppe Reviewers:
Report Upstream: N/A Work issues:
Branch: u/mkoeppe/use__matmul__operator____ (Commits, GitHub, GitLab) Commit: da5104cc5c3b00d2ccc68b76b65c2d1c2c56fa98
Dependencies: Stopgaps:

Status badges

Description (last modified by mkoeppe)

#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

Change History (34)

comment:1 follow-up: 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: Changed 15 months ago by mkoeppe

Replying to nbruin:

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: 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

Replying to 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]).

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

comment:5 Changed 15 months ago by egourgoulhon

  • Cc egourgoulhon added

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: 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

Replying to 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:

ea74b6fAdd support for __matmul__ in the coercion model
d939b2cUpdate doctests for py3
e8d7924Merge branch 't/22760/add_support_for___matmul___in_the_coercion_model' into t/30244/use__matmul__operator____
a982dd8sage.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: 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

Replying to mkoeppe:

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

Replying to 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.

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

Replying to 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]).

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:

c08e5b5sage.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:

d263307IdentityMorphism: Override __matmul__, not __mul__
b9808f8PoorManMap: Add __matmul__, let __mul__ just call __matmul__
3efa472ActionEndomorphism: Override __matmul__, not __mul__
430dad5ConstructionFunctor: Switch to __matmul__ for functor composition; delegate from __mul__ to __matmul__
14ad8cdWordMorphism: Also define __matmul__
9c85802MatrixSpace._get_action_: Handle matmul for matrix-matrix actions
9d1dccdMatrix.__matmul__: New
e84a7dcMatrixMorphism_abstract: Define __matmul__, not __mul__; use @ in examples
66206f9FanMorphism: Define __matmul__, not __mul__
63b8892HyperbolicIsometry: 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:

9588f52FiniteSetEndoMap*: 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:

f626394TensorWithIndices: 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:

3006833FiniteSetEndoMap*: Define __matmul__, delegate to it from __mul__
de14c34TensorWithIndices: Make __matmul__ an alias of __mul__
451eb95Map.__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:

348f680TensorWithIndices: 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:

878e923CoercionModel.verify_coercion_maps: Use @ instead of * for composition
da5104cPrecomposedAction.__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.