Re: [eclipse-clp-users] Large integers and delayed goals

Hi Joachim,

Thank you for detailed response.

There is a solution for A^4 + B^4 + C^4 = D^4, but it's probably too
large to find it in feasible amount of time:

a = 95800, b = 217519, c = 414560, d = 422481

As for your first suggestion to solve with ic in floating point
arithmetics and then test with exact arithmetic. Of course we will
weed out false positives with this approach, but isn't it possible
that a solution in exact arithmetic will not be solution in floating
point arithmetic?

Also for this particular problem stating C #< D instead of C #=< D
helps to weed out at least some false positives.

Sergey.

On Sat, Mar 17, 2012 at 1:29 AM, Joachim Schimpf <jschimpf_at_coninfer.com> wrote:
> Sergey Dymchenko wrote:
>> Hi,
>>
>> I want to find positive natural numbers A, B, C, D, such that A^4 +
>> B^4 + C^4 = D^4.
>> My program do this:
>>
>>     A #> 0, B #> 0, C #> 0, D #> 0,
>>     A #=< B, B #=< C, C #=< D,
>>     A*A*A*A + B*B*B*B + C*C*C*C #= D*D*D*D,
>>     labeling([A, B, C, D])
>>
>> But I get incorrect result [1, 1, 9742, 9742] with 10 delayed goals.
>> As far as I understand, the system finds that the result is imprecise,
>> but correct enough because 9742^4 is a large number.
>> Probably floating arithmetic is used, not big integers...
>> Is there a way to force precise arithmetic for large integers?
>
> Hi Sergey,
>
> you are right, library(ic) uses floating point arithmetic, which
> gives you 53 bits of precision.  Now, 9742^4 is just above 2^53,
> that's why you get the first pseudo-result for this magnitude.
>
> A simple, brute force way to avoid the pseudo-results is to add
> a test after search:
>
>     A #> 0, B #> 0, C #> 0, D #> 0, A #=< B, B #=< C, C #=< D,
>     A^4 + B^4 + C^4 #= D^4,
>     labeling([A, B, C, D]),
>     A^4 + B^4 + C^4 =:= D^4.   %%% HERE %%%
>
> A somewhat more elegant way is the following:
>
>     A #> 0, B #> 0, C #> 0, D #> 0, A #=< B, B #=< C, C #=< D,
>     [ic,suspend]:(A^4 + B^4 + C^4 #= D^4),  %%% HERE %%%
>     labeling([A, B, C, D]).
>
> By annotating the critical constraint with [ic,suspend], we set up
> two versions of the constraint: the 'ic' version (as before), plus
> the 'suspend' version (which delays until ground and then tests
> the equality using integer arithmetic).  If you have multiple
> constraints, you only need to do this for those that are at
> risk of imprecision due to the large numbers.
>
>
> While this solves the precision problem, I'm afraid I can't help
> with actually finding solutions - are there any?
> There is no magic in the propagation that finds [A,B,C] triplets
> whose 4th powers sum up to a 4th power, so it's mostly search.
> With your labeling order of A,B,C,D, you get stuck at
>
> ?- A #> 0, B #> 0, C #> 0, D #> 0, A #=< B, B #=< C, C #=< D,
>    [ic, suspend] : (A ^ 4 + B ^ 4 + C ^ 4 #= D ^ 4),
>    A = 1,
>    B = 1.
>
> A = 1
> B = 1
> C = C{85 .. 1.0Inf}
> D = D{85 .. 1.0Inf}
> There are 9 delayed goals.
> Yes (0.00s cpu)
>
> and then you keep trying values for C from 85 up to infinity.
> Unless there really is a solution with A=B=1, you will look for
> larger Cs forever, never going beyond the initial values for A and B.
>
> It may be better to label D first, beccause this immediately leads
> to bounds on the other variables:
>
> ?- A #> 0, B #> 0, C #> 0, D #> 0, A #=< B, B #=< C, C #=< D,
>    [ic, suspend] : (A ^ 4 + B ^ 4 + C ^ 4 #= D ^ 4),
>    D = 6.
>
> A = A{3 .. 5}
> B = B{3 .. 5}
> C = C{3 .. 5}
> D = 6
> There are 11 delayed goals.
> Yes (0.00s cpu)
>
> Now the remaining search space for [A,B,C] is bounded and can
> be exhausted in finite time.  So you can do:
>
> ?- A #> 0, B #> 0, C #> 0, D #> 0, A #=< B, B #=< C, C #=< D,
>    [ic, suspend] : (A ^ 4 + B ^ 4 + C ^ 4 #= D ^ 4),
>    labeling([D]),
>    writeln(d = D),
>    labeling([A, B, C]).
>
> d = 6
> d = 7
> d = 8
> d = 9
> d = 10
> ...
>
> and at least have the satisfaction of seeing something happening ;)
>
> [By the way, I have used X^4 instead of X*X*X*X in all the
> constraints, because that generally gives better propagation]
>
> Cheers,
> Joachim
>
> ------------------------------------------------------------------------------
> This SF email is sponsosred by:
> Try Windows Azure free for 90 days Click Here
> http://p.sf.net/sfu/sfd2d-msazure
> _______________________________________________
> ECLiPSe-CLP-Users mailing list
> ECLiPSe-CLP-Users_at_lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/eclipse-clp-users