Opened 2 years ago

Last modified 3 months ago

#29886 needs_work enhancement

Add generator of symmetric matrices

Reported by: gh-Ivo-Maffei Owned by:
Priority: major Milestone: sage-9.8
Component: linear algebra Keywords: symmetric_matrices matrix_spaces
Cc: Matthias Köppe, Travis Scrimshaw Merged in:
Authors: Ivo Maffei Reviewers: Samuel Lelièvre, ...
Report Upstream: N/A Work issues:
Branch: u/gh-Ivo-Maffei/symmetric_matrices (Commits, GitHub, GitLab) Commit: 08fc8931f350c764b29aa1f76ac2a504de7cac5b
Dependencies: Stopgaps:

Status badges

Description (last modified by Samuel Lelièvre)

We add a method to the MatrixSpace class that returns a generator over the symmetric matrices in the space. The "symmetry" is defined by the user.

Attachments (1)

traceback.txt (9.6 KB) - added by gh-Ivo-Maffei 2 years ago.
mutable matrices issue full traceback

Download all attachments as: .zip

Change History (40)

comment:1 Changed 2 years ago by Samuel Lelièvre

Remove space before closing parenthesis:

- def symmetric_matrices(self, f, g=None ):
+ def symmetric_matrices(self, f, g=None):

Typo: doesn't contains -> doesn't contain.

Use a space after each comma in code, and around + signs.

Use double backticks, e.g. ``f`` or ``ValueError`` for code formatting in docstrings.

Typo "funcition" -> function.

Use code formatting for None: if it is None -> if it is ``None``

- skew-symmetric matrices:
+ Skew-symmetric matrices:

Typo: MetricSpace -> MatrixSpace.

Hint: after making changes, rebuild Sage and test with

sage -t src/sage/matrix/matrix_space.py

before pushing. That would catch things such as MetricSpace.

Replace this doctest whose output is hard to parse:

sage: gen = MS.symmetric_matrices(f)
sage: for M in gen:
....:     print(M)

by this one whose output takes up less space and is easier to parse:

sage: list(MS.symmetric_matrices(f))

Likewise, replace this:

sage: gen = MS.symmetric_matrices(f)
sage: i = 0
sage: for M in gen:
....:     i+= 1
....:     print(M)
....:     if i == 3: break

with this more concise:

sage: gen = MS.symmetric_matrices(f)
sage: [next(gen) for _ in range(3)]

Remove spaces around M in def make_symmetric( M ):.

Simplify "can't have symmetric matrices if they are not square" to "symmetric matrices must be square" or "only square matrices can be symmetric".

Use a space after # for comments, for instance:

- #When the base ring is not finite,
+ # When the base ring is not finite,

For inline comments, double space before # and single space after. For instance:

matrix_entries = []  # entries of matrix with lower half 0

No need to comment #make M symmetric before make_symmetric(M), function name says it all.

Rather than importing

import sage.combinat.integer_vector

and using sage.combinat.integer_vector.IntegerVectors, instead import

from sage.combinat.integer_vector import IntegerVectors

and then use IntegerVectors.

Last edited 2 years ago by Samuel Lelièvre (previous) (diff)

comment:2 Changed 2 years ago by Samuel Lelièvre

The PEP-0008 Python style guide and PEP-0257 Docstring conventions are recommended reading.

I should read them again too, I realise I need a refresher on some points there.

comment:3 Changed 2 years ago by git

Commit: a24e6a07971f44715896523d15b8c61633ef9e529b5145b551a552ece2f31c0cf6a42b7c10d08f89

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

1b38de3fixed formatting
9b5145bremoved trailing whitespaces

comment:4 Changed 2 years ago by gh-Ivo-Maffei

Authors: gh-Ivo-MaffeiIvo Maffei
Status: newneeds_review

comment:5 in reply to:  1 Changed 2 years ago by gh-Ivo-Maffei

Thanks for the detailed review. I think I fixed all formatting issues now. I also took the liberty to change a bit the __iter__ method.

comment:6 Changed 2 years ago by Samuel Lelièvre

