4
$\begingroup$

Mathematica can solve $v^2- d u^2=4$ quickly, even for nonsmall $d$:

d = 400004;
Reduce[v^2 - d u^2 == 4, Integers]
d = 400012;
Reduce[v^2 - d u^2 == 4, Integers]

My difficulty is that I need the smallest solution with positive $u,v$, preferably in the form {u, v}. FindInstance gives one solution, but not necessarily the smallest.

Reduce[v^2 - d u^2 ==4 && 0 < u && 0 < v, Integers]

is also fast, but is there a fast & elegant way?

$\endgroup$
  • 1
    $\begingroup$ Is this post about the Project Euler Problem 66? $\endgroup$ – Αλέξανδρος Ζεγγ Jan 7 '18 at 11:13
  • 1
    $\begingroup$ [Taking the liberty to respond for Kevin] @AlexanderZeng No $\endgroup$ – Daniel Lichtblau Jan 7 '18 at 23:35
  • $\begingroup$ I need to examine a range of $d$ around a million. But now that you've pointed it out, I'll submit the solution to #66, too. $\endgroup$ – Kevin O'Bryant Jan 9 '18 at 17:59
  • $\begingroup$ For searching purposes: this is a generalized Pell equation. $\endgroup$ – J. M. will be back soon Mar 5 '18 at 0:08
7
$\begingroup$

You can obtain the smallest solution using the following

d = 400004;
b = Convergents[Sqrt[d]][[-2]];
{vS, uS} = 2*{Numerator[b], Denominator[b]}

{639176293975850025902542507002869107405272920896001249322708894389485\ 74385666998528790229868783788066961187502, \ 1010621404580133400642412209621549153797917953531796195330954316636633\ 23751876289858094985610380884355499500}

Check:

vS^2 - d uS^2

4

Reference: Olds, Carl Douglas. Continued fractions. Vol. 18. New York: Random House, 1963.

$\endgroup$
5
$\begingroup$

The answer by @AnjanKumar is very concise and very fast for $d=400004$; however, it may not be correct for all values of $d$. Some values of $d$ have more than one fundamental solution, that is, there are different classes of solution built on different fundamental solutions.

Consider $d=60$ and $n=4$ in $v^2-d*u^2=n$. The previous solution via convergents is $\{v,u\}=\{62,8\}$.

