Sunday, January 20, 2013

Hard Problems, Part 2

Some problems are hard nuts to crack. We have discussed this one before. Other problems are simply statically complex. As we will see, when things get hairy, there are still choices to make and techniques which work every time.

A statically complex problem has so many cases that it becomes daunting just to imagine a solution. This comes to mind because, of course, I am working on one of these problems right now. But I can't discuss it. Instead, we will look at another such problem and what it took to solve it.

In the 1980 I was working on a three-dimensional problem: rendering solid fractals. First Benoit Mandelbrot and then Loren Carpenter made famous the construction of mountains using fractal techniques. These were so-called plane-to-line functions, where the domain was a plane and the range was a line.

How to make nice fractals

In plainer English, you could create a mountain by making the height (z) a function of x and y. So, x and y are subdivided, and each time you subdivide, you add some random amount of z that depends upon the size of the x,y patch. In this way your mountain has a predictable fractal dimension. But now I would create such a thing by using Fourier transforms in the frequency domain: create a circular spot, transform it from the spatial domain to the frequency domain, randomize its phase, then transform it back into the spatial domain. This is shown in my blog post Textures, Part 2.

Once this technique is demonstrated, it becomes tempting to implement it. Which, of course, I did! I remember going to SIGgraph in 1980 and bringing a poster created at Calma on a Versatec raster plotter of some fractal mountains. Each facet was shaded using dither patterns, which I coded specifically for the poster. I presented it to Benoit Mandelbrot.

I'm sure Mandelbrot could've cared less about my rather amateurish fractal forgeries. There were so many people there, also. But I was undaunted and continued to play with fractals.

And so here comes the problem: what if you make a space-to-line function?

A space-to-line function is a function where w can be evaluated as a function of x, y, and z. Imagine that every point in space has a w value. And that w varies fractally. This construction can make clouds and rocks and all sorts of things!

So I decided to create a renderer that could display such a thing. The first way it could be displayed was simply to facet the zeroset of the surface, where w = 0. In a plane-to-line function, the zeroset z = 0 is simply a coastline. In a space-to-line function, the zeroset w = 0 is an actual surface of the object.

The problem can then be refined to this: how can you render the zeroset of the space-to-line function?

My answer to the problem was an early implementation of what is now called Marching Cubes. This technique was invented later by Lorensen and Cline and published in the 1987 SIGgraph proceedings. But I was using it in 1984. I won't say I invented it. Actually, it's a rather obvious solution.

My solution to the problem was to split the volume into small cubes. By evaluating the function for w at each corner of the cube, it will be found to be less than or equal to zero (a 0 bit) or greater than zero (a 1 bit). Since each cube has 8 corners, this was conveniently an 8-bit number.

Thus there were 256 cases to consider. That's a lot of cases, hence it was statically complex.

This is the first time I ever wrote a program to write the program to solve the problem. The program that generated the code for marching cubes was cloudtst.c, written in C on February 6, 1984, designed to run on a DEC VAX 11-780.

When a problem is this complex, it's hard to get every case right and consistent with every other case. This leads to the need to make sure each case is correct by having each case be generated by a program.

I categorized the eight bits into topological cases. And then I used the program to decide, for one of the 256 cases, which topological case it fell into. Then I mapped the solution for the topological case onto the 8-bit case.

Actually, the way I framed the problem was even more clever, since I allowed myself the luxury of subdividing a cube into 8 sub-cubes and faceting them. This meant that I could recursively subdivide the volume I was rendering, only subdividing it where it needed to be subdivided. This would minimize the total number of facets output and thus cut down on the time to render the fractal volume.

Code that writes code

It sounds a bit recursive, doesn't it? But there are several advantages to writing code that writes code. I will discuss them now.

The first and most obvious thing is this: if a program writes the code, then you won't have to. This is advantageous because there might be a lot of code to write. With a lot of complicated cases. And a program can make sure that each case is correct. So you won't have to. All you have to make sure of is that the program that analyzes the cases will get the right answer in each case.

The second advantage happens when you are generating a lower-level language. The program can see to all the details. This is useful when optimizing something. Lower-level languages are inherently more efficient because they are closer to the machine they run on. There is no hidden complexity. Each instruction is atomic. If you have to break each operation into smaller operations, the program can do it slavishly for you. This is the reason that people write compilers.

The third and most important advantage is also useful to you. To fix bugs, even those occurring in several cases, you only have to change the program that writes the code, not all the cases. This turns out to be characteristic of the code-writing-code process. Each bug in the generated code can be attacked one-by-one in the generating code. And eventually a perfectly correct program is the result. This can be done without tediously fixing each case in the generated code one at a time. Actually, many cases in the generated code can come from one part of the program, so a single bug fix in the generating program can fix several bugs in the generated program. This illustrates the power of this approach.

What language should you use? You have two choices to make. First, you choose a language to write the generating program in. This is usually a high-level language. One that you can work in with extreme freedom and ease. This can be something as utilitarian as C or ripe in data structures like C++. It can also be Python, a common choice nowadays. The second choice is the language used for the program you are generating. This should generally be chosen for its efficiency. Sometimes you do not have a choice. For instance, I have written code that generates programs in various forms of assembler, in C, in fragment shaders, and in OpenCL.



1 comment: