curv icon indicating copy to clipboard operation
curv copied to clipboard

How to construct polyhedra

Open gsohler opened this issue 7 years ago • 16 comments
trafficstars

Hi Group,

Right now I am quite excited about curv, this is why i wanted to port all the solids from Mercury.

I have:

===============

let phi =(sqrt(5)+1)/2; testvec=identity(3); GDFVectors = [ normalize([1, 0, 0]), normalize([0, 1, 0]), normalize([0, 0, 1]),

    normalize([1, 1, 1 ]),
    normalize([-1, 1, 1]),
    normalize([1, -1, 1]),
    normalize([1, 1, -1]),

    normalize([0, 1, phi+1]),
    normalize([0, -1, phi+1]),
    normalize([phi+1, 0, 1]),
    normalize([-phi-1, 0, 1]),
    normalize([1, phi+1, 0]),
    normalize([-1, phi+1, 0]),

    normalize([0, phi, 1]),
    normalize([0, -phi, 1]),
    normalize([1, 0, phi]),
    normalize([-1, 0, phi]),
    normalize([phi, 1, 0]),
    normalize([-phi, 1, 0])

];

c0=GDFVectors[0]; c1=GDFVectors[1]; c2=GDFVectors[2]; c3=GDFVectors[3]; c4=GDFVectors[4]; c5=GDFVectors[5]; c6=GDFVectors[6]; c7=GDFVectors[7];

my_object = make_shape { dist(x,y,z,t) = let v=max( abs(dot([x,y,z],c0)), abs(dot([x,y,z],c1)), abs(dot([x,y,z],c2)), abs(dot([x,y,z],c3)), abs(dot([x,y,z],c4)), abs(dot([x,y,z],c5)), abs(dot([x,y,z],c6)), );

// tmp1=max( // for (c in [c3,c4,c5,c6]) // abs(dot([x,y,z],c)), // );

    in v- 1;
    is_3d = true;
    bbox = [[-1,-1,-1],[1,1,1]];

};

    tmp2=max(
    for (c in [c3,c4,c5,c6])
            abs(dot([1,1,1],c)),
     );

in

my_object

================

But I am facing small problems and i am not sure if curv yet support this. I see that such functionality is already coded, but i'd like to further automate this.

  • Is it possible to move the for loop inside the make_shape and further inside the max function ? I get an ERROR: at field .dist: this expression is not supported by the geometry compiler: curv::For_Op

  • Is it possible to directly access GDFVectors[] without using temporary variables c0,c1,c2,c3 ? I get an ERROR: syntax error in delimited phrase

I tried several syntax constructs, but none of them worked, I would appreciate any help here.

best regards Guenther

gsohler avatar Aug 16 '18 09:08 gsohler

The geometry compiler is limited to a subset of Curv, and many things that look like they should work, do not. In part, this is unavoidable, because in a fragment shader, you cannot allocate memory, there are no recursive functions, there are no pointers. But the geometry compiler has more limitations than are necessary. I'd like to redesign it at some point, and remove as many limitations as possible.

The good news is that you can successfully transliterate most GLSL programs into Curv colour and dist functions. You need to write code that is very close to the original GLSL, for the most part. I do support the 'for' loop, but only in a restricted form. for (i in 0..10 by 2) is compiled into for (float i=0.0; i <10.0; i += 2.0).

Two big limitations for translating GLSL code are: no arrays, and no matrices. So there is no direct translation for GDFVectors in the Mercury library. I did look at the Mercury code when I was implementing the platonic solids: see cube, tetrahedron, octahedron, dodecahedron, icosahedron in std.curv, and you may see some resemblances. I think I may have copied some data from the GDFVectors table, even if the code is different.

Fixing these limitations (eg, adding support for GLSL arrays and matrices) is important to me: I'd like that done for the 1.0 release.

doug-moen avatar Aug 16 '18 15:08 doug-moen

Hi Doug,

Thank you for the exact for loop translation pattern.

Actually in Mercury/GLSL they model solids with a for loop. their code is:

