Opened 2 years ago

Last modified 12 months ago

#19112 needs_review enhancement

Add a function "isometry" to the quadratic forms package.

Reported by: tgaona Owned by:
Priority: minor Milestone: sage-7.1
Component: quadratic forms Keywords: isometry
Cc: annahaensch, tgaona Merged in:
Authors: Tyler Gaona Reviewers:
Report Upstream: N/A Work issues:
Branch: u/tgaona/ticket/19112 (Commits) Commit: 2b394f41bd18eb22c6558a1441444692b79caec7
Dependencies: Stopgaps:

Description

Adds a function "isometry" that returns an isometry from one quadratic form to another, provided that it exists.

Change History (52)

comment:1 Changed 22 months ago by tgaona

  • Branch set to u/tgaona/ticket/19112

comment:2 Changed 22 months ago by tgaona

  • Commit set to 9e211fecd4b525ad30276d37d442a68ce9bec8f0
  • Status changed from new to needs_review

New commits:

0d32f4bImplements a method to compute isometries between positive definite quadratic forms.
9e211feAdds a flag parameter to short_vector_list_up_to_length to pass to PARI's

comment:3 Changed 22 months ago by jdemeyer

  • Milestone changed from sage-6.9 to sage-6.10
  • Status changed from needs_review to needs_info
  • Work issues needs implementation deleted

This looks like it reinvents the already-existing method is_globally_equivalent_to. What does your code add?

comment:4 Changed 22 months ago by tgaona

Hi Jeroen,

Sorry for the slow response. This method should be able to return isometries for forms over fields, where I believe is_globally_equivalent_to just works over the ring of integers. For example:

sage: Q = DiagonalQuadraticForm(QQ, [1, 5])
sage: F = QuadraticForm(QQ, 2, [1, 12, 81])
sage: Q.is_globally_equivalent_to(F)
False
sage: Q.is_rationally_isometric(F)
True
sage: T = Q.isometry(F)
sage: T
[  1  -2]
[  0 1/3]
sage: Q.Gram_matrix() == T.transpose() * F.Gram_matrix() * T
True

So Q is equivalent to F over the rationals, but is_globally_equivalent_to doesn't recognize this. This method was intended to complement the is_rationally_isometric method, so perhaps it should be refactored to work in a similar manner as is_globally_equivalent_to, i.e, adding an optional flag to is_rationally_isometric to return the transition matrix for the isometry.

comment:5 Changed 22 months ago by jdemeyer

Right, I think it's better to mention this more clearly in the documentation (I always get confused when quadratic-form people talk about isomorphic/isometric/equivalent, I never know what they mean).

comment:6 Changed 22 months ago by jdemeyer

That being said, I think you are really over-complicating the algorithm.

Suppose one of the forms is given by the diagonal 5*x0^2 + .... Then really all you need to do in the first step is to find a vector v of "length" 5. This can be done with some algorithm based on PARI's qfsolve(): the only issue is that qfsolve() finds a projective point, but you really want an affine point.

So I suggest to do the following:

  1. implement an affine version of qfsolve() to solve equations of the form Q(x) = C where Q is a quadratic form and C is a constant. Do this in a separate Sage ticket.
  2. use this new function to implement isometry on this ticket.

This will have several major advantages:

  1. your code will be a lot faster.
  2. it will no longer run forever if there is isometry (which is unacceptable).
  3. you can remove the assumption that your quadratic form is positive-definite.

comment:7 Changed 22 months ago by jdemeyer

Also, your branch is based on an old version of Sage. You should really develop from the latest beta version (currently 6.10.beta2).

comment:8 Changed 22 months ago by jdemeyer

Please follow http://doc.sagemath.org/html/en/developer/coding_basics.html#documentation-strings for the correct formatting of docstrings.

comment:9 Changed 22 months ago by jdemeyer

  • Status changed from needs_info to needs_work

comment:10 Changed 22 months ago by tgaona

I will make it more clear in the documentation that the method returns a transformation matrix, as opposed to saying it returns an isometry.

