Sage: Ticket #15245: Pfaffian of a skew-symmetric matrix
https://trac.sagemath.org/ticket/15245
<p>
I couldn't believe my eyes when I saw we don't have the Pfaffian implemented in Sage.
</p>
<p>
Attached is an implementation that computes it over any ring in the combinatorial way using perfect matchings. This is probably not an optimal algorithm (and I feel like a lot could already be gained by improving up the <a class="missing wiki">PerfectMatchings?</a>(n) iterator without even changing the algorithm) but it's enough for my combinatorial needs.
</p>
<p>
Other than this, the patch makes perfect matchings iterable (no, they weren't) and improves some docstrings related to rook polynomials. The <a class="closed ticket" href="https://trac.sagemath.org/ticket/14117" title="#14117: defect: Jordan normal form not computed for nilpotent matrix over rational ... (closed: fixed)">#14117</a> dependency is because both patches edit matrix2.py and would probably cause some fuzz.
</p>
<p>
Apply:
</p>
<ul><li><a class="attachment" href="https://trac.sagemath.org/attachment/ticket/15245/trac_15245-pfaffian-dg.patch" title="Attachment 'trac_15245-pfaffian-dg.patch' in Ticket #15245">trac_15245-pfaffian-dg.patch</a><a class="trac-rawlink" href="https://trac.sagemath.org/raw-attachment/ticket/15245/trac_15245-pfaffian-dg.patch" title="Download"></a>
</li></ul>en-usSagehttps://trac.sagemath.org/chrome/site/logo_sagemath_trac.png
https://trac.sagemath.org/ticket/15245
Trac 1.2Darij GrinbergTue, 01 Oct 2013 04:13:10 GMTstatus changed
https://trac.sagemath.org/ticket/15245#comment:1
https://trac.sagemath.org/ticket/15245#comment:1
<ul>
<li><strong>status</strong>
changed from <em>new</em> to <em>needs_review</em>
</li>
</ul>
TicketDarij GrinbergTue, 01 Oct 2013 04:17:12 GMTdependencies set
https://trac.sagemath.org/ticket/15245#comment:2
https://trac.sagemath.org/ticket/15245#comment:2
<ul>
<li><strong>dependencies</strong>
set to <em>#14117</em>
</li>
</ul>
TicketTravis ScrimshawWed, 02 Oct 2013 15:33:40 GMT
https://trac.sagemath.org/ticket/15245#comment:3
https://trac.sagemath.org/ticket/15245#comment:3
<p>
Hey Darij,
</p>
<p>
A few things I'd like to see changed:
</p>
<ul><li>Could you call <code>self.is_skew_symmetric()</code> instead of implementing your own test?
</li><li>Instead of doing things such as <code>3x3</code> in the examples block, could you put it in full latex <code>`3 \times 3`</code>?
</li><li>For the <code>algorithm</code> input, I'd rather see something like
<pre class="wiki">- ``algorithm`` -- string, the algorithm to use; currently the following
algorithms have been implemented:
* ``'definition'`` - using the definition given by perfect matchings
</pre>or a variant with a better name for the algorithm.
</li></ul><p>
Also since this file is so big and often subject to changes, it would be better to not include as many whitespace changes.
</p>
<p>
Thanks,<br />
Travis
</p>
TicketDarij GrinbergWed, 02 Oct 2013 17:24:12 GMT
https://trac.sagemath.org/ticket/15245#comment:4
https://trac.sagemath.org/ticket/15245#comment:4
<p>
Hi Travis,
</p>
<p>
thanks once again for the reviews! I fixed all of your issues apart from not using <code>is_skew_symmetric()</code> because that method doesn't check diagonal entries to be 0 (it only checks them to satisfy x = -x, which is not the same in characteristic 2):
</p>
<pre class="wiki">sage: M = Matrix(Zmod(2), [[0,1],[1,1]])
sage: M.is_skew_symmetric()
True
</pre><p>
If you ask me, this is a bug, but I don't want to grasp into yet another hornet's nest. But thanks for having me look at my check code again; it contained an ugly indentation error...
</p>
<p>
I was pretty sure that matrix algorithms are the least busy part of the code, given how long <a class="closed ticket" href="https://trac.sagemath.org/ticket/14117" title="#14117: defect: Jordan normal form not computed for nilpotent matrix over rational ... (closed: fixed)">#14117</a> took to be reviewed. In hindsight I shouldn't have done unrelated docstring fixes, but the rook stuff was just pretty close and caught my eyes.
</p>
<p>
Best regards,<br />
Darij
</p>
TicketDarij GrinbergWed, 02 Oct 2013 17:24:52 GMTstatus, description changed
https://trac.sagemath.org/ticket/15245#comment:5
https://trac.sagemath.org/ticket/15245#comment:5
<ul>
<li><strong>status</strong>
changed from <em>needs_review</em> to <em>positive_review</em>
</li>
<li><strong>description</strong>
modified (<a href="/ticket/15245?action=diff&version=5">diff</a>)
</li>
</ul>
TicketDarij GrinbergWed, 02 Oct 2013 17:25:14 GMT
https://trac.sagemath.org/ticket/15245#comment:6
https://trac.sagemath.org/ticket/15245#comment:6
<p>
for the patchbot:
</p>
<p>
apply trac_15245-pfaffian-dg.patch
</p>
TicketNils BruinWed, 02 Oct 2013 18:02:47 GMT
https://trac.sagemath.org/ticket/15245#comment:7
https://trac.sagemath.org/ticket/15245#comment:7
<p>
Replying to <a class="ticket" href="https://trac.sagemath.org/ticket/15245#comment:4" title="Comment 4">darij</a>:
</p>
<blockquote class="citation">
<p>
thanks once again for the reviews! I fixed all of your issues apart from not using <code>is_skew_symmetric()</code> because that method doesn't check diagonal entries to be 0 (it only checks them to satisfy x = -x, which is not the same in characteristic 2):
</p>
</blockquote>
<p>
I think that's the definition of skew symmetric/antisymmetric: that A<sup>T</sup>=-A. It just happens to be the case that the concepts "symmetric" and "antisymmetric" coincide in characteristic 2.
</p>
<p>
In other words: skew symmetric matrices don't have to have 0 on their diagonal. Note that the terminology comes from bilinear forms, where alternating means (v,v)=0, antisymmetric means (v,w)=-(w,v) and symmetric means (v,w)=(w,v). "alternating" is not "antisymmetric" in characteristic 2.
</p>
TicketDarij GrinbergWed, 02 Oct 2013 18:08:25 GMT
https://trac.sagemath.org/ticket/15245#comment:8
https://trac.sagemath.org/ticket/15245#comment:8
<p>
That's the hornet's nest I was talking about. With the definition you give, the Pfaffian lacks many of its nice properties like squaring to the determinant. It is the definition used on <a class="ext-link" href="http://en.wikipedia.org/wiki/Skew-symmetric_matrix"><span class="icon"></span>http://en.wikipedia.org/wiki/Skew-symmetric_matrix</a> but not the definition used on <a class="ext-link" href="http://en.wikipedia.org/wiki/Pfaffian"><span class="icon"></span>http://en.wikipedia.org/wiki/Pfaffian</a> . Since mathematicians can't agree, I figured it is easier to check for the kind of skew-symmetry needed in the definition of the Pfaffian rather than push that definition into the rest of the code. Is there a better solution?
</p>
TicketJeroen DemeyerSun, 06 Oct 2013 09:43:18 GMTstatus changed
https://trac.sagemath.org/ticket/15245#comment:9
https://trac.sagemath.org/ticket/15245#comment:9
<ul>
<li><strong>status</strong>
changed from <em>positive_review</em> to <em>needs_work</em>
</li>
</ul>
TicketJeroen DemeyerSun, 06 Oct 2013 09:43:42 GMTstatus changed
https://trac.sagemath.org/ticket/15245#comment:10
https://trac.sagemath.org/ticket/15245#comment:10
<ul>
<li><strong>status</strong>
changed from <em>needs_work</em> to <em>needs_review</em>
</li>
</ul>
<p>
darij: it looks like you reviewed your own patch. I'd say this still needs a "formal" review.
</p>
TicketDarij GrinbergSun, 06 Oct 2013 13:38:14 GMT
https://trac.sagemath.org/ticket/15245#comment:11
https://trac.sagemath.org/ticket/15245#comment:11
<p>
Oops -- I just realized what Travis gave were comments, not a review. Sorry!
</p>
TicketTravis ScrimshawSun, 06 Oct 2013 15:02:29 GMT
https://trac.sagemath.org/ticket/15245#comment:12
https://trac.sagemath.org/ticket/15245#comment:12
<p>
Hey Darij,
</p>
<p>
Sorry for letting this slip away.
</p>
<p>
For the skew-symmetric, how about instead calling <code>is_skew_symmetric()</code> and then checking that the diagonal entries are 0 with
</p>
<pre class="wiki">all(d == 0 for d in self.diagonal())
</pre><p>
One more minor thing: the <code>AUTHORS:</code> block shouldn't be indented.
</p>
<p>
Best,<br />
Travis
</p>
TicketDarij GrinbergSun, 06 Oct 2013 16:58:27 GMT
https://trac.sagemath.org/ticket/15245#comment:13
https://trac.sagemath.org/ticket/15245#comment:13
<p>
I've given the <code>is_skew_symmetric</code> method an additional keyword variable now. All the rest is OK?
</p>
TicketNils BruinSun, 06 Oct 2013 17:24:12 GMT
https://trac.sagemath.org/ticket/15245#comment:14
https://trac.sagemath.org/ticket/15245#comment:14
<p>
Replying to <a class="ticket" href="https://trac.sagemath.org/ticket/15245#comment:13" title="Comment 13">darij</a>:
</p>
<blockquote class="citation">
<p>
I've given the <code>is_skew_symmetric</code> method an additional keyword variable now. All the rest is OK?
</p>
</blockquote>
<p>
That keyword should really be <code>is_alternating</code>, so the more appropriate thing would be to supply a method <code>is_alternating</code>.
</p>
<p>
If pfaffians need alternating rather than skew-symmetric (I haven't checked) then the confusion in terminology just comes from the fact that people looking at pfaffians haven't considered characteristic 2. That kind of thing happens all the time, and it's the kind of thing that computer algebra systems need to be a little more pedantic about than math literature, since you don't get to say "in this paper, with <THIS> we mean <SOMETHING ELSE>".
</p>
<p>
It may well be that pfaffians are simply not all that useful in characteristic 2, so that people didn't bother with them (Cayley certainly wouldn't have).
</p>
TicketDarij GrinbergSun, 06 Oct 2013 22:12:43 GMTattachment set
https://trac.sagemath.org/ticket/15245
https://trac.sagemath.org/ticket/15245
<ul>
<li><strong>attachment</strong>
set to <em>trac_15245-pfaffian-dg.patch</em>
</li>
</ul>
<p>
new version, separating skew-symmetry from alternatingness systematically
</p>
TicketDarij GrinbergSun, 06 Oct 2013 22:14:39 GMT
https://trac.sagemath.org/ticket/15245#comment:15
https://trac.sagemath.org/ticket/15245#comment:15
<p>
Done. A number of people seem to be lax about characteristic 2, among them Knuth from whom I had expected this the least.
</p>
TicketTravis ScrimshawMon, 07 Oct 2013 16:52:42 GMT
https://trac.sagemath.org/ticket/15245#comment:16
https://trac.sagemath.org/ticket/15245#comment:16
<p>
I'm happy with it (having an <code>is_alternating()</code> method). Nils?
</p>
TicketNils BruinMon, 07 Oct 2013 22:41:03 GMT
https://trac.sagemath.org/ticket/15245#comment:17
https://trac.sagemath.org/ticket/15245#comment:17
<p>
Replying to <a class="ticket" href="https://trac.sagemath.org/ticket/15245#comment:16" title="Comment 16">tscrim</a>:
</p>
<blockquote class="citation">
<p>
I'm happy with it (having an <code>is_alternating()</code> method). Nils?
</p>
</blockquote>
<p>
Yep. It's a pedantic difference, but since this gets exposed in our official matrix API I think it's worth being precise.
</p>
TicketTravis ScrimshawTue, 08 Oct 2013 14:31:44 GMTstatus changed; reviewer set
https://trac.sagemath.org/ticket/15245#comment:18
https://trac.sagemath.org/ticket/15245#comment:18
<ul>
<li><strong>status</strong>
changed from <em>needs_review</em> to <em>positive_review</em>
</li>
<li><strong>reviewer</strong>
set to <em>Travis Scrimshaw</em>
</li>
</ul>
<p>
Then it's a positive review.
</p>
TicketNils BruinTue, 08 Oct 2013 15:59:59 GMT
https://trac.sagemath.org/ticket/15245#comment:19
https://trac.sagemath.org/ticket/15245#comment:19
<p>
Oh, a comment that may be worthwhile for future work:
</p>
<p>
The generic implementations of is_skewsymmetric and is_alternating are probably horribly slow compared to what one can do on specific classes (see, e.g. <a class="closed ticket" href="https://trac.sagemath.org/ticket/15104" title="#15104: enhancement: Special case modn_dense matrix operations to improve performance (closed: fixed)">#15104</a>). Since the difference between is_alternating and is_skewsymmetric only is apparent in characteristic 2, we'd probably get better performance if one of the two calls the other if required. If the specific implementation only applies to characteristic not equal to 2, it only needs to implement one fast method and the generic other one will call it.
</p>
<p>
Probably <code>is_skewsymmetric</code> is the better choice for being the "main" method, because it's the more widely used term. We'd get something like:
</p>
<pre class="wiki"> def is_alternating(self):
if self.base_ring().characteristic() !=2:
return self.is_skewsymmetric()
<rest of code>
</pre>
TicketDarij GrinbergTue, 08 Oct 2013 16:38:17 GMT
https://trac.sagemath.org/ticket/15245#comment:20
https://trac.sagemath.org/ticket/15245#comment:20
<p>
Not sure about it. <code>is_skewsymmetric</code> should be a tad slower than <code>is_alternating</code> (not seriously so -- it just calls diagonal elements twice rather than once). And I never understood what a characteristic of a ring is; chances are this is another thing not consistently understood in Sage. What we probably cannot do is ask whether the characteristic is not 2; if anything, we should ask for 2 to be invertible. But even then, I fear that excepting the <code>NotImplementedError</code> errors will ruin the speed benefits we get from unifying the code.
</p>
<p>
Thanks, Travis and Nils, for the review and the helpful comments!
</p>
TicketNils BruinTue, 08 Oct 2013 16:55:35 GMT
https://trac.sagemath.org/ticket/15245#comment:21
https://trac.sagemath.org/ticket/15245#comment:21
<p>
Replying to <a class="ticket" href="https://trac.sagemath.org/ticket/15245#comment:20" title="Comment 20">darij</a>:
</p>
<blockquote class="citation">
<p>
Not sure about it. <code>is_skewsymmetric</code> should be a tad slower than <code>is_alternating</code> (not seriously so -- it just calls diagonal elements twice rather than once).
</p>
</blockquote>
<p>
If that matters people can always implement both.
</p>
<blockquote class="citation">
<p>
What we probably cannot do is ask whether the characteristic is not 2;
</p>
</blockquote>
<p>
uh ...
</p>
<pre class="wiki">sage: GF(2).characteristic() !=2
False
sage: ZZ.characteristic() != 2
True
</pre><blockquote class="citation">
<p>
if anything, we should ask for 2 to be invertible.
</p>
</blockquote>
<p>
No, that's not the correct check. In ZZ and (ZZ/6ZZ), 2 is not invertible and yet skew symmetric is the same as alternating.
</p>
<p>
You could check whether 1+1==0, but that's more expensive.
</p>
TicketTravis ScrimshawTue, 08 Oct 2013 18:03:16 GMT
https://trac.sagemath.org/ticket/15245#comment:22
https://trac.sagemath.org/ticket/15245#comment:22
<p>
I think what we really should be checking is if it has positive even characteristic (see the ZZ/4 example in the patch). For example doing something like
</p>
<div class="wiki-code"><div class="code"><pre><span class="k">def</span> <span class="nf">is_alternating</span><span class="p">(</span><span class="bp">self</span><span class="p">):</span>
<span class="k">if</span> <span class="ow">not</span> <span class="bp">self</span><span class="o">.</span>is_skew_symmetric<span class="p">():</span>
<span class="k">return</span> <span class="bp">False</span>
c <span class="o">=</span> <span class="bp">self</span><span class="o">.</span>base_ring<span class="p">()</span><span class="o">.</span>characteristic<span class="p">()</span>
<span class="k">return</span> c <span class="o">!=</span> <span class="mi">0</span> <span class="ow">or</span> c <span class="o">%</span> <span class="mi">2</span> <span class="o">!=</span> <span class="mi">0</span> \ <span class="c1"># If past here, we have pos. even char.</span>
<span class="ow">or</span> <span class="nb">all</span><span class="p">(</span><span class="bp">self</span><span class="o">.</span>get_unsafe<span class="p">(</span>i<span class="p">,</span>i<span class="p">)</span> <span class="o">==</span> <span class="mi">0</span> <span class="k">for</span> i <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span><span class="bp">self</span><span class="o">.</span>_ncols<span class="p">))</span>
</pre></div></div><p>
I do agree that <code>is_skew_symmetric()</code> should be the "main" method, also because it is a larger class of matrices in positive even characteristic. However, I don't see a way to speed up <code>is_skew_symmetric()</code> for the mod case, so I'm happy with the two ([more] optimized) implementations.
</p>
TicketDarij GrinbergTue, 08 Oct 2013 18:35:12 GMT
https://trac.sagemath.org/ticket/15245#comment:23
https://trac.sagemath.org/ticket/15245#comment:23
<p>
Skew symmetric is NOT the same as alternating in ZZ/(6 ZZ); think of a 3 on the main diagonal. Just having 1 + 1 != 0 is not enough. Infinite characteristic in the sense of "ZZ embeds into the ring" is not enough since we can have things like ZZ[X] / (2X).
</p>
TicketNils BruinTue, 08 Oct 2013 19:03:35 GMT
https://trac.sagemath.org/ticket/15245#comment:24
https://trac.sagemath.org/ticket/15245#comment:24
<p>
Replying to <a class="ticket" href="https://trac.sagemath.org/ticket/15245#comment:23" title="Comment 23">darij</a>:
</p>
<blockquote class="citation">
<p>
Skew symmetric is NOT the same as alternating in ZZ/(6 ZZ); think of a 3 on the main diagonal. Just having 1 + 1 != 0 is not enough. Infinite characteristic in the sense of "ZZ embeds into the ring" is not enough since we can have things like ZZ[X] / (2X).
</p>
</blockquote>
<p>
Oops, sorry. I should use non-integral domains a little more. I withdraw my proposal then. Deciding whether alternating is the same as skew symmetric is not so easy after all, so we should just have two separate methods and if people want the tests to be quicker they will have to implement both themselves.
</p>
TicketJeroen DemeyerSat, 12 Oct 2013 09:47:54 GMTstatus changed; resolution, merged set
https://trac.sagemath.org/ticket/15245#comment:25
https://trac.sagemath.org/ticket/15245#comment:25
<ul>
<li><strong>status</strong>
changed from <em>positive_review</em> to <em>closed</em>
</li>
<li><strong>resolution</strong>
set to <em>fixed</em>
</li>
<li><strong>merged</strong>
set to <em>sage-5.13.beta1</em>
</li>
</ul>
Ticket