Coding: Fleeting Thoughts

[phlip]

↶
Maybe even more importantly, any operation on a float has to return another float, since the value may be outside the range representable by ints.

What's intfloor(2147483649.5)?

Well, we're talking about a language with bignums, so that's not a concern here.

The lowlevel C floor function, that makes sense for. And the lowlowlevel assembly floor function it makes even more sense, because it means that opcode can be implemented entirely by the floating point components and it doesn't need to connect to the integer registers.

But a higher-level function doesn't need to be beholden to the lowlevel API...

[ahammel]

↶
Maybe even more importantly, any operation on a float has to return another float, since the value may be outside the range representable by ints.

What's intfloor(2147483649.5)?

Well, we're talking about a language with bignums, so that's not a concern here.

Even python doesn't have infinite and NaN bigints, though (right?).

I guess intfloor(NaN) could just throw an exception.

[Yakk]

For -2¹⁰²⁴ < x < 2¹⁰²⁴ (or thereabouts, I CBA to calculate the exact endpoints), every integer is exactly representable as a float, so you shouldn't get any rounding errors at all. Outside that range, every possible float value is an exact integer, so any rounding errors will take you to ints, not away from them (and in practise it just means that floor(x) == x, since x is always going to be integral).

Pigeonhole says that range is too large?

[phlip]

Yeah, had a brainfart there. 2¹⁰²⁴ is just the maximum number you can fit in a double at all, not range of integers in doubles.

That cutover is at 2⁵² to 2⁵³. Every integer with absolute value less than 2⁵³ can be represented exactly as a double, and every double bigger than 2⁵² is exactly an integer.

[ucim]

Ok, I didn't realize that floats could hold exact integers up to that big, and was confusing it with the results of float operations that should be integers not always coming out that way. For example, while 6/3 is always 2, is the float (or double) result of 6.0/3.0 guaranteed to be equal to the float representation of exactly 2?

If 6.0/3.0 is 1.99...9913, then floor(6.0/3.0) would be either 1 or 1.0; wrong either way. However, floor(6.0)/floor(3.0) would always be two if floor returned an int.

Jose

[Demki]

Seeing as floating point division is done by dividing the mantissa and subtracting the exponents, in the case of 6.0/3.0, you should get exactly 2.0, as the mantissa are 1.1 and 1.1, and the exponents would be 2 and 1( so we have (1.1 × 2^(2))/(1.1 × 2^(1)) so we get (1.1/1.1) × 2^(1), and since the mantissa are represented as integers, we have 1 × 2^(1) which is 2

9/3 would be (1.001 × 2^(3))/(1.1 × 2^(1)), so we have (1.001)/(1.1) × 2^(2), and applying integer division to (1.001/1.1) (which is just 1001/11 × 2^(-2)) gives 11 × 2^(-2) × 2^(2) = 11 × 2^(0), and since floating point values need to be normalized(ie of the form 1.mantissa × 2^(exponent)) we have 1.1 × 2^(1)

So yes, you get the expected result for integers. (Note the mantissa is basically a fixed point number, where the point is before the leftmost position and 1 is prepended, since all numbers start with 1. In normalized form we save a bit by not representing it in the mantissa)

[ucim]

So yes, you get the expected result for integers.

So then, why is

Code: Select all

float length;
float width;
length = 10.0;
for (width=1.0 ; width=width+.01 ; length==width)
   print ("still too long");
print ("Finally, perfect!");


wrong?

Jose

[korona]

So yes, you get the expected result for integers.

So then, why is

Code: Select all

float length;
float width;
length = 10.0;
for (width=1.0 ; width=width+.01 ; length==width)
   print ("still too long");
print ("Finally, perfect!");


wrong?

Jose

Because floats cannot represent 1/100 exactly. Just like the decimal base cannot represent 1/3 exactly (if you only allow a finite number of digits). If you replace 0.01 by 0.25, 0.625 or similar (i.e. any number that is a sum of inverses of powers of 2) it works perfectly fine.

[Demki]

So yes, you get the expected result for integers.

So then, why is

Code: Select all

float length;
float width;
length = 10.0;
for (width=1.0 ; width=width+.01 ; length==width)
   print ("still too long");
print ("Finally, perfect!");


wrong?

Jose

As Korona said, floating point values cannot accurately represent any fraction that does not have a power of 2 as its denominator.
Just curious, what language is that? First time I see a for loop in the order (initialization; increment; condition) instead of the order (initialization; condition; increment)

[ucim]

Because floats cannot represent 1/100 exactly.

Ah, ok. Makes sense. So if I did
floor(length)==floor(width)
then it wouldn't be wrong? (at least it wouldn't slither past the condition)