Description: modified (diff)
Summary: Generator of symmetric matrices.Add generator of symmetric matrices

Use is None to test if something is None:

- if g == None:
+ if g is None:

I would also use f=None to mean f is the identity -- giving the standard "symmetric" matrices:

if f is None:
    f = lambda x: x

I would write the input as:

- ``f`` -- function or ``None`` (optional; default: ``None``);
  if ``None``, use the identity function

- ``g`` -- function or ``None`` (optional; default: ``None``);
  if ``None``, use the same function as for ``f``

Need to decide what symmetric should mean

  • Among possible definitions:
    • subdiagonal entries are images by f of the corresponding superdiagonal entries, i.e. for any i < j we have A[j,i] = f(A[i,j])
    • nondiagonal entries are images by f of their mirror entries, i.e. for any i != j we have A[j,i] = f(A[i,j])
  • These are the same if f is an involution; if that is assumed it needs to be said
  • Accordingly, do all or only half the checks

I would avoid formulas in the docstring's first sentence and clarify in a second sentence:

Return a generator of all matrices in this matrix space
that are "symmetric" as defined by `f` and `g`.

This means subdiagonal entries are images by ``f``
of the corresponding superdiagonal entries, i.e.
for any ``i < j`` we have ``A[j,i] = f(A[i,j])``,
and diagonal elements are stable by ``g`` ,
i.e. for any ``i`` we have ``A[i,i] = g(A[i,i])``.

The default value of ``None`` for ``f`` means
``f`` is the identity, and the default value
of ``None`` for ``g`` means ``g`` is ``f``.

Need :: here to start the doctest block:

Skew-symmetric matrices::

I would use v rather than iv for the loop variable.

Comments could be more concise:

- # If the number of entries is zero, then just
- # yield the empty matrix in that case and return
+ # If no entries, yield empty matrix and return
  if number_of_entries == 0:
      yield self(0)
      return
- # When the base ring is not finite, then we should go
- # through and yield the matrices by "weight", which is
- # the total number of iterations that need to be done
- # on the base ring to reach the matrix.
+ # For infinite rings, yield matrices by "weight",
+ # i.e. by number of iterations needed to reach them.
- # In the finite case, we do a similar thing except that
- # instead of checking if the diagonal is correct after creating the vector
- # we can select all possible diagonal elements a priori
+ For finite rings, allow all ring elements on the diagonal a priori.

comment:7 Changed 2 years ago by git

Commit: 9b5145b551a552ece2f31c0cf6a42b7c10d08f8908fc8931f350c764b29aa1f76ac2a504de7cac5b

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

08fc893fixed more code formatting; allow f=None for nomal symmetric matrices

comment:8 Changed 2 years ago by gh-Ivo-Maffei

Thanks for the review!! I apologise for the delay. I got caught in other tickets.

comment:9 Changed 2 years ago by Dima Pasechnik

reviewer name?

comment:10 Changed 2 years ago by gh-Ivo-Maffei

@slelievre (Samuel Lelièvre) did some reviewing, but never green lighted this. They probably never felt sure about the code.

comment:11 Changed 2 years ago by Dima Pasechnik

Cc: Matthias Köppe added
Reviewers: Samuel Lelièvre

Matthias, does it fit into the tensors stuff you are working on?

comment:12 Changed 2 years ago by Dima Pasechnik

Reviewers: Samuel LelièvreSamuel Lelièvre, ...

comment:13 Changed 2 years ago by Matthias Köppe

Yes, I think it does. In #30229 I define free modules of tensors with prescribed index symmetries. This covers the case of symmetric and antisymmetric matrices, of course. There is a basis and so if the base ring is an enumerated set, so is the module.

comment:14 Changed 2 years ago by Dima Pasechnik

I think that skew-symmetric matrices are requested often enough for them to warrant a named function (which would just forward naturally to what you altready implemented).

comment:15 Changed 2 years ago by Matthias Köppe

