Improving the solving process using Cornacchia
Hi All, It has been sometime since I last wrote about my project. I spent most of my time improving the solving methods used in binary quadratic forms. I also did a fair bit of refactoring of the earlier work and bug fixing. The implementation of ternary quadratic forms is now almost done and once the pull request by Pernici gets merged, we can use his function for finding the solutions for the quadratic congruence and the time taken to find the solutions for ternary / binary quadratic forms will be reduced considerably.
Tweaks to the algorithms used in solving quadratic binary forms
I had a misunderstanding earlier that the binary quadratic form can be converted to the form only in the case and is not a perfect square. Here . It occurred to me only when I read  more carefully that any binary quadratic with and can be transformed to the latter form. With this new discovery, our treatment of the binary quadratic becomes more general and efficient. Now in both the cases is not a perfect square and , we transform the equation to the form . The former case was earlier addressed under solving generalized Pell equation as the equation reduces to a generalized Pell equation in that case. My earlier approach with the latter case was to use brute force. I managed to implement the Cornacchia’s algorithm for this case as suggested by Pernici. Cornacchia’s algorithm can be used to solve the equation with and . In our form under the case , and , so these two conditions are automatically fulfilled. Below is a rough outline of the Cornacchia’s algorithm. Refer  and  for more details.
Solving by Cornacchia’s method
1. Let be such that .
2. Construct a set such that . We can ommit solutions greater than .
3. For each set .
4. Define the sequences, and for .
5. Find the terms until the condition is met.
6. Set .
7. If and is a perfect square we have the solution and .
8. Add to a set .
9. When we are done with all contains all the solutions to the equation .
Cornacchia’s method finds solutions such that . So it won’t find solutions to the equation since the solutions of this equation is and . So we have to extract out square factors in and find solutions to the equation with this new value and reconstruct the solutions. For example in the above example, 2 is a square factor of 20 (i.e square of two divides 20) so, we find solutions to , i.e for the equation and then reconstruct the solutions for the original equation.
diophantine() and check_solutions()
I replaced the main routine of the module by a new method `diophantine()`. What the old method,`diop_solve()` did was, when given with an input equation, it tried to find the type of the equation using `classify_diop()` and called the appropriate solving function accordingly. So when given an equation like , diop_solve will identify this as a binary quadratic and try to solve this using methods for binary quadratics. But the best thing to do in this case is to factorize the above equation as and solve and separately. After that, we can construct the general solution by adding parameters. The general solution will be constructed after solving and will be added to the solution set for . This was done by introducing the new routine `diophantine()`. This is now the main routine for the module and when we call `diophantine(eq)`, first `diophantine()` tries to factor `eq` and then solve each factor using `diop_solve()`. Then it creates the general solution and that will be added to the solution set for `eq`. This is a very good approach and it is also a very common method used when we try to solve Diophantine equations manually.
I also unified various kinds of functions in `test_diophantine.py` which were used to find whether the solutions returned by a function satisfies the input equation. I wrote a new method `check_solutions()`. Aaron should be thanked for proposing such a method and now the `test_diophantine.py` looks nicer and whoever trying to hack it doesn’t need to study various kinds of solution checking functions. To be honest, I was running out of names for these functions because a new one had to be created for each equation type.
I spent a considerable time fixing the bugs in the module, especially towards the end of the last week. I had done some mistakes when transforming general ternary quadratic equation into the form . Also, in `diophantine()` method, I hadn’t check for the case when None is returned to indicate the insolubility of the equation. I also found a bug in the `factor_list()` function under Python3. Here is the issue I created for it: https://code.google.com/p/sympy/issues/detail?id=3968
Hmm, That is pretty much everything I did during the past weeks. I hope in the upcoming week I can start on solving extended Pythagorean equation.
 Solving the equation ax^2 + bxy + cy^2 + dx + ey + f = 0, John P.Robertson, May 8, 2003. http://www.jpr2718.org/ax2p.pdf
 Solving ax^2 + by^2 = m by Cornacchia’s method, Number theory web, http://www.numbertheory.org/php/cornacchia.html
 A. Nitaj, L’algorithme de Cornacchia, Expositiones Mathematicae 13 (1995), 358-365