Just curious, what language is that? First time I see a for loop in the order (initialization; increment; condition) instead of the order (initialization; condition; increment)

It's in C-mustard. It's an unintentional language that isn't taught much any more (though it is used more often than one might suspect).

Jose

[Demki]

Because floats cannot represent 1/100 exactly.

Ah, ok. Makes sense. So if I did
floor(length)==floor(width)
then it wouldn't be wrong? (at least it wouldn't slither past the condition)

If you did that you would indeed get it to terminate, it probably will not be after the expected number of iterations. You could actually just use(in this case) floor(width)==length
It is still better to order comparisons on floating points rather than equality.

Just curious, what language is that? First time I see a for loop in the order (initialization; increment; condition) instead of the order (initialization; condition; increment)

It's in C-mustard. It's an unintentional language that isn't taught much any more (though it is used more often than one might suspect).

Jose

Well, looking for it I found there are way too many programming languages.

[ucim]

Well, looking for it I found there are way too many programming languages.

Yes, there are.

"Mustard" is OTTish for "screwup". C-mustard is the language I inadvertently created a few posts ago. It's identical to C except that in some cases the for statement has its conditions in a different order. Alas, I never created proper documentation for it.

I guess I could have also named it after C-sharp, calling it C-flat-on-your-face.

Jose

[Xanthir]