I agree that it would be nice to have facilities in matrix space for symmetric and skew-symmetric matrices, but I think it would be better to define the subspaces/modules of these matrices than to define enumeration methods. Enumeration comes for free from our general facilities for enumerated sets.

comment:16 Changed 2 years ago by gh-Ivo-Maffei

So you suggest having a general SymmetricMatrixSpace class which extends MatrixSpace or having the "standard" SymmetricMatrixSpace and SkewSymmetricMatrixSpace classes and leaving everything else to an enumeration method?

comment:17 Changed 2 years ago by Matthias Köppe

Cc: Travis Scrimshaw added

I don't think new classes are really needed. Instead, just methods that construct submodules.

Unfortunately, currently MatrixSpace does not fully implement the category of modules with bases (see (sage.categories.modules_with_basis) that it claims to belong to.

sage: M = MatrixSpace(QQ, 3, 3)
sage: M in ModulesWithBasis(QQ)
True
sage: A = M([[0, 0, 1], [0, 0, 0], [0, 0, 0]])
sage: A
[0 0 1]
[0 0 0]
[0 0 0]
sage: M.submodule([A])
AttributeError: 'MatrixSpace_with_category' object has no attribute 'get_order'

I would suggest to start with fixing the implementation of this functionality of MatrixSpace.

Then in the next step, you can define the methods that construct the specific modules that you need.

comment:18 Changed 2 years ago by Dima Pasechnik

I am inclined to wave this through, as it's needed for Ivo's GSoC project, and the time is tight.

comment:19 Changed 2 years ago by Matthias Köppe

I don't object to adding these specific enumeration methods, but I think adding this interface with the generality of this f and g thing is not a good idea, because a future implementation along the lines that we just discussed will not be able to support general f and g.

Last edited 2 years ago by Matthias Köppe (previous) (diff)

comment:20 Changed 2 years ago by Matthias Köppe

See #30276 for the type of symmetry that we hope to support.

comment:21 Changed 2 years ago by Travis Scrimshaw

Status: needs_reviewneeds_work

I usually thing of the sub/superdiagonal matrices as just being the entries immediately off the diagonal, not in all of the lower/upper triangular part.

A lot of what this method does is really confusing to me. The algorithm should be explained in the docstring, in particular the order which it is iterating through everything.

I also somewhat agree with Matthias here. I would first expect this to return a subspace, which doesn't work for generic f and g. Plus the word "symmetric" seems misleading to me for generic f and g. Perhaps as a middle ground, we call this symmetric_matrices_iterator().

Comments from a quick look at the code:

I don't see why you need to do v2.clone() since you never modify it. You simply get a sublist that to reassign to v.

self(0) -> self.zero().

This is very strange: for _ in range(nrows):. I don't understand this code block.

comment:22 Changed 2 years ago by Dima Pasechnik

I suggest to implement here only the minimum needed for the GSoC tickets to be done (on a different branch, to preserve this, I suppose).

Changed 2 years ago by gh-Ivo-Maffei

Attachment: traceback.txt added

mutable matrices issue full traceback

comment:23 Changed 2 years ago by gh-Ivo-Maffei

I spent some time looking at the submodule method to see if there was some easy fix. I think I got 90% of the way. There were some issues in the method sage.categories.finite_dimensional_modules_with_basis.echelon_form.

  • src/sage/categories/finite_dimensional_modules_with_basis.py

    diff --git a/src/sage/categories/finite_dimensional_modules_with_basis.py b/src/sage/categories/finite_dimensional_modules_with_basis.py
    
    a b class FiniteDimensionalModulesWithBasis(CategoryWithAxiom_over_base_ring): 
    337337                sage: C.echelon_form([x[0] - x[1], 2*x[1] - 2*x[2], x[0] - x[2]])
    338338                [x[0] - x[2], x[1] - x[2]]
    339339            """
    340340            if order is not None:
    341                 order = self._compute_support_order(order)
     341                order = self._compute_support_order(elements, order)
     342
     343            # order may be an ordering of the support of elements
     344            # hence not an ordering of the whole basis as needed below
     345            # so we extend order to an ordering of the whole basis
     346            orderList = list(order)
     347            for k in self.basis().keys():
     348                if k not in orderList:
     349                    orderList.append(k)
     350
     351            full_order = tuple(orderList)
     352
     353           
    342354            from sage.matrix.constructor import matrix
    343             mat = matrix(self.base_ring(), [g._vector_(order=order) for g in elements])
     355            mat = matrix(self.base_ring(), [g._vector_(order=full_order) for g in elements])
    344356            # Echelonizing a matrix over a field returned the rref
    345357            if row_reduced and self.base_ring() not in Fields():
    346358                try:

After applying that patch there is an issue with mutable matrices not being hashable. The issue seems to lie in sage.misc.cachefunc and I don't have enough knowledge of python or the inner working of Sage to proceed.

Honestly, I had a huge amount of issues with matrices/vectors not being hashable, so I wonder why aren't they immutable by default.

comment:24 in reply to:  22 Changed 2 years ago by Matthias Köppe

Replying to dimpase:

I suggest to implement here only the minimum needed for the GSoC tickets to be done (on a different branch, to preserve this, I suppose).

I have opened #30334 for the MatrixSpace.submodule issue.

comment:25 Changed 2 years ago by Matthias Köppe

#30334 is not urgent.

Given that submodule does not work, my suggestion for the present ticket would be the following. First implement a method named perhaps MatrixSpace.symmetric_submodule_basis(phase=1). As in #30276, for phase=1, return a basis of the submodule of the symmetric matrices, and for phase=-1, skew-symmetric matrices. Error for all other values of phase.

comment:26 Changed 2 years ago by gh-Ivo-Maffei

I don't fully understand #30334. What I'm using this method for, is to generate skew-symmetric matrices with diagonal entries 0 (the standard definition A^T = -A would allow diagonal entries 1 in GF(2)) and Hermitian matrices. Can these be achieved for some phase?

comment:27 Changed 2 years ago by Matthias Köppe

That's important information - from looking at the code it was not clear to me what generality you need.

For Hermitian matrices - over what field?

comment:28 in reply to:  27 Changed 2 years ago by gh-Ivo-Maffei

Replying to mkoeppe:

For Hermitian matrices - over what field?

Finite fields GF(r^2). Practically only GF(4) and GF(9).

comment:29 Changed 2 years ago by Matthias Köppe

Unfortunately the Hermitian matrices over a field k do not form a subspace (over k)...

comment:30 Changed 2 years ago by gh-Ivo-Maffei

They form a subgroup and they seem quite common... If you think there is no place for a generator of some sort here, then I'll just have the hermitian matrices created where I use them, so that nothing is exposed to the user.

comment:31 Changed 2 years ago by Dima Pasechnik

well, Hermitian matrices form a subspace over the fixed (by the automorphism) subfield. Probably our matrices aren't flexible enough for this, though.

comment:32 in reply to:  22 Changed 2 years ago by Matthias Köppe

Replying to dimpase:

I suggest to implement here only the minimum needed for the GSoC tickets to be done (on a different branch, to preserve this, I suppose).

To move forward, I'd suggest to move this methods to the code where they are needed. It can be be refactored later when we figure out the right framework.

comment:33 Changed 2 years ago by Matthias Köppe

Milestone: sage-9.2sage-9.3

comment:34 Changed 22 months ago by Matthias Köppe

Milestone: sage-9.3sage-9.4

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

comment:35 Changed 17 months ago by Samuel Lelièvre

Submodules of matrix spaces are taken care of at #31995 (superseding #30334).

comment:36 Changed 17 months ago by Matthias Köppe

Milestone: sage-9.4sage-9.5

Setting a new milestone for this ticket based on a cursory review.

comment:37 Changed 12 months ago by Matthias Köppe

Milestone: sage-9.5sage-9.6

comment:38 Changed 8 months ago by Matthias Köppe

Milestone: sage-9.6sage-9.7

comment:39 Changed 3 months ago by Matthias Köppe

Milestone: sage-9.7sage-9.8
Note: See TracTickets for help on using tickets.