2*{Numerator[#], Denominator[#]} &[Convergents[Sqrt[60]][[-2]]]

{62,8}

However, the second class of solutions for $d=60$ has fundamental solution {8,1}, which is smaller.

8^2 - 60*1^2

4

For $d=400012$, my imperfect code (below) says there are 2 classes of solution (the first two of three returned values are identical), and the convergents method gives the largest of the corresponding fundamental solutions. But in fact, the smallest solution is

{797961680157890791240296184184442534003247415937056231677647802477835\ 5033476938301119995768929905092055262086898273030683968332913381379325\ 36366968953450845322267394081990452, 1261669272519709656748944247157199207977946870608786859077399143331190\ 3335367323625818174514114436565606506960785400148895680270988581555924\ 43412490263634202021762111820545}

I am not a mathematician and do not know how to fix the code based on Convergents. For large $d$, I have horrible code based on Reduce.

FermatReduce[d_Integer, n_Integer] :=
   Block[{x, y, a, b, m},
      If[VectorQ[#], {x -> #[[1]], y -> #[[2]]}, 
         Map[{x -> #[[1]], y -> #[[2]]} &, #]] &[
            DeleteCases[
               Simplify[
                  Apply[List, Reduce[{x^2 - d y^2 == n, x > 0, y > 0}, {x, y}, Integers] /.
                    {x == a_Integer && y == b_Integer -> {a, b},
                     C[1] \[Element] Integers && C[1] >= 0 && x == a_ && y == b_ :> ({a, b} /. C[1] -> m),
                     C[1] \[Element] Integers && C[1] >= 1 && x == a_ && y == b_ :> ({a, b} /. C[1] -> m + 1)}]],
            {x_Integer, y_Integer}]]]

FermatReduceN[d_Integer, n_Integer, mmax_: 1] :=
   With[{s = FermatReduce[d, n]},
        If[s =!= False, If[VectorQ[s], Transpose[Simplify[({x, y} /. s) /. m -> Range[0, mmax - 1]]],
           Sort[Flatten[Map[Transpose, Simplify[({x, y} /. s) /. m -> Range[0, mmax - 1]]], 1]]], {}]]

Testing with $d=60$ returns two fundamental solutions.

FermatReduceN[60, 4, 1]

{{8, 1}, {62, 8}}

The first case of $d=400004$ has one fundamental solution.

FermatReduceN[400004, 4, 1]

{{63917629397585002590254250700286910740527292089600124932270889438948\ 574385666998528790229868783788066961187502, 10106214045801334006424122096215491537979179535317961953309543166366\ 3323751876289858094985610380884355499500}}

The second case of $d=400012$ has two distinct fundamental solutions, the first and duplicate second of which are the smallest.

{{79796168015789079124029618418444253400324741593705623167764780247783\ 5503347693830111999576892990509205526208689827303068396833291338137932\ 536366968953450845322267394081990452, 12616692725197096567489442471571992079779468706087868590773991433311\ 9033353673236258181745141144365656065069607854001488956802709885815559\ 2443412490263634202021762111820545}, {797961680157890791240296184184442534003247415937056231677647802477835\ 5033476938301119995768929905092055262086898273030683968332913381379325\ 36366968953450845322267394081990452, 12616692725197096567489442471571992079779468706087868590773991433311\ 9033353673236258181745141144365656065069607854001488956802709885815559\ 2443412490263634202021762111820545}, {636742843000404002094201555617301178841445766262822122120505990753779\ 5717151103883018329373667345793746180265783230617455632746373301328580\ 3340913711943568284600677189686290013138415404806484947896724925142846\ 2606112523010581813076268617866384855712991589283925886210162447065345\ 346667709872705161453248309833515799078936771060814604566610219164302, 1006763732503411310533316722371570972693591699672489670428809573524\ 1743199532188167623411725516527477574571258078601173084338076716798459\ 4809316096277820789345143098516384464615848188948686356800775221809589\ 9235029680122446424203249605735651518329420622133464056744059240122938\ 436704113746518789382368654651027417415826893061292222194371027436340}\ }

$\endgroup$
  • $\begingroup$ Thanks. Does your code assumes that Reduce will return C[1] >= 0 or C[1] >= 1, but other constraints are sometimes generated, like C[1] >= -1. I'm not clear on the role of "m" in your code. $\endgroup$ – Kevin O'Bryant Jan 7 '18 at 3:40
  • 2
    $\begingroup$ I think one can solve the system (original quadratic and convergent ratio) for all convergents, and just keep the integer solutions.This example would be d = 400004; frax = Convergents[ContinuedFraction[Sqrt[d]]]; sol = Table[ Solve[{u^2 - d*v^2 == 4, u/v == frax[[j]], u >= 1}, {u, v}], {j, Length[frax] - 1}]; Cases[Flatten[{u, v} /. sol, 1], {_Integer, _Integer}]. $\endgroup$ – Daniel Lichtblau Jan 7 '18 at 16:33
  • $\begingroup$ I have not seen a case of C[1]>=-1 returned by Reduce, which does not mean it never happens. I recommend the code in the comment by @DanielLichtblau, which seems to be the most efficient solution, and it avoids cases of C[1]. The parameter m defaults to 1, and is the number of solutions to return for each class. $\endgroup$ – KennyColnago Jan 8 '18 at 16:51
1
$\begingroup$

Here's the code I'm using now, but it isn't nearly fast enough for what I need.

FirstPositiveSolution[d_] := 
  Module[{fi, u, v, u0, v0},
    {u, v} = ({u0, v0} /. 
       First[
       FindInstance[v0^2 - d u0^2 == 4 && v0 > 0 && u0 > 0, {u0, v0}, Integers]]);

    While[(fi = 
       FindInstance[v0^2 - d u0^2 == 4 && v0 > 0 && 0 < u0 < u, {u0, v0}, Integers]) != {},
       {u, v} = ({u0, v0} /. First[fi])];
    {u, v}
  ]

I find a solution with FindInstance, and then insist the next solution is smaller, repeat until FindInstance fails to find an instance.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.