The important thing to remember, as said above, is that floats are not random. They don't produce weird results "just because". They're a base-2 value, with a movable "window" of precision. As a base-2 value, they can only represent values divisible by a power of 2. This is largely incompatible with base-10 fractions, tho, and so most base-10 values turn into an inexact base-2 value (just like anything not divisible by powers of 2 or 5 can't be written as an exact decimal in base-10), and then you get further slight divergence when you convert them back into base-10 for printing.

[Qaanol]

The important thing to remember, as said above, is that floats are not random. They don't produce weird results "just because". They're a base-2 value, with a movable "window" of precision. As a base-2 value, they can only represent values divisible by a power of 2. This is largely incompatible with base-10 fractions, tho, and so most base-10 values turn into an inexact base-2 value (just like anything not divisible by powers of 2 or 5 can't be written as an exact decimal in base-10), and then you get further slight divergence when you convert them back into base-10 for printing.

In addition to the accuracy problem you describe, there is also an issue of precision loss when adding (or subtracting) floating-point values of disparate magnitude:

a = x + y
b = a - y // algebraically equals x, but won't always in practice
c = b - x // algebraically equals 0, but won't always in practice

This may come up when summing a large number of values. As the total gets bigger, the low-order bits of the values “fall off the end” and get discarded instead of being added to the sum.

You can correct for it via compensated summation, by keeping a tally of how much has been lost and adding it to the next value before summing.

[raudorn]

Assume initialization and bounds checking in the following.

Code: Select all

class Grid {
private:
    double* _data;
   
    double Cell(int pos) {
        return _data[pos];
    }
   
public:
    void crunch() {
        for(int i = 0; i < largeNumber; i++) {
            // Some arithmetic expression like
            float foo = (this->Cell(i) + this->Cell(i+1)) * this->Cell(i-1) + this->Cell(i);
        }
    }
}

Runtime: x seconds

Code: Select all

class Grid {
private:
    double* _data;
   
public:
    void crunch() {
        for(int i = 0; i < largeNumber; i++) {
            // Some arithmetic expression like
            float foo = (_data[i] +_data[i+1]) * _data[i-1] + _data[i];
        }
    }
}

Runtime: 3 * x seconds

[incoherent, angry screaming]
[breaking computer noises]
[sobbing]

[korona]

And that, kids, is why we pass -O2 on the command line.
(Runtime: 0 seconds because the dead store and the loop is eliminated )

Seriously though, it does look quite strange. What compiler / options were you using? Can you paste the real code (the real calculation and initialization and bounds checking?

It could be a cache thrashing artifact.

[raudorn]

I'm afraid the real code is a bit too lengthy to be properly discussed here. Basically we have a data grid that calculates various derivatives like

Code: Select all

real_t Grid::DC_vdv_y(const Iterator &it, const real_t &alpha) const {
  real_t ft = pow((this->Cell(it) + this->Cell(it.Top())), 2.0)
    - pow((this->Cell(it.Down()) + this->Cell(it)), 2.0);
  real_t st = fabs(this->Cell(it) + this->Cell(it.Top())) * (this->Cell(it) - this->Cell(it.Top()))
    - fabs(this->Cell(it.Down()) + this->Cell(it)) * (this->Cell(it.Down()) - this->Cell(it));
  return (0.25 * (ft + alpha * st)) / _geom->Mesh()[1];
}

and I thought "Hm, all those this->Cell(it) calls are superfluous, since they just access the data that is stored in the 'this' instance anyway." So I replaced them with direct calls to _data[it] and suddenly performance is three times slower. Default flags are "-Wall", "-Wextra", "-pedantic", "-std=c++11" and "-O3". Honestly, I'm not sure what they all do, since this is a class project and we use a given build system.

[Flumble]

Default flags are "-Wall", "-Wextra", "-pedantic", "-std=c++11" and "-O3". Honestly, I'm not sure what they all do, since this is a class project and we use a given build system.

The first three are for printing out every possible warning (any sign that you might not strictly follow the C++ standards), the std flag is for compiling like it's C++11 (rather than an older or newer version of the language) and -O3 is for optimizing the executable for speed (rather than memory usage, executable size or whatever).

Could it somehow be because of the int argument of Cell? I'd try running it with size_t as argument type and see if it runs thrice as slow.

[Yakk]

Another possivility is that the conversion to value means you don't have a `double&` but instead a `double`; ie, more locality. Simpler locality (working with temporary copies, not references to data in the array directly) could make a difference. You can emulate this with a `static_cast<double>( m_data[it] )`.

So cast in-to-`int` and out-to-`double` separately and together and see if any get faster.

[raudorn]

So, today I tried to replicate the problem... and couldn't. The same changes to the code now had no measurable impact on performance at all. I would've liked to look further into it, but we had more pressing concerns.

Another possivility is that the conversion to value means you don't have a `double&` but instead a `double`; ie, more locality.

I assumed something along those lines, but since the artifact went away, I doubt it again.

I guess the moral here is that how you measure performance is just as important and non-trivial than improving it.

[liberonscien]

Something occurred to me a while ago:
Is it possible to code a program that randomly generates new code, which the program then tests?
Like assign a number to every code piece, generate a string of pseudorandom numbers, translate that string into a generated code which the program then runs.
There would have to be failsafes in place that check and see if the generated code is generating errors. If it detects errors in the debug session it would start over and generate a new code. Also, if it somehow generates a code that loops, it would have to have a system in place that restarts the generation after a specified period of time.
For this to have any value at all, it would need to be able to save its generated codes.

I find myself wondering if something like this is possible, theoretically speaking, or is it impossible to create something like this?

[raudorn]

It's certainly possible, though I can't point to an existing implementation at the moment. That said, I can think of a few problems from the musing I did on the same topic:

1.) The halting problem means that we can't have a function that tells us whether the generated program will halt. What we can do is have a bunch of heuristics that will check against unreasonably long execution time or maybe check the loop conditions, but this will result in false negatives, e.g. programs that are useful and will halt, but not in the arbitrary time limit we set.
2.) The space of syntactically correct programs is vastly larger than the space of syntactically correct and useful programs, so we can't really expect to get a useful result anytime soon, where soon can take values like the lifespan of the universe.
3.) Checking the result is non-trivial in two ways:
a) If we check only the return value, we ignore useful programs that do not return anything, but have useful side effects. I guess in some languages all side effects are considered evil, so maybe in those it's not a problem.
b) We either can only check for classes of test parameters (e.g. test for -1, 0, 1, 2, 99999 and NaN) or we restrict them to a necessarily narrow space and test all combinations. The first is not guarranteed to properly check against invalid behaviour, since it can be hidden in parts of the parameter space that are not covered by the test. The second one is impossible in practice, unless you have a problem with an extremely narrow parameter space.

[Tub]

The good news: programs are just finite-length strings of bits, so you can enumerate them if you want.

The bad news: proper testing is impossible (not just due to the halting problem), and heuristic testing will pass too many broken programs.

Say, you're testing a sorting algorithm on [], [1,2,3,4,5] and [3,2,5,4,1]. That might otherwise qualify as reasonable test cases. But which of the following two programs would you expect to be randomly generated first?