As to your second point, I agree that a cleaner, more efficient replacement for the first step of the algorithm is very desirable. However, I have looked at the source for PARI's qfsolve() and it isn't clear to me how to adapt it to return an affine vector x such that Q(x) = C. I would appreciate it if you could offer more guidance here, otherwise, I'm not sure how to proceed.

  1. it will no longer run forever if there is isometry (which is unacceptable).

I can fix this by throwing an exception when is_rationally_isometric() returns false for the two forms. The reason I didn't include this initially is that I was working off of Sage 6.8 but is_rationally_isometric() was added in 6.9.

I will fix the formatting in the docstring.

comment:11 Changed 22 months ago by jdemeyer

  • Dependencies set to #18669

Given that your new functionality relies on rational_diagonal_form(), it would be a lot better in fact to base your patch on top of #18669.

comment:12 Changed 22 months ago by jdemeyer

About solving quadratic forms:

Suppose you want to solve Q(x) = c. Consider the quadratic form Q(x) - c*z^2 = 0, where z is an extra variable. Find a solution (x,z) to this quadratic form using qfsolve().

Case 1: If z != 0, then Q(x/z) = c.

Case 2: We found a solution Q(x) = 0. Let e be any vector such that B(x,e) != 0, where B is the bilinear form corresponding to Q. To find e, just try all unit vectors (0,..0,1,0...0). Let a = (c - Q(e))/B(x,e) and let y = e + a*x. Then Q(y) = c.

It would be great if you could implement this on a separate ticket.

comment:13 Changed 22 months ago by tgaona

  • Dependencies changed from #18669 to #18669 #19533

Thanks for the clear explanation. I've opened up a separate ticket for this, and will address the various changes you've suggested.

comment:14 Changed 21 months ago by git

  • Commit changed from 9e211fecd4b525ad30276d37d442a68ce9bec8f0 to afcd21b014a9c543b02c0b957c6d2c028d60b17f

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

51ba9f8Adds isometry method for quadratic forms.
afcd21bMerge branch 'research' into ticket/19112

comment:15 Changed 21 months ago by git

  • Commit changed from afcd21b014a9c543b02c0b957c6d2c028d60b17f to 8a040f3e2bdc5b312bd14b1fe36e85dda87ead91

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

8a040f3Merge branch 'develop' into ticket/19112

comment:16 Changed 21 months ago by git

  • Commit changed from 8a040f3e2bdc5b312bd14b1fe36e85dda87ead91 to 8f8ce1575dcfb53fee524a9e546d297d22aa82a0

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

3ee588dUse PARI to diagonalize quadratic forms
2689325Merge tag '6.10.beta1' into t/18669/use_pari_to_compute_rational_diagonal_form__
45cbc2dFurther fixes to rational diagonal form
ee778c8Implements a method to compute isometries between positive definite quadratic forms.
27b2190Adds a flag parameter to short_vector_list_up_to_length to pass to PARI's
92de9efAdds isometry method for quadratic forms.
1b1a574Merge branch 'ticket/19533' into ticket/19112
2f8ebc7Bug fixes in 'isometry' method and refactoring to utilize the 'solve' method for quadratic forms
8f8ce15Merge branch 'u/tgaona/ticket/19112' of git://trac.sagemath.org/sage into ticket/19112

comment:17 Changed 21 months ago by tgaona

  • Status changed from needs_work to needs_review

comment:18 Changed 21 months ago by jdemeyer

Can you review #18669?

comment:19 Changed 21 months ago by jdemeyer

The patchbot complains about a TAB character.

comment:20 Changed 21 months ago by git

  • Commit changed from 8f8ce1575dcfb53fee524a9e546d297d22aa82a0 to d2b3dde11edb22ba199d4666b99ad753865c5322

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

d2b3ddeRemoved errant TAB character in quadratic_form.py

comment:21 Changed 21 months ago by jdemeyer

  • Dependencies #18669 #19533 deleted

comment:22 Changed 21 months ago by jdemeyer