// Version with without exponent, creates objects with sharp edges and flat faces float fGDF(vec3 p, float r, int begin, int end) { float d = 0; for (int i = begin; i <= end; ++i) d = max(d, abs(dot(p, GDFVectors[i]))); return d - r; }

when I try this in my code:

cube_solid = make_shape { dist(x,y,z,t) = let d=0; for(i in 0..2) d=max(d, abs(dot([x,y,z],GDFVectors[i])));

    in d - 1;
    is_3d = true;
    bbox = [[-1,-1,-1],[1,1,1]];

};

===

All the vector should eliminate to constants, but it tells: Error: not an operation

Where is the error ? can i look into the generated GLSL code ?

gsohler avatar Aug 17 '18 06:08 gsohler

Hi Guenther. You are now trying to use parts of Curv that are not documented. So I'll try to fill in the missing details.

for (int i = begin; i <= end; ++i)
    d = max(d, abs(dot(p, GDFVectors[i])));

This is a syntax error. The body of your for loop is d = ..., which is a definition, not an assignment statement. The body of a for loop must be some kind of statement, or a parenthesized statement sequence, and a definition is not a statement.

Curv has assignment statements, but they are not documented. The syntax is d := ....

An assignment statement assigns a new value to a mutable variable. Mutable variables are not documented. To use a mutable variable, you must first define it, using a definition of the form

   var d := <initial value>;

Mutable variables and assignment statements can only be used inside of a do clause. Curv is a pure functional language, originally inspired by OpenSCAD and Haskell. It is not an imperative language like Python, Javascript, etc. But Curv contains a simple embedded imperative language whose variable definitions and statements can only be used inside a do clause. The syntax is

    do
        <mutable variable definitions and statements>
    in
        <expression that references mutable variables defined above>

I use this embedded imperative language when I am translating code from GLSL.

Here are some examples:

  • Uses a while loop: https://github.com/doug-moen/curv/blob/master/examples/mandelbrot.curv
  • Notice how I chain let clauses and do clauses: https://github.com/doug-moen/curv/blob/master/examples/smoke.curv

doug-moen avatar Aug 17 '18 21:08 doug-moen

Also, you are using an array, GDFVectors, inside of a dist function. That won't work, because arrays are not currently supported by the geometry compiler. It's a limitation I want to fix for release 1.0.

doug-moen avatar Aug 17 '18 21:08 doug-moen

Hi, beeing able to use the new "do - in" syntax i was able to implement the for loops as intended - see code below - There are just two issues left:

  • Using function basing on function "gdf_index" the program is very perormant anymore
  • looking at model "dod_wire" there are artefacts, where only "half" of model is drawn, the rest is not. I could not yet spot the problem

P.S. is there a better way in curv to calculate max(values) and 2nd least max(values) ???

========== CODE starts here ==================

let phi =(sqrt(5)+1)/2;

gdf_index(i) =

if(i == 0) normalize([1, 0, 0])
else if(i == 1) normalize([0, 1, 0])
else if(i == 2) normalize([0, 0, 1])

else if(i == 3) normalize([1, 1, 1 ])
else if(i == 4) normalize([-1, 1, 1])
else if(i == 5) normalize([1, 1, -1])
else if(i == 6) normalize([1, -1, 1])

else if(i == 7) normalize([0, 1, phi+1])
else if(i == 8) normalize([0, -1, phi+1])
else if(i == 9) normalize([phi+1, 0, 1])
else if(i == 10) normalize([-phi-1, 0, 1])
else if(i == 11) normalize([1, phi+1, 0])
else if(i == 12) normalize([-1, phi+1, 0])

else if(i == 13) normalize([0, phi, 1])
else if(i == 14) normalize([0, -phi, 1])
else if(i == 15) normalize([1, 0, phi])
else if(i == 16) normalize([-1, 0, phi])
else if(i == 17) normalize([phi, 1, 0])
else  normalize([-phi, 1, 0])

;