Code: Select all

for i = 1 to length(input)
  for j = 2 to length(input)
    if input[i] < input[i-1]
      t = input[i]
      input[i] = input[i-1]
      input[i-1] = t
return input

Code: Select all

if (input == [])
  return []
return [1,2,3,4,5]

Any course on theoretical computer science will teach you why brute-force programming is a bad idea. Even if you bypass the halting problem with a timeout, and even if you bypass proper testing with a handful of testcases, your runtime is still exponential against code length, and any code of value is too long for a computer to randomly generate within any reasonable timespan.
Seriously, there are hard upper bounds on the computing capacity of the whole observable universe, and if you do the math, then you're running out of universe before you've tested programs with maybe 50 characters. In some languages, that's not even enough for "hello world".

[phlip]

It's also worth name-dropping Rice's theorem, a generalisation of the halting problem that any question like "is this a sorting algorithm" is undecidable.

[Tub]

RRDTool seems like a pretty useful thing. It can store regular measurements (like CPU temperature, system load etc) in a way that avoids the overhead of RDBMS et al, and it'll automatically provide lower-resolution series with different cutoffs. For example, store the last 24 hours in 1-minute accuracy, the last week with 1 hour accuracy, the last year with daily accuracy.

Neat.

There's a whole ecosystem of tools that'll collect your data and insert it into a rrdtool database, making it pretty easy to setup comprehensive monitoring for your servers, your network devices, your home automation etc. collectd seems lightweight, feature-rich and modular.

Great.

It'll also produce little graphs from the stored data, With just a few lines of code, you can program a dashboard for your servers' status (or pick one of several existing implementations).

Awesome.

To produce the graphs, it uses cairo. When I want to install rrdtool, the dependencies include cairo, mesa, libdrm, several X libraries and wayland. On a headless server.

What the fuck. Who considered that to be a good idea?

[korona]

