Hi Ondrej, Thank you very much for the quick response.
I copied these two functions from my current diophantine.py and I forgot to copy the lines containing the imports. Sorry about that. I think both the outputs [2] and [3] are correct. In [3], since we called `descent(2, 1)` the equation is `w**2 = 2*x**2 + y**2`, and the output `(x, y, w) = (0, 1, 1)` satisfies the above equation. As I noticed, outputs generated by `descent()` when either `A` or `B` equals 1 are correct. Things go wrong when both A and B are different from 1. For example: >>>descent(2, 7) # equation w**2 = 2*x**2 + 7*y**2 (1, 7, 3) Here,` (x, y, w) = (1, 7, 3)` and it is wrong because `3**2 != 2 * (1**2) + 7 * (7**2)`. I'll send you a notebook as soon as possible. Thank you and best regards, Thilina. On Tue, Jul 16, 2013 at 2:06 AM, Ondřej Čertík <ondrej.cer...@gmail.com>wrote: > On Mon, Jul 15, 2013 at 2:31 PM, Ondřej Čertík <ondrej.cer...@gmail.com> > wrote: > > Thilina, > > > > Thanks. You can use ipython notebook for exactly this exploration (let > > me know if you have any problems setting it up). Here it is: > > > > http://nbviewer.ipython.org/6003151 > > P.S. One has to put > > from math import floor > > on top, and also the "divisors" function is not defined. Thilina, in > the meantime, you can send > a notebook that produces the wrong results and is executable and we > can then debug it. > > IPython notebook is very useful for such exploratory debugging. > > Ondrej > > > > > So the output [2] is correct, isn't it? But the output [3] seems > > incorrect to me, isn't it? Let me know which other output you found > > incorrect. > > > > I'll have a look at the book later. > > > > Ondrej > > > > On Mon, Jul 15, 2013 at 1:22 PM, Thilina Rathnayake > > <thilina.r...@gmail.com> wrote: > >> I attached the code herewith as a .py file so you can easily run it. > >> > >> > >> On Tue, Jul 16, 2013 at 12:47 AM, Thilina Rathnayake > >> <thilina.r...@gmail.com> wrote: > >>> > >>> Hi Ondrej, > >>> > >>> Thank you very much for the words of encouragement. I am implementing > the > >>> descent method. > >>> I have used the solver you pointed out to validate answers for > quadratic > >>> Diophantine equations. > >>> We need something which solve the the ternary quadratic equations. > >>> > >>> I will continue implementing the algorithm as you suggested. I think > >>> that's the best way. > >>> I am confident that we will find some way to validate the results. I > >>> couldn't completely check PARI/GP > >>> and sage. They might contain a solver for ternary quadratic forms. > >>> > >>> I implemented the descent algorithm but it doesn't return the correct > >>> results. Here is the code. > >>> > >>> def descent(A, B): > >>> """ > >>> Uses Lagrange's method to finda a solution to $w^2 = Ax^2 + By^2$. > >>> Output a tuple $(x_0, y_0, z_0)$ which a solution to the above > >>> equation. > >>> This solution is used as a base to calculate all the solutions. > >>> """ > >>> if abs(A) > abs(B): > >>> y, x, w = descent(B, A) > >>> return (x, y, w) > >>> > >>> if A == 1: > >>> return (1, 0, 1) > >>> > >>> if B == 1: > >>> return (0, 1, 1) > >>> > >>> if quadratic_congruence(A, B, 0.5) != None: > >>> r = quadratic_congruence(A, B, 0.5) > >>> Q = (r**2 - A)//B > >>> div = divisors(Q) > >>> > >>> B_0 = 0 > >>> > >>> for i in div: > >>> if isinstance(sqrt(Q // i), Integer): > >>> B_0, d = i, sqrt(Q // i) > >>> break > >>> > >>> if B_0 > 0: > >>> X, Y, W = descent(A, B_0) > >>> return (r*X + W, B*d*Y, A*X + r*W) > >>> > >>> def quadratic_congruence(a, m, k=1.0): > >>> """ > >>> Solves the quadratic congruence $x^2 \equiv a \ (mod \ m)$. Returns > >>> the > >>> first solution in the range $0 .. \lfloor k*m \rfloor$. > >>> Return None if solutions do not exist. Currently uses bruteforce. > >>> Good enough for ``m`` sufficiently small. > >>> > >>> TODO: An efficient algorithm should be implemented. > >>> """ > >>> upper_bound = floor(abs(k*m)) + 1 > >>> > >>> for i in range(upper_bound): > >>> if (i**2 - a)% abs(m) == 0: > >>> return i > >>> > >>> return None > >>> > >>> Perhaps you could take a look at it and tell me where I have done > wrong. I > >>> tried to follow > >>> almost the same symbols as in the algorithm. Here B_0 corresponds to > B' in > >>> the algorithm. > >>> Please take a look at it when you have some free time. > >>> > >>> Thank you and best regards, > >>> Thilina > >>> > >>> > >>> On Mon, Jul 15, 2013 at 9:01 PM, Ondřej Čertík < > ondrej.cer...@gmail.com> > >>> wrote: > >>>> > >>>> I see, it's the Chapter IV. > >>>> Are you going to implement the "Sieving method" (IV.3.1)? Or the > >>>> "descent method"? > >>>> > >>>> So, if I take the Pythagorean formula a^2 + b^2 = c^2, your new solver > >>>> will be able to give me all integer solutions to this? That would be > >>>> really cool. > >>>> > >>>> As far as a solver, I found this one > >>>> (http://www.alpertron.com.ar/QUAD.HTM). > >>>> But it can only do equations for "x" an "y". You need equations for > >>>> "x", "y" and "z". > >>>> > >>>> I would suggest that you implement the solver, and then we just check > >>>> the solutions numerically, that should be quite easy. We will not know > >>>> for sure that you got all the solutions, but just cite the algorithm > >>>> reference in the docstring, so that users know what it does. That > >>>> would be good start. If somebody later finds some other software that > >>>> can do this, or some other reference with examples and results, we can > >>>> make sure that it returns all the solutions. > >>>> > >>>> Ondrej > >>>> > >>>> On Sun, Jul 14, 2013 at 11:21 PM, Thilina Rathnayake > >>>> <thilina.r...@gmail.com> wrote: > >>>> > This is what I use. > >>>> > > >>>> > “The Algorithmic resolution of Diophantine Equations”, Nigel P. > Smart, > >>>> > Chapter IV, Quadratic ternary forms. > >>>> > > >>>> > Here is the book: > >>>> > > >>>> > > >>>> > > >>>> > On Mon, Jul 15, 2013 at 2:20 AM, Aaron Meurer <asmeu...@gmail.com> > >>>> > wrote: > >>>> >> > >>>> >> On Sun, Jul 14, 2013 at 12:09 AM, Thilina Rathnayake > >>>> >> <thilina.r...@gmail.com> wrote: > >>>> >> > Hi All, > >>>> >> > > >>>> >> > I am planning to implement ternary quadratic forms, i.e > equations of > >>>> >> > the > >>>> >> > form, > >>>> >> > a*x**2 + by**2 + cz**2 + fxy + gyz + hzx = 0. It would be better > If > >>>> >> > I > >>>> >> > can > >>>> >> > find > >>>> >> > a system which currently implement this so I can validate my > >>>> >> > results. If > >>>> >> > you > >>>> >> > know > >>>> >> > of any system which solves this or a source which have good > >>>> >> > literature > >>>> >> > on > >>>> >> > the problem, > >>>> >> > please let me know. > >>>> >> > >>>> >> What is the literature that you plan to use to get the algorithm in > >>>> >> the first place? > >>>> >> > >>>> >> Aaron Meurer > >>>> >> > >>>> >> > > >>>> >> > Also I have to solve the quadratic congruence x**2 = D (mod m) > as a > >>>> >> > sub > >>>> >> > problem to implement the algorithm I found. I found a few > algorithms > >>>> >> > on > >>>> >> > this > >>>> >> > but > >>>> >> > none of them explain precisely how to solve the general equation, > >>>> >> > all > >>>> >> > they > >>>> >> > do is > >>>> >> > solve the equation for m prime and gcd(D, m) = 1 and just ask to > use > >>>> >> > Chinese > >>>> >> > remainder theorem to combine the results to solve for the general > >>>> >> > case > >>>> >> > where > >>>> >> > m > >>>> >> > is not prime and gcd(D, m) not necessarily equal to 1. Any help > on > >>>> >> > this > >>>> >> > would > >>>> >> > be highly appreciated. > >>>> >> > > >>>> >> > Regards, > >>>> >> > Thilina > >>>> >> > > >>>> >> > -- > >>>> >> > You received this message because you are subscribed to the > Google > >>>> >> > Groups > >>>> >> > "sympy" group. > >>>> >> > To unsubscribe from this group and stop receiving emails from it, > >>>> >> > send > >>>> >> > an > >>>> >> > email to sympy+unsubscr...@googlegroups.com. > >>>> >> > To post to this group, send email to sympy@googlegroups.com. > >>>> >> > Visit this group at http://groups.google.com/group/sympy. > >>>> >> > For more options, visit https://groups.google.com/groups/opt_out > . > >>>> >> > > >>>> >> > > >>>> >> > >>>> >> -- > >>>> >> You received this message because you are subscribed to the Google > >>>> >> Groups > >>>> >> "sympy" group. > >>>> >> To unsubscribe from this group and stop receiving emails from it, > send > >>>> >> an > >>>> >> email to sympy+unsubscr...@googlegroups.com. > >>>> >> To post to this group, send email to sympy@googlegroups.com. > >>>> >> Visit this group at http://groups.google.com/group/sympy. > >>>> >> For more options, visit https://groups.google.com/groups/opt_out. > >>>> >> > >>>> >> > >>>> > > >>>> > -- > >>>> > You received this message because you are subscribed to the Google > >>>> > Groups > >>>> > "sympy" group. > >>>> > To unsubscribe from this group and stop receiving emails from it, > send > >>>> > an > >>>> > email to sympy+unsubscr...@googlegroups.com. > >>>> > To post to this group, send email to sympy@googlegroups.com. > >>>> > Visit this group at http://groups.google.com/group/sympy. > >>>> > For more options, visit https://groups.google.com/groups/opt_out. > >>>> > > >>>> > > >>>> > >>>> -- > >>>> You received this message because you are subscribed to the Google > Groups > >>>> "sympy" group. > >>>> To unsubscribe from this group and stop receiving emails from it, > send an > >>>> email to sympy+unsubscr...@googlegroups.com. > >>>> To post to this group, send email to sympy@googlegroups.com. > >>>> Visit this group at http://groups.google.com/group/sympy. > >>>> For more options, visit https://groups.google.com/groups/opt_out. > >>>> > >>>> > >>> > >> > >> -- > >> You received this message because you are subscribed to the Google > Groups > >> "sympy" group. > >> To unsubscribe from this group and stop receiving emails from it, send > an > >> email to sympy+unsubscr...@googlegroups.com. > >> To post to this group, send email to sympy@googlegroups.com. > >> Visit this group at http://groups.google.com/group/sympy. > >> For more options, visit https://groups.google.com/groups/opt_out. > >> > >> > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to sympy+unsubscr...@googlegroups.com. > To post to this group, send email to sympy@googlegroups.com. > Visit this group at http://groups.google.com/group/sympy. > For more options, visit https://groups.google.com/groups/opt_out. > > > -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to sympy+unsubscr...@googlegroups.com. To post to this group, send email to sympy@googlegroups.com. Visit this group at http://groups.google.com/group/sympy. For more options, visit https://groups.google.com/groups/opt_out.