solid(i1,i2) = make_shape { dist(x,y,z,t): do var d := 0; var i := i1; while(i <= i2) ( d := max(d,abs(dot([x,y,z],gdf_index(i)))); i := i+1; )

in d - 1;
    is_3d : true;
bbox : [[-1,-1,-1],[1,1,1]];

};

wire(i1,i2) = make_shape { dist(x,y,z,t) = do var d := 0; var i := 0;

var imax := 0;

var max1 := -inf;
i := i1;
while(i <= i2)
(
	d := abs(dot([x,y,z],gdf_index(i)))-1;
	if(d>max1) (
		max1 := d;
		imax := i;
	);

	i := i+1;
);

var max2 := -inf;
i := i1;
while(i <= i2)
(
	d := abs(dot([x,y,z],gdf_index(i)))-1;
	if(d>max2 && i != imax) (
		max2 := d;
	);

	i := i+1;
);

in mag(max1,max2)-0.05;
    is_3d = true;
bbox = [[-2,-2,-2],[2,2,2]];

};

cube_solid = solid(0,2); oct_solid = solid(3,6); ico_solid = solid(3,12); dod_solid = solid(13,18);

cube_wire = wire(0,2); oct_wire = wire(3,6); ico_wire = wire(3,12); dod_wire = wire(13,18);

== END of Code

gsohler avatar Aug 20 '18 12:08 gsohler

if statements are very slow on a GPU. That long chain of if statements in gdf_index looks like a performance killer.

for loops have a better chance of being optimized by the GPU driver than while loops. I suggest converting

var i := 0;
i := i1;
while(i <= i2) (
	...
	i := i+1;
);

to

for (i in i1..i2) (
   ...
);

The Nvidia driver has a very aggressive compiler. It might decide to unroll the for loop, which might then allow it to eliminate the if chain in each iteration of the loop.

And if that doesn't work, then I would need to do some work on the Curv GL compiler so that it can produce code that is easier for Nvidia to optimize. Adding array support is the obvious thing, but there might be easier things to be tried first.

doug-moen avatar Aug 20 '18 15:08 doug-moen

Actually, I'm getting 60 FPS on the nvidia gtx 1050 without any code changes.

doug-moen avatar Aug 20 '18 15:08 doug-moen

I put in the for loop now and its reasonably faster. However it does not display an FPS count in my place anymore. Do you also see the artefacts (wires missing) in your place ? is it the the X-Ray spheres dont fit thorugh the wireframes in the back to exposure the front wireframe contours ? Or do the X-Rays start from the front ? Is it possible to configure the initial sphere size ?

gsohler avatar Aug 20 '18 17:08 gsohler

I'm out of town, so I just have my 2011 laptop now. I get 20-30 FPS with dod_solid and 20 FSP with dod_wire. Changing while loop -> for loop doesn't change rendering time. The FPS is being displayed. My laptop GPU is really old and relatively slow, but it has the closed source Nvidia driver, so maybe that's why it is apparently beating your system.

The artefacts are caused by a bad distance function. Not Lipschitz-1. If you try

   dod_wire >> lipschitz 2

then the problem goes away. This affects rendering speed, so you might want to experiment and find the lowest Lipschitz number > 1 that fixes the problem.

Nice work, by the way.

doug-moen avatar Aug 20 '18 20:08 doug-moen

General suggestions for writing efficient distance functions:

  • Avoid if statements and if expressions.
  • Avoid loops, but simple for loops over a numeric range are faster than while loops.
  • Use 'vectorized' operations where possible. The fundamental vector types are vec2, vec3 and vec4. For example, [x,y,z] is a vec3. [x,y,z]+[a,b,c] is vectorized addition, and it runs as fast as adding two numbers, due to vector hardware in the GPU.

For this specific task, you are creating wireframes of polyhedra. I would suggest looking at Knighty's "fold and cut polyhedra" technique. It's more general than the Mercury polyhedra code, and should be fast.

  • Syntopia's blog post, with some pictures and a link: http://blog.hvidtfeldts.net/index.php/2012/01/knots-and-polyhedra/
  • The Fragmentarium project on github contains multiple examples of this technique: https://github.com/3Dickulus/FragM/blob/master/Fragmentarium-Source/Examples/Knighty%20Collection/Fold%20and%20Cut%20Polyhedra%20-%20Original.frag