I haven't looked closely at the code, but some clean-up should be done:

  1. the function diagonal_isometry should have doctests
  2. remove commented code like #print ...

comment:23 follow-up: Changed 21 months ago by jdemeyer

Why does the final term need to be a special case? I don't like special cases unless they are justified.

comment:24 Changed 21 months ago by jdemeyer

Are you sure it's worth to check is_diagonal? Is there significant gain?

Are you sure it's worth to check is_rationally_isometric? I would just try to compute the isometry and let it fail if it doesn't exist.

comment:25 Changed 21 months ago by jdemeyer

Do you really need the copy in the function isometry?

comment:26 Changed 21 months ago by jdemeyer

Can you explain this block of code:

+            # Find a vector w such that Q(v) = F(w) where v = [1, ..., 0]
+            # v, w = vectors_of_common_length_dev(Q, F, q_basis, f_basis, i)
+            v = vector([0] * (n - i))
+            index = 0;
+            while True:
+                v[index] = v[index] + 1
+                index = (index + 1) % (n - i)
+                c = Q(v)
+                try:
+                    w = F.solve(c)
+                    #print("Find vectors {0} and {1} such that Q(v) = F(w)".format(v, w))        
+                    if not zero_row(f_basis, w, i) and not zero_row(q_basis, v, i):
+                        break
+                except ArithmeticError:
+                    # No solution found, try another vector.
+                    pass

comment:27 Changed 21 months ago by jdemeyer

  • Status changed from needs_review to needs_work

comment:28 Changed 21 months ago by jdemeyer

Can you move the helper functions out of the diagonal_isometry function and also doctest them?

comment:29 Changed 21 months ago by git

  • Commit changed from d2b3dde11edb22ba199d4666b99ad753865c5322 to 4b22c5155459c95fe2fac874c49d8bc7843f08c5

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

4b22c51Added doctesting to helper methods of isometry

comment:30 Changed 21 months ago by git

  • Commit changed from 4b22c5155459c95fe2fac874c49d8bc7843f08c5 to 582ae4500eb3652446c6f829310298953f8dbd02

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

582ae45Minor bugfixes/refactoring on isometry

comment:31 in reply to: ↑ 23 ; follow-up: Changed 21 months ago by tgaona

Replying to jdemeyer:

Why does the final term need to be a special case? I don't like special cases unless they are justified.

It was originally necessary because the short_vector_list_up_to_length() function didn't work for 1-dimensional quadratic forms, and the final iteration of the algorithm operates on two 1-dimensional forms. However, since that's been replaced by the PARI method, it's no longer necessary to have it as a special case, so thanks for pointing that out.

Are you sure it's worth to check is_diagonal? Is there significant gain?

I suppose not, and simplifies the code to remove it, so I will.

Are you sure it's worth to check is_rationally_isometric?

I believe this is necessary. The loop that looks for two vectors such that the modified bases including them will be non-singular will try vectors indefinitely until it finds a satisfactory pair. If the forms aren't equivalent there's no guarantee such a pair will be found.

Do you really need the copy in the function isometry?

Nope. That's a relic from an early version when isometry and diagonal_isometry were one function. I'll remove it, thanks for catching that.

Can you explain this block of code:

+            # Find a vector w such that Q(v) = F(w) where v = [1, ..., 0]
+            # v, w = vectors_of_common_length_dev(Q, F, q_basis, f_basis, i)
+            v = vector([0] * (n - i))
+            index = 0;
+            while True:
+                v[index] = v[index] + 1
+                index = (index + 1) % (n - i)
+                c = Q(v)
+                try:
+                    w = F.solve(c)
+                    #print("Find vectors {0} and {1} such that Q(v) = F(w)".format(v, w))        
+                    if not zero_row(f_basis, w, i) and not zero_row(q_basis, v, i):
+                        break
+                except ArithmeticError:
+                    # No solution found, try another vector.
+                    pass

