Talk:Sequence of non-squares

I chose this as a way to show how easy it is to investigate functions in a programming language. --Paddy3118 08:42, 24 August 2008 (UTC)

Stability, Accuracy

The formula need to be investigated for numeric stability. Calculation of sqrt is inexact, it need to be shown that for all n in question, floor(1/2 + sqrt(n)) yields the exact result. A minimal requirement for this (though insufficient, I guess) is that sqrt has an error below 0.5. Note that the addition following to squaring will normalize for big n. Also conversion of those to floating point becomes quickly inexact when 32-bit floats are used. --Dmitry-kazakov 09:51, 24 August 2008 (UTC)

I suggest they choose/use the precision that gives correct results over the stated range, or reduce the range to suit the precision they have and add a comment. I think its true that most languages on Rosetta supporting floating point also support more than 32 bit FP values. --Paddy3118 11:28, 24 August 2008 (UTC)

Maybe it would be better to require an implementation not to use floating-point at all. I mean that m=floor(1/2 + sqrt(n)) should be solved iteratively using, say Newton method, in integers. That would guaranty that rounding error do not mangle the result. The case arises when sqrt(n) lands into the interval m.49999..m.50001. This is instable if the precision of floating-point numbers is 0.0001, because one cannot decide whether the floor of + 0.5 were m or m+1. Mathematically, for any given precision of floating-point numbers, there exist an infinite number of n, such that sqrt(n) is in the interval like above. There exist algorithms for integer sqrt(n) defined as the integer lower bound of real sqrt(n). One could modify them for 0.5+sqrt(n), that would be the cleanest way.

I have just noticed that there is a typo error in the task name: "sequance" instead of "sequence". --Dmitry-kazakov 13:30, 24 August 2008 (UTC)

Whoops! I have asked. --Paddy3118 13:56, 24 August 2008 (UTC)

The best way to avoid numerical instability is of course to remain integer all the time. My thoughts to this:

Let's assume that k is the largest integer so that k^2 <= n. Then it's obvious that sqrt(n) = k + x, where 0 <= x < 1.

Therefore:

• floor(1/2 + sqrt(n)) = k if x < 1/2, i.e. n < (k + 1/2)^2 = k^2 + k + 1/4. Since n and k are integers, this simplifies to n <= k^2 + k = k(k+1).
• floor(1/2 + sqrt(n)) = k+1 otherwise.

Maintaining k is quite simple: Since our n grows one by one, we just have to increment k as soon as n==(k+1)^2. Looking again at it, one sees that this can be further optimized: Since we need k+1 as soon as we are larger than k(k+1), we can just use that condition (which we would have to check anyway) for incrementing k. Therefore the algorithm would go like this (untested)

n := 1, k := 1
proposition_is_right := true
while (n <= 1000000)
  result := n + k
  if (result is square)
    proposition_is_right := false
  n := n + 1
  if (n > k*(k+1))
    k := k + 1

This doesn't use any floating point, and therefore isn't susceptible to rounding errors. Note that the largest number occuring in the calculation (namely k*(k+1)) isn't too much larger than n, so unless you get close to the upper end of your integer type range, you shouldn't get overflow. --Ce 15:33, 2 September 2008 (UTC)

You can rewrite it as n - k*k > k without risking integer overflow. Then k*k from the former step can be stored in order to reduce the number of multiplications down to O(sqrt(n)). So if this task is really about a sequence, i.e. n+floor(...) is not required to be programmed as a closed function, then this is looks like the most efficient solution. --Dmitry-kazakov 17:47, 2 September 2008 (UTC)

Well, using that (k+1)^2 = k^2 + 2k + 1, you can even completely eliminate all multiplications. However the stated goal of the task (see top of discussion page) was "to show how easy it is to investigate functions in a programming language." So while this algorithm is certainly a solution of the task (well, printing out some values is missing) avoiding the numeric problems, I'm not sure if it isn't against the idea of the task. Maybe a clarification of the task would be a good idea. --Ce 07:36, 3 September 2008 (UTC)

Try and implement something close to the given function that is being investigated, at least at first. --Paddy3118 18:03, 3 September 2008 (UTC)

Zero

I thought we were in the realm of Number Theory where Wikipedia states they don't include zero as a natural number, so went with the term Natural Number, but the term Positive Integer might be more exact. --Paddy3118 11:36, 24 August 2008 (UTC)

OCaml: superfluous truncate?

In this function for OCaml:

 let nonsqr n = n + truncate (floor (0.5 +. sqrt (float n)));;

Why do you need truncate as well as floor? Oh wait, does that change a floating point number with no fractional part into an integer? You may have guessed, I don't know OCaml. --Paddy3118 07:38, 25 August 2008 (UTC)

bc problem

Doesn't follow the task. Should print the given, limited range, then assert the condition for a million. Not print a million. --Paddy3118 13:05, 5 March 2009 (UTC)

To Bc makes no sense since there's no limit. No "assert" checked needed. You can put 10^6 or 10^100 as you prefer... (only time and memory (for arbitrary precision) are the limit. --ShinTakezou 13:09, 5 March 2009 (UTC)

Ah sorry, understood the sense now! Instead of reading the task, I've looked at sources, and get wrong understanding of the task, sorry! Going to fix it! --ShinTakezou 13:10, 5 March 2009 (UTC)

Back on my steps: the assertion exists to check the validity of the formula upto 10^6, or to check the validity of the implementations (because of precision) upto 10^6? I think the formula has no such limit; so the check is there because of limits that bc, by design, has not... so to bc makes not too much sense the "assertion" part. Anyway, it is there now, if I've got it right this time! --ShinTakezou 13:27, 5 March 2009 (UTC)