I recommend Syntopia's blog and the Fragmentarium project as a rich source of techniques and GLSL code that can be ported to Curv.

doug-moen avatar Aug 21 '18 12:08 doug-moen

Here's a dodecahedron wireframe, using Knighty's code:

let
Type = 5; // symmetry group type, range 3..5
U = 0; // barycentric coordinate for the principle node, range 0..1
V = 0; // range 0..1
W = 1; // range 0..1
SRadius = 0.1; // segments radius

nc = [-0.5,-cospin,scospin];
cospin = cos(pi/Type);
scospin = sqrt(0.75-cospin*cospin);

p = normalize(U*pab+V*pbc+W*pca);
pab = [0,0,1];
pbc = normalize[scospin,0,0.5];
pca = normalize[0,scospin,cospin];

fold pt =
    do  var pos := pt;
        for (i in 0..Type) (
            pos := [abs(pos[X]), abs(pos[Y]), pos[Z]];
            let t = -2*min(0,dot(pos,nc));
            in pos := pos + t*nc;
        );
    in pos;

D2Segments pt =
    let pos = pt - p;
        dla = mag(pos-min(0,pos[X])*(1,0,0));
        dlb = mag(pos-min(0,pos[Y])*(0,1,0));
        dlc = mag(pos-min(0,dot(pos,nc))*nc);
    in min(min(dla,dlb),dlc)-SRadius;//max(max(dla,dlb),max(dlc,dlp))-SRadius;

in
make_shape {
    dist pt = D2Segments(fold(pt[[X,Y,Z]]));
    is_3d = true;
    bbox = [[-1,-1,-1]-SRadius, [1,1,1]+SRadius];
} 

doug-moen avatar Aug 21 '18 13:08 doug-moen

works perfect for me! I am wondering if curv can also implement some graphical sliders like FragM does ? Appearently curv can already show "video", this is picture dependant on time. Can the picture also be dependant on the value of a slider ?

gsohler avatar Sep 19 '18 08:09 gsohler

My two main goals between now and end of 2018† are:

  • A hierarchical namespace for standard libraries. I want to add lots of new libraries, and provide a mechanism for user-contributed standard libraries. This is close to being ready.
  • Graphical sliders. This will involve new language syntax (for specifying parameters that can be controlled by sliders), a new GLSL compiler (that can output GLSL "uniform variables" for slider parameters), and a GUI component (graphical sliders in the Viewer window). This is a big job, and the coding hasn't started yet.

You are correct, the image will animate in real time as the slider is manipulated.

(†) I work slowly, and all dates are aspirational.

doug-moen avatar Sep 19 '18 14:09 doug-moen

Just a passing comment because I've been playing around with this stuff.

Knighty's "fold and cut polyhedra" technique.

Is this pretty much the same as https://en.wikipedia.org/wiki/Wythoff_construction ? (aside from - I suspect - only supporting convex polyhedra i.e. integer Schwarz triangles)

andybak avatar Sep 23 '18 16:09 andybak

@andybak: For a deep dive into the math behind Knighty's method, see: https://syntopia.github.io/Polytopia/polytopes.html

Syntopia doesn't mention the Wythoff construction in his essay, and I don't know how the techniques are related. Let us know what you figure out. And consider contributing a Curv library for building polyhedra using this math.

doug-moen avatar Sep 23 '18 16:09 doug-moen

There is now a way to contribute built-in libraries using the new lib mechanism.

At some point, I'd like to have a lib.polyhedron library that provides a flexible set of primitives for constructing symmetrical polyhedrons. Possibly based on Knighty's "fold and cut polyhedra" method.

I mentioned that there is no support for static arrays of data in a shape dist function. I plan to fix this for release 1.0. However, I don't think this feature is needed for constructing polyhedra.

doug-moen avatar Sep 30 '18 22:09 doug-moen