This block finds a pair of vectors v and w such that Q(v) == F(v). These vectors will represent a linear combination of the vectors in the basis for each quadratic form. It's necessary that modifying the bases to include these vectors not produce a basis whose matrix is singular. So essentially this loop just looks for a pair of vectors satsifying these properties. It starts with v = [1, 0, 0] (for a 3-dimensional form) and finds w by calling F.solve(Q(v)). If this pair doesn't work, it finds a new v and starts over. I get new v's by incrementing each term in the vector so the first few vectors that are generated are: [1, 0, 0], [1, 1, 0], [1, 1, 1], [2, 1, 1]...

Also, I realized that matrices have an is_singular function, so I'm getting rid of the zero_row function.


New commits:

582ae45Minor bugfixes/refactoring on isometry

comment:32 Changed 21 months ago by tgaona

  • Status changed from needs_work to needs_review

comment:33 Changed 20 months ago by jdemeyer

More comments later, but already this: you make changes to short_vector_list_up_to_length() which I think do not belong to this ticket. In fact, it looks more likely that they belong to #14868. So I suggest you undo those changes here and move them to a new branch on #14868.

comment:34 in reply to: ↑ 31 Changed 20 months ago by jdemeyer

Replying to tgaona:

This block finds a pair of vectors v and w such that Q(v) == F(v). These vectors will represent a linear combination of the vectors in the basis for each quadratic form. It's necessary that modifying the bases to include these vectors not produce a basis whose matrix is singular. So essentially this loop just looks for a pair of vectors satsifying these properties. It starts with v = [1, 0, 0] (for a 3-dimensional form) and finds w by calling F.solve(Q(v)). If this pair doesn't work, it finds a new v and starts over. I get new v's by incrementing each term in the vector so the first few vectors that are generated are: [1, 0, 0], [1, 1, 0], [1, 1, 1], [2, 1, 1]...

Sorry, I still don't understand the algorithm. You really need more comments in the code. In particular, you need to explain the purpose of the variables i, q_basis, f_basis, qb and fb.