At least on Debian RRDTool just pulls in cairo. This is fine: cairo is a generic vector graphics drawing library and fine for headless systems; it's the right tool for the job. However cairo is built with support for output to X11 surfaces and thus pulls in X11 libraries (in wayland based systems X11 is replaced by a combination of wayland and OpenGL, i.e. mesa + libdrm). That makes sense if you're running a Desktop (GTK+ draws via cairo; I'm not sure about Qt) but not so much if you're on a headless server. I do not know if cairo's upstream supports building the X11 part as a separate module; in any case this seems to be cairo's or the maintainer's fault and not a problem of RRDTool.

[Wildcard]

Something occurred to me a while ago:
Is it possible to code a program that randomly generates new code, which the program then tests?
Like assign a number to every code piece, generate a string of pseudorandom numbers, translate that string into a generated code which the program then runs.

Consider the equivalent for English:

Generate strings of random words. Long, long strings. A whole book. Do this millions of times.

Now sift through the output and try to find the books which make sense and express some idea clearly or at least entertainingly.

The vast majority will be garbage.

Even if we throw out all books for which any page is garbage, the vast, vast majority of what remains will still be garbage.

Go to the library; pull 500 independently and randomly selected books off the shelves; tear a random page out of each and stack the result and call it a book. That type of result is exponentially more likely than finding the actual contents of one of the books at the library from our random generation.

Jonathan Swift had a section along this line in "Gulliver's Travels." Not a new idea by any means.

Further philosophical thoughts on this subject spoilered for somewhat off-topicness:

Spoiler:

The interpretation of the result of this thought experiment is usually missed by mathematicians, who like to number everything; and that is: Ideas are not enumerable. Ideas are created.

It's entirely possible to create a new idea that's never been thought of before. In order to communicate the idea, you will have to either (a) compare it to other ideas that have been thought of before, so that you can talk about it using words that other people know, or (b) just directly cause others to experience the new idea themselves (perhaps a new kind of perception like smelling in color), at which time you can invent new words to talk about this idea which you both know.

(By the way, if you're not convinced, try to describe color to a man who was born blind. Then tell me again that all ideas are enumerable. The symbols might be enumerable, but the ideas certainly aren't, and likewise "s/.\(.\)$/\1/" didn't mean anything until we assigned it meaning—despite the potential existence of the symbols.)

[Link]

(I'm not even sure where to ask this, so I figured the Coding Fleeting Thoughts thread would be a good place to start.) I'm working on a numeric research project using Python (with Numpy, SciPy, etc.), and I have the misfortune of being faced with a tensor equation of the form T^a_bc x^b x^c+M^a_b x^b=y^a (in Einstein notation) that I need to solve for the vector x. Does anyone know if there are any pre-existing packages to tackle this sort of monstrosity?

[Thesh]

Do you use Linux, by any chance? You can load SageMath libraries from Python, or just write python inside of it (easier since it already loads most of the packages you want), which can solve equations, plot graphs, do a ton of mathstuffs.

http://www.sagemath.org/

[Link]

Yep, and I know about Sage. AFAIK it doesn't really do anything you can't do with Numpy and SciPy, though. There's also the issue of speed: I need this to be as efficient as possible, and anything I've done in Sage has thus far been slower than the equivalent in pure Numpy.

[Thesh]

Have you tried GraviPy ("Tensor Calculus Package for General Relativity")? That's what google came up with. Don't know the performance, although it is based on SciPy. If speed is really that critical, I would look into a library for a compiled language, like C++. Python is just not that fast.

[Link]

↶
Have you tried GraviPy ("Tensor Calculus Package for General Relativity")? That's what google came up with. Don't know the performance, although it is based on SciPy. If speed is really that critical, I would look into a library for a compiled language, like C++. Python is just not that fast.

Nope, that's symbolic; I need a numeric package. I basically want numpy.linalg.solve, but for quadratic tensor equations instead of linear matrix ones.

I agree that C or C++ would be a bit better if I could find a highly optimised library, but Python is the go-to language in my group, and the idea is that I write the code in a way that others could use it with very little introduction, so switching languages isn't an option.

[ahammel]

I agree that C or C++ would be a bit better if I could find a highly optimised library, but Python is the go-to language in my group, and the idea is that I write the code in a way that others could use it with very little introduction, so switching languages isn't an option.

If you do find a C(++) library that does what you want, you could always write Python bindings for it.

[Link]

That's true. In fact, I recently discovered cffi, which makes accessing C functions incredibly easy. Still, it would be nice to have something that can interface with Numpy without too much glue code.

[Xeio]

Sometimes I wonder if I'm actually terrible but just manage to luck my way through everything.

Like, say I'm not quite sure an easy way to test XML submission to a particular service endpoint? Why not just debug the IIS process, put a breakpoint on the BackgroundWorker's AssignWork and just create the work object wholesale during runtime! Oh it has to read an encrypted but I only have the unencrypted file? I'll just step over all that decryption and do a "x = File.ReadAllText()".

I mean, as long as I don't care about tesing the submission process... right?

TL;DR: I have this nagging feeling VS's awesome debug tools are teaching me bad habits.

[>-)]

i've been playing around with brute force testing some open but obscure number theory conjectures on a high end GPU recently. the speedup is nice about 10x ... but not as much as expected. i think this is due to the relatively low int64 compute performance on a gpu compared to the float32 performance, it's maybe only 10%. however, i recently heard of this neat trick for speeding up integer computations by casting them to floats and proving hoping that nothing goes wrong when you do the same operations, so i'm about to try that out now.

i also have to wait until the int64's are small enough before casting them to floats or else i lose too much precision. the algorithm does make the numbers smaller at each iteration, so that's good news. the bad news is that because of the extremely limited shared memory on the GPU, this will probably involve some tricky manipulation in order to cast int64's to float32's *in place*

edit: managed to get 3x speedup by doing everything with floats

[Xenomortis]

*Has a pull request on github with a bunch of commits*
*Goes on holiday for nearly three weeks*
*Comes back to see someone added some more commits*

git reset --hard abc123
git push -f

I was productive today.

[Yakk]

Considered using floats as 17 bit ints, and composing 4 of them to make a 64 bit int?

[>-)]

That sounds pretty horrifying, and i'm not sure it would run faster after all the overhead. Luckily, I was dealing with prime gaps, so even though the primes themselves had to be put in int64's, the gaps were nice and small enough to fit in a float

edit:
hmm, so according to nvidia's cuda documentation, int32 add/subtraction is as fast as float32 operations. only int32 multiplcation is expensive, and all int64 instructions are very expensive. maybe just casting to int32 would have been sufficient after all.

[Xanthir]

Considered using floats as 17 bit ints, and composing 4 of them to make a 64 bit int? :)

Why 17? You've got 24 bits of integer precision in a float (including the hidden bit), so you only need three of them.