My feeling is that the current algorithm is too complicated (I don't think you need to loop for v), but I cannot really say how to improve it since I don't understand it.

Last edited 20 months ago by jdemeyer (previous) (diff)

comment:35 Changed 20 months ago by jdemeyer

  • Status changed from needs_review to needs_info

comment:36 Changed 19 months ago by annahaensch

I have looked over tgaona's algorithm, and it makes sense to me. The block of code mentioned in comment 26 is a bit confusing, but I believe that it is necessary since its eventually necessary to invert, and therefore avoid the possibility of having a singular change of basis matrix.

Another source of confusion, I believe, is in the description of the helper function modify basis. I do believe that it does precisely what it needs to do in the context of the code, but the description is not quite correct, and a little bit misleading, since it is always a function acting on a submatrix of the original basis matrix. A better description might read "Given a lattice L with basis matrix M and a vector, v=(v_1,...,v_n) of length n, this function extends the basis {b_1,...,b_n} of an underlying nxn orthogonal component of L to contain the vector v_1b_1+...+v_nb_n."

The functions diagonal_isometry, compute_gram_matrix_from_basis, modify_basis, and graham_schmidt are all very particular to this parent function isometry, so (I believe) the standard Sage convention is to begin those with an underscore. Also, I believe a better name for the parent function would be explicit_isometry just to differentiate this from the binary function returning True or False.

Other than that, I'll concede that there may very well be a faster way to run this algorithm, but as tgaona has it, it's certainly correct.

comment:37 Changed 19 months ago by git

  • Commit changed from 582ae4500eb3652446c6f829310298953f8dbd02 to 17d12c6faa5f18a9f822e30239bf201e2fe76c5f

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

6d9be4fremoved changes to 'short_vector_list_up_to_length()'
17d12c6Prefixed helper functions with an underscore and reworded documentation.

comment:38 Changed 19 months ago by tgaona

  • Status changed from needs_info to needs_review

comment:39 Changed 19 months ago by jdemeyer

Sorry, but 34 hasn't been addressed yet. The code is simply too hard to understand and therefore impossible to review for me.

And why did you rename isometry -> explicit_isometry? That's much harder to find.

comment:40 Changed 19 months ago by jdemeyer

  • Status changed from needs_review to needs_work

Merge conflicts with Sage 7.0

comment:41 Changed 19 months ago by git

  • Commit changed from 17d12c6faa5f18a9f822e30239bf201e2fe76c5f to 56b5225a8caec33e4b4366a90d1686fc0a5f6768

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

56b5225Merge branch 'develop' of git://trac.sagemath.org/sage into ticket/19112

comment:42 Changed 18 months ago by git

  • Commit changed from 56b5225a8caec33e4b4366a90d1686fc0a5f6768 to a56cb960a859fd2cb318d00bc06eacebb8cccb91

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

132da46attempt to simplify algorithm
24e01feSimplified algorithm and code for 'isometry' and its helpers
a56cb96Merge branch 'u/tgaona/ticket/19112' of git://trac.sagemath.org/sage into isometry2

comment:43 Changed 18 months ago by git

  • Commit changed from a56cb960a859fd2cb318d00bc06eacebb8cccb91 to 642df1071e10e3a771abd7d02ceff017c0a3ef30

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

642df10Fixed typo in documentation example for _diagonal_isometry

comment:44 Changed 18 months ago by tgaona

  • Status changed from needs_work to needs_review

comment:45 Changed 18 months ago by jdemeyer

  • Milestone changed from sage-6.10 to sage-7.1
  • Status changed from needs_review to needs_work

Why this added line?

+            v = parilist[i]

comment:46 Changed 18 months ago by jdemeyer

I am still missing some high-level explanation of the algorithm. For example, what is FM?

comment:47 follow-up: Changed 17 months ago by annahaensch

I understand the algorithm, I don't think it's so complicated. It just looks for a vector w in W of length corresponding to, V[0][0], the upper left most entry in V. Then it simply extends the basis of W to contain w, preforms Gram Schmidt (fixing w), and then picks off <Q(w)> as an orthogonal component of W and <V[0][0]> as an orthogonal component of V. I feel satisfied the the algorithm does what it claims to do, and that it is mathematically correct.

A few notes on the documentation:

  1. In your documentation for _diagonal_isometry() you call the function isometry() (of course it works out to be the same thing here, since the forms are diagonal to begin with.
  1. Your documentation for _gram_schmidt() is a bit misleading, since it looks like you are treating the columns of a Gram matrix as a set of vectors to be orthogonalized, rather, do the process on some sort of (preferably basis) matrix.

In response to comment 45, I'm also not sure why v=parilist[i] got added. Is there a reason you're touching short_vector_list_up_to_length()? I think something may have happened on the merge.

comment:48 in reply to: ↑ 47 Changed 17 months ago by jdemeyer

Replying to annahaensch:

I understand the algorithm

The fact that you understand it, is not sufficient. It should be documented in the code such that everybody can understand it.

comment:49 Changed 17 months ago by git

  • Commit changed from 642df1071e10e3a771abd7d02ceff017c0a3ef30 to 70ed4293b680ec8859683a8e1c3b7791c1c67e73

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

2b8f694removed old line from quadratic_form__automorphisms
70ed429Fixed documentation errors and added more comments to _diagonal_isometry

comment:50 Changed 17 months ago by tgaona

  • Status changed from needs_work to needs_review

comment:51 Changed 12 months ago by git

  • Commit changed from 70ed4293b680ec8859683a8e1c3b7791c1c67e73 to 0b8841492ca2e5714879cadacb6076f595927bd7

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

44c29efMerge branch 'master' into isometry2
0b88414Changed 'from quadratic_form' to 'from sage.quadratic_forms.quadratic_form'

comment:52 Changed 12 months ago by git

  • Commit changed from 0b8841492ca2e5714879cadacb6076f595927bd7 to 2b394f41bd18eb22c6558a1441444692b79caec7

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

2b394f4Removed trailing whitespace.
Note: See TracTickets for help on using tickets.