florian
/
triangle_group_limit_set
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Compare commits

base: florian:public
Branches Tags
florian:public
florian:master
florian:asymmetric_distance
..
compare: florian:master
Branches Tags
florian:public
florian:master
florian:asymmetric_distance

3 Commits
public ... master

Author SHA1 Message Date
Florian Stecker 78769593a7 movie command line arguments 2023-02-05 15:10:16 -05:00
Florian Stecker a8b4bb7c2c implement rotation group and play around with it a bit 2022-07-29 14:38:48 +09:00
Florian Stecker 35932782e9 new .gitignore 2022-03-02 10:53:19 -06:00
8 changed files with 1040 additions and 526 deletions
Show Stats Expand all files Collapse all files

3
.gitignore vendored Normal file
View File

@ -0,0 +1,3 @@
*.o
limit_set
*.pdf

50
README.md
View File

@ -1,50 +0,0 @@
# triangle_group_limit_set - visualizer for triangle group limit sets in the projective plane
This program visualizes the fractal limit set of triangle groups in SL(3,R) in the projective plane. For the mathematical details, see [our paper].
## Installation & Prerequisites ##
This program is written C and only runs on Linux with the X window system (because that's what I happen to be using). It uses [Cairo] for drawing and X11 for input. We also require the [Gnu Scientific Library (GSL)] including CBLAS.
To build, just make sure GCC as well as these libraries are installed, and run
make
Then we can run it for example like this:
./limit_set 5 5 5 1 1 1 1.0
The arguments are 6 integers (p1, p2, p3, q1, q2, q3) and a floating point numbers t. Together they discribe the triangle group under consideration (t can be changed with the arrow keys).
## Key bindings
Drag the mouse to move the image, drag with Shift pressed to rotate.
| Key | function |
|----------|-----------------------------------------------------------------------------------------------------|
| PageUp | increase t by 2% |
| PageDown | decrease t by 2% |
| Right | increase t by 0.2% |
| Left | decreaes t by 0.2% |
| Up | increase t by 0.002% |
| Down | decrease t by 0.002% |
| Space | (reset t to original value) |
| R | cycle through affine charts |
| l | show limit curve |
| L | show limit curve as line or as points |
| f | generate limit curve using attracting / repelling fixed points of conjugates of the Coxeter element |
| d | show dual limit curve |
| r | show fixed points and lines of generating reflections |
| a | show fixed points and lines of Coxeter elements |
| c | show the orbit of an arbitrary point (point can be changed with Shift+Click) |
| b | show certain conics touching the limit curve |
| B | show the "boxes" used to approximate the limit curve |
| p | save a screenshot of the current image (in PDF format) |
| t | toggle info text |
| i | (print some info to terminal?) |
| x | (show rotated reflectors?) |
| M | (make a movie) |
[Cairo]: https://cairographics.org/
[Gnu Scientific Library (GSL)]: https://www.gnu.org/software/gsl/
[our paper]: https://arxiv.org/abs/2106.11349

644
draw.c
View File

@ -46,21 +46,22 @@ int fixedPoints(DrawingContext *ctx, const char *word, vector_t *out)
{
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix *ev = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(tmp);
for(int i = 0; i < strlen(word); i++) {
if(word[i] == ' ')
continue;
multiply_right(tmp, gen[word[i]-'a'], ctx->ws);
if(word[i] >= 'a' && word[i] <= 'c')
multiply_right(tmp, gen[word[i]-'a'], ctx->ws);
else if(word[i] >= 'A' && word[i] <= 'C')
multiply_right(tmp, gen[word[i]-'A'+3], ctx->ws);
}
int count = real_eigenvectors(tmp, ev, ctx->ws);
LOOP(i) LOOP(j) out[i].x[j] = gsl_matrix_get(ev, j, i);
releaseTempMatrices(ctx->ws, 5);
releaseTempMatrices(ctx->ws, 8);
return count;
}
@ -69,9 +70,9 @@ int wordEigenvalues(DrawingContext *ctx, const char *word, double *out)
{
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_vector *ev = getTempVector(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(tmp);
for(int i = 0; i < strlen(word); i++) {
@ -83,7 +84,7 @@ int wordEigenvalues(DrawingContext *ctx, const char *word, double *out)
LOOP(i) out[i] = gsl_vector_get(ev, i);
releaseTempMatrices(ctx->ws, 4);
releaseTempMatrices(ctx->ws, 7);
releaseTempVectors(ctx->ws, 1);
return count;
@ -96,10 +97,20 @@ void drawPoint(DrawingContext *ctx, point_t p)
cairo_t *C = ctx->cairo;
cairo_save(C);
cairo_arc(C, p.x, p.y, 4.0/ctx->dim->scalefactor, 0, 2*M_PI);
cairo_arc(C, p.x, p.y, 3.0/ctx->dim->scalefactor, 0, 2*M_PI);
cairo_close_path(C);
cairo_fill(C);
cairo_restore(C);
/*
cairo_save(C);
cairo_move_to(C, p.x, p.y);
cairo_close_path(C);
cairo_set_line_cap(C, CAIRO_LINE_CAP_ROUND);
cairo_set_line_width(C, 10.0/ctx->dim->scalefactor);
cairo_stroke(C);
cairo_restore(C);
*/
}
void drawSegment2d(DrawingContext *ctx, point_t a, point_t b)
@ -183,13 +194,22 @@ void drawSegmentWith(DrawingContext *ctx, vector_t a, vector_t b, vector_t line,
LOOP(i) x[i] = 0.0;
LOOP(i) LOOP(j) {
cofactor = gsl_matrix_get(ctx->cob, (i+1)%3, (j+1)%3) * gsl_matrix_get(ctx->cob, (i+2)%3, (j+2)%3)
- gsl_matrix_get(ctx->cob, (i+1)%3, (j+2)%3) * gsl_matrix_get(ctx->cob, (i+2)%3, (j+1)%3);
- gsl_matrix_get(ctx->cob, (i+1)%3, (j+2)%3) * gsl_matrix_get(ctx->cob, (i+2)%3, (j+1)%3);
x[i] += cofactor * line.x[j];
}
// t = parameter on segment of intersection with line, s(t) = a + (b-a)*t
tline = -(a_.x*x[0] + a_.y*x[1] + x[2])/((b_.x - a_.x)*x[0] + (b_.y - a_.y)*x[1]);
/*
printf("tline: %f\n", tline);
point_t s;
s.x = a_.x - (b_.x-a_.x)*tline;
s.y = a_.y - (b_.y-a_.y)*tline;
drawPoint(ctx, s);
*/
if((tline < 0 || tline > 1) != contains) {
drawSegment2d(ctx, a_, b_);
} else { // need to draw complementary segment
@ -197,7 +217,6 @@ void drawSegmentWith(DrawingContext *ctx, vector_t a, vector_t b, vector_t line,
m.x = ctx->dim->center_x;
m.y = ctx->dim->center_y;
r = ctx->dim->radius;
// equation is coeff2 t^2 + 2 coeff1 t + coeff0 = 0
coeff0 = (a_.x - m.x)*(a_.x - m.x) + (a_.y - m.y)*(a_.y - m.y) - r*r;
coeff1 = (a_.x - m.x)*(b_.x - a_.x) + (a_.y - m.y)*(b_.y - a_.y);
@ -218,6 +237,7 @@ void drawSegmentWith(DrawingContext *ctx, vector_t a, vector_t b, vector_t line,
}
// level 3: boxes and polygons
void drawPolygon(DrawingContext *ctx, int segments, int sides, ...)
{
va_list args;
@ -308,13 +328,13 @@ void drawRotationOrbitFrame(DrawingContext *ctx, gsl_matrix *frame, vector_t sta
vector_t v[3], w;
point_t p;
double parameter, startangle;
int iterations = 200;
int iterations = 2000;
gsl_matrix *inverse = getTempMatrix(ctx->ws);
gsl_vector *start_v = getTempVector(ctx->ws);
gsl_vector *start_in_frame = getTempVector(ctx->ws);
cairo_t *C = ctx->cairo;
// computeRotationMatrix(ctx, frame, word);
// computeRotationMatrixFrame(ctx, frame, word);
LOOP(i) LOOP(j) v[i].x[j] = gsl_matrix_get(frame, j, i);
LOOP(i) gsl_vector_set(start_v, i, start.x[i]);
@ -350,7 +370,7 @@ void drawRotationOrbit(DrawingContext *ctx, const char *word, vector_t start)
{
gsl_matrix *frame = getTempMatrix(ctx->ws);
computeRotationMatrix(ctx, frame, word);
computeRotationMatrixFrame(ctx, frame, word);
drawRotationOrbitFrame(ctx, frame, start);
releaseTempMatrices(ctx->ws, 1);
@ -368,7 +388,7 @@ void drawDualRotationOrbit(DrawingContext *ctx, const char *word, vector_t start
gsl_vector *start_in_frame = getTempVector(ctx->ws);
cairo_t *C = ctx->cairo;
computeRotationMatrix(ctx, frame, word);
computeRotationMatrixFrame(ctx, frame, word);
LOOP(i) LOOP(j) v[i].x[j] = gsl_matrix_get(frame, j, i);
LOOP(i) gsl_vector_set(start_v, i, start.x[i]);
@ -413,7 +433,7 @@ void drawArcWithOutput(DrawingContext *ctx, const char *word, vector_t start, ve
gsl_vector *vector_in_frame = getTempVector(ctx->ws);
cairo_t *C = ctx->cairo;
computeRotationMatrix(ctx, frame, word);
computeRotationMatrixFrame(ctx, frame, word);
LOOP(i) LOOP(j) v[i].x[j] = gsl_matrix_get(frame, j, i);
LOOP(i) gsl_vector_set(vector, i, start.x[i]);
@ -555,11 +575,10 @@ void drawAttractors(DrawingContext *ctx)
vector_t p[6][3];
vector_t l[6][3];
fixedPoints(ctx, "abc", p[0]);
fixedPoints(ctx, "bca", p[1]);
fixedPoints(ctx, "cab", p[2]);
// fixedPoints(ctx, "abacacbabc", p[3]);
// fixedPoints(ctx, "abcabcabcabcab", p[3]);
fixedPoints(ctx, "cba", p[0]);
fixedPoints(ctx, "bac", p[1]);
fixedPoints(ctx, "acb", p[2]);
fixedPoints(ctx, "a cab a", p[3]);
fixedPoints(ctx, "b abc b", p[4]);
fixedPoints(ctx, "c bca c", p[5]);
@ -633,6 +652,12 @@ void drawCurvedBox(DrawingContext *ctx, int base, const char *conj, int style)
fixedPoints(ctx, word, p[10]);
// conjugate_word("bca b abc b acb", modifier, conj, word);
// fixedPoints(ctx, word, p[6]);
// conjugate_word("bca baca cab acab acb", modifier, conj, word);
// fixedPoints(ctx, word, p[7]);
LOOP(j) l[0][j] = cross(p[0][(3-j)%3], p[0][(4-j)%3]);
LOOP(j) l[1][j] = cross(p[1][(3-j)%3], p[1][(4-j)%3]);
LOOP(j) l[10][j] = cross(p[10][(3-j)%3], p[10][(4-j)%3]);
@ -679,68 +704,233 @@ void drawCurvedBox(DrawingContext *ctx, int base, const char *conj, int style)
}
}
groupelement_t *left(const char *word, groupelement_t *g)
{
int n = strlen(word);
for(int i = n-1; i >= 0; i--) {
if(word[i] == ' ')
continue;
g = g->adj[word[i]-'a'];
if(!g)
break;
}
return g;
}
void drawBoxes(DrawingContext *ctx)
{
gsl_matrix *rot = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
gsl_matrix *frame = getTempMatrix(ctx->ws);
gsl_matrix *frame2 = getTempMatrix(ctx->ws);
gsl_vector *startpoint_drawbasis = getTempVector(ctx->ws);
gsl_vector *startpoint_globalbasis = getTempVector(ctx->ws);
gsl_matrix **elements = getTempMatrices(ctx->ws, ctx->n_group_elements);
cairo_t *C = ctx->cairo;
cairo_save(C);
vector_t p[22][3];
vector_t l[22][3];
vector_t alpha[6];
vector_t ptmp[3];
char word[100], word2[100];
fixedPoints(ctx, "abc", p[0]);
fixedPoints(ctx, "bca", p[1]);
fixedPoints(ctx, "cab", p[2]);
fixedPoints(ctx, "bacabab", p[3]);
fixedPoints(ctx, "bcacabacb", p[4]);
cairo_set_line_width(C, 2.0/ctx->dim->scalefactor);
cairo_set_source_rgb(C, 0.6, 0.6, 0.6);
drawRotationOrbit(ctx, "ab", p[0][0]);
drawRotationOrbit(ctx, "bc", p[0][0]);
drawRotationOrbit(ctx, "ca", p[0][0]);
vector_t p[22][3];
vector_t fp[3];
vector_t l[22][3];
vector_t alpha[6];
vector_t ptmp[3];
vector_t start;
vector_t start2;
char word[100], word2[100];
if(ctx->mode >= 2) {
drawRotationOrbit(ctx, "bcabcb", p[1][0]); // bcC
drawRotationOrbit(ctx, "abcabcba", p[0][0]); // abcC
fixedPoints(ctx, "cba", p[0]);
fixedPoints(ctx, "acb", p[1]);
fixedPoints(ctx, "bac", p[2]);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(elements[0]);
for(int i = 1; i < ctx->n_group_elements; i++) {
if(ctx->group[i].length % 2)
continue;
int letter = ROTATION_LETTER(ctx->group[i].letter, ctx->group[i].parent->letter);
multiply(gen[letter], elements[ctx->group[i].parent->parent->id], elements[i]);
}
cairo_set_source_rgb(C, 0, 0, 0);
fixedPoints(ctx, "abababcbaba", p[3]);
gsl_vector_set(startpoint_drawbasis, 0, ctx->marking.x);
gsl_vector_set(startpoint_drawbasis, 1, ctx->marking.y);
gsl_vector_set(startpoint_drawbasis, 2, 1);
solve(ctx->cob, startpoint_drawbasis, startpoint_globalbasis, ctx->ws);
LOOP(i) start.x[i] = gsl_vector_get(startpoint_globalbasis, i);
cairo_set_source_rgb(C, 0, 0, 1);
fixedPoints(ctx, "bc bca cb", p[3]);
/*
fixedPoints(ctx, "ABAB", fp);
drawCovector(ctx, cross(fp[0], fp[2]));
drawVector(ctx, fp[1]);
computeRotationMatrixFrame(ctx, frame, "c");
drawRotationOrbit(ctx, "c", start); // ba cb ab ac = C A c B
drawVector(ctx, start);
computeMatrix(ctx, tmp, "ABABc");
start2 = apply(tmp, start);
computeRotationMatrixFrame(ctx, frame, "c");
multiply(tmp, frame, frame2);
drawVector(ctx, start2);
drawRotationOrbitFrame(ctx, frame2, apply(tmp, start));
computeMatrix(ctx, tmp, "ABABcABABC");
start2 = apply(tmp, start);
computeRotationMatrixFrame(ctx, frame, "c");
multiply(tmp, frame, frame2);
drawVector(ctx, start2);
drawRotationOrbitFrame(ctx, frame2, apply(tmp, start));
*/
for(groupelement_t *cur = &ctx->group[180]; cur->parent; cur = cur->parent)
fputc('a'+cur->letter, stdout);
fputc('\n', stdout);
queue_t queue;
queue_init(&queue);
queue_put(&queue, 0);
int current;
groupelement_t *cur, *next;
for(int i = 0; i < ctx->n_group_elements_combinatorial; i++)
ctx->group[i].visited = 0;
cur = &ctx->group[0];
computeRotationMatrixFrame(ctx, frame, "c");
for(int i = 0; i < 1000; i++) {
if(ctx->group[i].length % 2 != 0)
continue;
multiply(elements[i], frame, frame2);
drawRotationOrbitFrame(ctx, frame2, apply(elements[i], start));
}
/*
while((current = queue_get(&queue)) != -1) {
cur = &ctx->group[current];
if(cur->visited > 4)
continue;
if(cur->id < ctx->n_group_elements) {
multiply(elements[cur->id], frame, frame2);
drawRotationOrbitFrame(ctx, frame2, apply(elements[cur->id], start));
}
next = left("ab ab", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
next = left("cbac cbac", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
next = left("cacbca cacbca", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
next = left("cabcacbcac cabcacbcac", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
next = left("acbcacba acbcacba", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
next = left("bcacbc bcacbc", cur);
if(next && next->visited == 0) {
queue_put(&queue, next->id);
next->visited = cur->visited+1;
}
}
*/
// drawRotationOrbit(ctx, "a", p[1][0]);
// drawRotationOrbit(ctx, "b", p[2][0]);
cairo_restore(C);
releaseTempMatrices(ctx->ws, 5);
releaseTempMatrices(ctx->ws, 9 + ctx->n_group_elements);
releaseTempVectors(ctx->ws, 2);
}
void drawBoxes2(DrawingContext *ctx)
{
gsl_matrix *rot = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
cairo_t *C = ctx->cairo;
cairo_save(C);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
vector_t p[4][3];
fixedPoints(ctx, "abc", p[0]);
fixedPoints(ctx, "bca", p[1]);
fixedPoints(ctx, "cab", p[2]);
fixedPoints(ctx, "cba", p[0]);
fixedPoints(ctx, "acb", p[1]);
fixedPoints(ctx, "bac", p[2]);
cairo_set_line_width(C, 2.5/ctx->dim->scalefactor);
/*
cairo_set_source_rgb(C, 0, 1, 1);
fixedPoints(ctx, "abc", p[0]);
drawRotationOrbit(ctx, "bc", p[0][0]); // (r^{-1}a)^2 C / cbr
fixedPoints(ctx, "bca", p[0]);
drawRotationOrbit(ctx, "abca", p[0][0]); // (r^{-1}a)^2 C / cbr
cairo_set_source_rgb(C, 1, 1, 0);
fixedPoints(ctx, "abc", p[0]);
drawRotationOrbit(ctx, "cbcabc", p[0][0]); // (r^{-1}a)^4 C / cbac r-
fixedPoints(ctx, "bca", p[0]);
drawRotationOrbit(ctx, "acbcabca", p[0][0]); // (r^{-1}a)^4 C / cbac r-
cairo_set_source_rgb(C, 1, 0, 1);
fixedPoints(ctx, "abc", p[0]);
drawRotationOrbit(ctx, "cbacabcabc", p[0][0]); // (r^{-1}a)^6 C / cbacba
fixedPoints(ctx, "bca", p[0]);
drawRotationOrbit(ctx, "acbacabcabca", p[0][0]); // (r^{-1}a)^6 C / cbacba
*/
/*
cairo_set_source_rgb(C, 1, 0, 0);
fixedPoints(ctx, "bca", p[3]);
drawRotationOrbit(ctx, "bcabcb", p[3][0]); // cb c a bc
fixedPoints(ctx, "abc", p[3]);
drawRotationOrbit(ctx, "abcabcba", p[3][0]); // cb c a bc
cairo_set_source_rgb(C, 0, 1, 0);
fixedPoints(ctx, "ababcba", p[3]);
drawRotationOrbit(ctx, "ababcababa", p[3][0]); // cb c a bc
fixedPoints(ctx, "abc", p[3]);
drawRotationOrbit(ctx, "abcabcabacba", p[3][0]); // cb c a bc
*/
// drawRotationOrbit(ctx, "bc", p[3][0]);
cairo_set_source_rgb(C, 0, 0, 0);
drawCurvedBox(ctx, 'A', "", 2);
cairo_set_source_rgb(C, 0, 0, 1);
// drawCurvedBox(ctx, 'C', "ab", 2);
// drawCurvedBox(ctx, 'C', "ab", 2);
// drawCurvedBox(ctx, 'B', "bca", 2);
// drawCurvedBox(ctx, 'B', "ba", 2);
// drawCurvedBox(ctx, 'B', "bca", 2);
// drawCurvedBox(ctx, 'B', "abca", 2);
// drawCurvedBox(ctx, 'C', "abab", 2);
// drawCurvedBox(ctx, 'C', "ababab", 2);
if(ctx->mode >= 3 && ctx->mode != 4) {
cairo_set_source_rgb(C, 0, 0, 0);
drawCurvedBox(ctx, 'B', "", 2);
@ -756,11 +946,15 @@ void drawBoxes2(DrawingContext *ctx)
drawCurvedBox(ctx, 'A', "bcbcbcbc", 2);
drawCurvedBox(ctx, 'A', "bcbcbcbcbc", 2);
drawCurvedBox(ctx, 'A', "bcbcbcbcbcbc", 2);
// drawCurvedBox(ctx, 'A', "b", 2);
// drawCurvedBox(ctx, 'A', "bcb", 2);
// drawCurvedBox(ctx, 'A', "bcbcb", 2);
// drawCurvedBox(ctx, 'A', "bcbcbcb", 2);
// drawCurvedBox(ctx, 'A', "bcbcbcbcb", 2);
}
if(ctx->mode == 6) {
cairo_set_source_rgb(C, 0, 0.8, 0.5);
drawCurvedBox(ctx, 'C', "ababab", 2);
drawCurvedBox(ctx, 'C', "abababab", 2);
drawCurvedBox(ctx, 'C', "ababababab", 2);
}
@ -769,45 +963,291 @@ void drawBoxes2(DrawingContext *ctx)
cairo_set_source_rgb(C, 0, 0.8, 0.5);
drawCurvedBox(ctx, 'C', "ab", 2);
drawCurvedBox(ctx, 'C', "abab", 2);
drawCurvedBox(ctx, 'C', "ababab", 2);
}
if(ctx->mode >= 8) {
cairo_set_source_rgb(C, 0, 0.8, 0.5);
drawCurvedBox(ctx, 'C', "b", 2);
drawCurvedBox(ctx, 'C', "bab", 2);
drawCurvedBox(ctx, 'C', "babab", 2);
}
if(ctx->mode >= 9) {
cairo_set_source_rgb(C, 0.8, 0, 0.5);
drawCurvedBox(ctx, 'B', "abca", 2);
drawCurvedBox(ctx, 'B', "abcaca", 2);
drawCurvedBox(ctx, 'B', "aba", 2);
drawCurvedBox(ctx, 'B', "abaca", 2);
cairo_set_source_rgb(C, 0.8, 0, 0.5);
drawCurvedBox(ctx, 'B', "abca", 2);
drawCurvedBox(ctx, 'B', "abcaca", 2);
drawCurvedBox(ctx, 'B', "abcacaca", 2);
drawCurvedBox(ctx, 'B', "aba", 2);
drawCurvedBox(ctx, 'B', "abaca", 2);
drawCurvedBox(ctx, 'B', "abacaca", 2);
drawCurvedBox(ctx, 'B', "ababca", 2);
drawCurvedBox(ctx, 'B', "ababcaca", 2);
drawCurvedBox(ctx, 'B', "ababa", 2);
drawCurvedBox(ctx, 'B', "ababaca", 2);
drawCurvedBox(ctx, 'B', "ababca", 2);
drawCurvedBox(ctx, 'B', "ababcaca", 2);
drawCurvedBox(ctx, 'B', "ababcacaca", 2);
drawCurvedBox(ctx, 'B', "ababa", 2);
drawCurvedBox(ctx, 'B', "ababaca", 2);
drawCurvedBox(ctx, 'B', "ababacaca", 2);
drawCurvedBox(ctx, 'B', "bca", 2);
drawCurvedBox(ctx, 'B', "bcaca", 2);
drawCurvedBox(ctx, 'B', "ba", 2);
drawCurvedBox(ctx, 'B', "baca", 2);
drawCurvedBox(ctx, 'B', "abababca", 2);
drawCurvedBox(ctx, 'B', "abababcaca", 2);
drawCurvedBox(ctx, 'B', "abababcacaca", 2);
drawCurvedBox(ctx, 'B', "abababa", 2);
drawCurvedBox(ctx, 'B', "abababaca", 2);
drawCurvedBox(ctx, 'B', "abababacaca", 2);
drawCurvedBox(ctx, 'B', "babca", 2);
drawCurvedBox(ctx, 'B', "babcaca", 2);
drawCurvedBox(ctx, 'B', "baba", 2);
drawCurvedBox(ctx, 'B', "babaca", 2);
drawCurvedBox(ctx, 'B', "bca", 2);
drawCurvedBox(ctx, 'B', "bcaca", 2);
drawCurvedBox(ctx, 'B', "bcacaca", 2);
drawCurvedBox(ctx, 'B', "ba", 2);
drawCurvedBox(ctx, 'B', "baca", 2);
drawCurvedBox(ctx, 'B', "bacaca", 2);
cairo_restore(C);
releaseTempMatrices(ctx->ws, 4);
drawCurvedBox(ctx, 'B', "babca", 2);
drawCurvedBox(ctx, 'B', "babcaca", 2);
drawCurvedBox(ctx, 'B', "babcacaca", 2);
drawCurvedBox(ctx, 'B', "baba", 2);
drawCurvedBox(ctx, 'B', "babaca", 2);
drawCurvedBox(ctx, 'B', "babacaca", 2);
drawCurvedBox(ctx, 'B', "bababca", 2);
drawCurvedBox(ctx, 'B', "bababcaca", 2);
drawCurvedBox(ctx, 'B', "bababcacaca", 2);
drawCurvedBox(ctx, 'B', "bababa", 2);
drawCurvedBox(ctx, 'B', "bababaca", 2);
drawCurvedBox(ctx, 'B', "bababacaca", 2);
}
// drawCurvedBox(ctx, 'B', "abab", 2);
// drawCurvedBox(ctx, 'B', "ababab", 2);
// drawCurvedBox(ctx, 'C', "b", 1);
// drawCurvedBox(ctx, 'C', "ababab");
// drawCurvedBox(ctx, 'C', "abababab");
// drawCurvedBox(ctx, 'C', "b");
// drawCurvedBox(ctx, 'C', "bab");
// drawCurvedBox(ctx, 'C', "babab");
// drawCurvedBox(ctx, 'C', "abababab");
// drawCurvedBox(ctx, 'C', "ababababab");
// drawCurvedBox(ctx, 'C', "abababababab");
cairo_set_source_rgb(C, 0, 0, 1);
// drawCurvedBox(ctx, 'B', "abca", 1);
// drawCurvedBox(ctx, 'B', "aba", 1);
// drawCurvedBox(ctx, 'B', "bca", 1);
// drawCurvedBox(ctx, 'B', "ba", 1);
// drawCurvedBox(ctx, 'B', "abaca");
// drawCurvedBox(ctx, 'B', "abcaca");
// drawCurvedBox(ctx, 'B', "abaca");
// drawCurvedBox(ctx, 'B', "aba");
// drawCurvedBox(ctx, 'B', "ababca");
// drawCurvedBox(ctx, 'B', "ababa");
// drawCurvedBox(ctx, 'B', "abababca");
// drawCurvedBox(ctx, 'B', "abababa");
// drawCurvedBox(ctx, 'B', "bca");
// drawCurvedBox(ctx, 'B', "ba");
// drawCurvedBox(ctx, 'B', "babca");
// drawCurvedBox(ctx, 'B', "baba");
// drawCurvedBox(ctx, 'B', "bababca");
// drawCurvedBox(ctx, 'B', "bababa");
cairo_set_source_rgb(C, 1, 0, 1);
// drawCurvedBox(ctx, 'A', "abcac");
// drawCurvedBox(ctx, 'A', "bcac");
// drawCurvedBox(ctx, 'A', "abcacbc");
// drawCurvedBox(ctx, 'A', "abcac");
// drawCurvedBox(ctx, 'A', "abcacbc");
// drawCurvedBox(ctx, 'A', "ababc");
// drawCurvedBox(ctx, 'A', "ababcbc");
/* drawCurvedBox(ctx, 'A', "abcac");
drawCurvedBox(ctx, 'A', "ababc");
drawCurvedBox(ctx, 'A', "abac");
drawCurvedBox(ctx, 'A', "ababcabc");
drawCurvedBox(ctx, 'A', "ababcac");
drawCurvedBox(ctx, 'A', "abababc");
drawCurvedBox(ctx, 'A', "ababac");
drawCurvedBox(ctx, 'A', "abababcabc");
drawCurvedBox(ctx, 'A', "abababcac");
drawCurvedBox(ctx, 'A', "ababababc");
drawCurvedBox(ctx, 'A', "abababac");
drawCurvedBox(ctx, 'A', "bcabc");
drawCurvedBox(ctx, 'A', "bcac");
drawCurvedBox(ctx, 'A', "babc");
drawCurvedBox(ctx, 'A', "bac");
drawCurvedBox(ctx, 'A', "babcabc");
drawCurvedBox(ctx, 'A', "babcac");
drawCurvedBox(ctx, 'A', "bababc");
drawCurvedBox(ctx, 'A', "babac");
drawCurvedBox(ctx, 'A', "bababcabc");
drawCurvedBox(ctx, 'A', "bababcac");
drawCurvedBox(ctx, 'A', "babababc");
drawCurvedBox(ctx, 'A', "bababac");
*/
// drawCurvedBox(ctx, 'A', "ababc");
cairo_set_source_rgb(C, 1, 0.5, 0);
// drawCurvedBox(ctx, 'C', "abcabcab");
// drawCurvedBox(ctx, 'B', "bca");
// drawCurvedBox(ctx, 'B', "abcacaca");
// drawCurvedBox(ctx, 'B', "abaca");
// drawCurvedBox(ctx, 'B', "aba");
// drawCurvedBox(ctx, 'B', "aba");
// drawCurvedBox(ctx, 'B', "abcaca");
// drawCurvedBox(ctx, 'B', "aba");
// drawCurvedBox(ctx, 'B', "abaca");
cairo_set_source_rgb(C, 0, 0, 1);
// drawCurvedBox(ctx, 'A', "abcabc");
// drawCurvedBox(ctx, 'A', "abcabcbc");
// drawCurvedBox(ctx, 'A', "abcacbc");
// drawCurvedBox(ctx, 'A', "abcac");
// drawCurvedBox(ctx, 'A', "bcabc");
// drawCurvedBox(ctx, 'A', "bcabcbc");
// drawCurvedBox(ctx, 'A', "bcacbc");
// drawCurvedBox(ctx, 'A', "bcac");
// drawCurvedBox(ctx, 'A', "abcac");
// drawCurvedBox(ctx, 'A', "abcac");
// drawCurvedBox(ctx, 'A', "abcacbc");
// drawCurvedBox(ctx, 'B', "abacaca");
// drawCurvedBox(ctx, 'B', "abacacaca");
cairo_set_source_rgb(C, 0, 1, 1);
fixedPoints(ctx, "ababcba", p[0]);
// drawRotationOrbit(ctx, "abcaba", p[0][0]);
/*
int p = 9;
vector_t fp[3][3],neutral_line[3],reflection_line[p],star[2*p],outer[2];
vector_t rotation_line = {0,0,1};
cairo_set_line_width(C, 1.5/ctx->dim->scalefactor);
cairo_set_source_rgb(C, 1, 0, 0);
multiply(gen[0], gen[1], rot);
fixedPoints(ctx, "abc", fp[0]);
fixedPoints(ctx, "bca", fp[1]);
fixedPoints(ctx, "cab", fp[2]);
LOOP(i) neutral_line[i] = cross(fp[i][0], fp[i][2]);
LOOP(j) reflection_line[0].x[j] = gsl_matrix_get(ctx->cartan, 0, j);
star[0] = cross(neutral_line[0],reflection_line[0]);
star[p] = cross(neutral_line[2],reflection_line[0]);
for(int j = 1; j < p; j++) {
reflection_line[j] = apply_transpose(rot, reflection_line[j-1]);
star[j] = apply(rot, star[j-1]);
star[j+p] = apply(rot, star[j+p-1]);
}
outer[0] = cross(neutral_line[0],reflection_line[5]);
outer[1] = cross(neutral_line[2],reflection_line[5]);
for(int j = 0; j < 8; j++) {
drawVector(ctx, star[j]);
drawVector(ctx, star[p+j]);
}
for(int j = 0; j < 7; j++) {
// if(j == 3)
// continue;
drawSegment(ctx, star[j%p], star[(j+1)%p+p]);
drawSegment(ctx, star[(j+1)%p+p], star[(j+1)%p]);
drawSegment(ctx, star[(j+1)%p], star[j%p+p]);
drawSegment(ctx, star[j%p+p], star[j%p]);
}
cairo_set_source_rgb(C, 0, 0, 1);
drawVector(ctx,outer[0]);
drawVector(ctx,outer[1]);
// drawCovector(ctx, neutral_line[0]);
// drawCovector(ctx, neutral_line[2]);
drawSegment(ctx, star[0], star[p]);
drawSegment(ctx, star[p], outer[1]);
drawSegment(ctx, outer[1], outer[0]);
drawSegment(ctx, outer[0], star[0]);
cairo_set_source_rgb(C, 0, 0, 0);
drawCovector(ctx, rotation_line);
*/
// for(int j = 0; j < 5; j++)
// drawCovector(ctx, reflection_line[j]);
// intersect attracting line with neutral line of the other element
// for(int j = 0; j < 2; j++)
// i[j] = cross(cross(p[j%2][0],p[j%2][1]),cross(p[(j+1)%2][0],p[(j+1)%2][2]));
// drawPolygon(ctx, 1, 4, p[0][0], i[0], p[1][0], i[1]);
// abc, ababcba, abababcbaba, ..., cab
// bca, acaba, abacababa, ..., babcb
/*
vector_t v[4][3];
vector_t i[4];
fixedPoints(ctx, "abc", v[0]);
fixedPoints(ctx, "bca", v[1]);
fixedPoints(ctx, "cab", v[2]);
fixedPoints(ctx, "acaba", v[3]);
i[0] = cross(cross(v[0][1],v[0][2]),cross(v[2][0],v[2][2]));
i[1] = cross(cross(v[1][1],v[1][2]),cross(v[3][0],v[3][2]));
i[2] = cross(cross(v[2][0],v[2][2]),cross(v[3][0],v[3][2]));
i[3] = cross(cross(v[0][0],v[0][2]),cross(v[1][0],v[1][2]));
// cairo_set_source_rgb(C, 0, 0, 1);
// drawPolygon(ctx, 1, 6, v[2][2], v[1][1], v[0][1], v[3][2], v[3][1], v[2][1]);
// drawPolygon(ctx, 1, 6, v[1][2], i[1], i[2], i[0], v[0][2], i[3]);
cairo_set_line_width(C, 1.5/ctx->dim->scalefactor);
cairo_set_source_rgb(C, 1, 0, 0);
// drawBox(ctx, "bca", "abc");
drawCurvedBox(ctx, 'A', "");
*/
/*
cairo_set_source_rgb(C, 1, 0.5, 0);
drawBox(ctx, "ab abc ba", "ab cab ba");
drawBox(ctx, "abab abc baba", "abab cab baba");
drawBox(ctx, "ababab abc bababa", "ababab cab bababa");
drawBox(ctx, "abababab abc babababa", "abababab cab babababa");
// drawBox(ctx, "ababababab abc bababababa", "ababababab cab bababababa");
drawBox(ctx, "b abc b", "b cab b");
drawBox(ctx, "bab abc bab", "bab cab bab");
drawBox(ctx, "babab abc babab", "babab cab babab");
drawBox(ctx, "bababab abc bababab", "bababab cab bababab");
// drawBox(ctx, "babababab abc babababab", "babababab cab babababab");
*/
cairo_restore(C);
releaseTempMatrices(ctx->ws, 7);
}
void drawRotatedReflectors(DrawingContext *ctx)
{
gsl_matrix *rot = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
cairo_t *C = ctx->cairo;
vector_t fp[3], fp2[3];
vector_t w;
@ -815,7 +1255,7 @@ void drawRotatedReflectors(DrawingContext *ctx)
cairo_set_source_rgb(C, 0.7, 0.7, 0.7);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
LOOP(i) LOOP(j) v[i].x[j] = gsl_matrix_get(ctx->cartan, i, j);
multiply(gen[0], gen[1], rot);
@ -836,7 +1276,7 @@ void drawRotatedReflectors(DrawingContext *ctx)
fixedPoints(ctx, "cacabac", fp2);
drawRotationOrbit(ctx, "ac", fp[0]);
releaseTempMatrices(ctx->ws, 4);
releaseTempMatrices(ctx->ws, 7);
}
void drawDualLimitCurve(DrawingContext *ctx)
@ -846,12 +1286,36 @@ void drawDualLimitCurve(DrawingContext *ctx)
cairo_save(C);
cairo_set_source_rgb(C, 0.5, 0.5, 1);
int n = 18;
vector_t p[n][3];
vector_t l[n][3];
vector_t ptmp[3], ltmp[3];
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
gsl_matrix **elements = getTempMatrices(ctx->ws, ctx->n_group_elements);
fixedPoints(ctx, "abc", p[0]);
// wordEigenvalues(ctx, "abc", ev);
// LOOP(i) LOOP(j) gsl_matrix_set(coxeter_fixedpoints, j, i, cox[0][i].x[j]);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(elements[0]);
for(int i = 1; i < ctx->n_group_elements; i++)
multiply(gen[ctx->group[i].letter], elements[ctx->group[i].parent->id], elements[i]);
vector_t p[3], l[3], v;
fixedPoints(ctx, "cba", p);
drawVector(ctx, p[0]);
drawVector(ctx, p[1]);
drawVector(ctx, p[2]);
LOOP(i) l[i] = cross(p[(i+1)%3], p[(i+2)%3]);
for(int i = 0; i < ctx->n_group_elements; i++) {
v = apply_transpose(elements[i], l[0]);
drawCovector(ctx, v);
}
releaseTempMatrices(ctx->ws, 3 + ctx->n_group_elements);
// releaseTempVectors(ctx->ws, 4);
/*
fixedPoints(ctx, "ab abc ba", p[1]);
fixedPoints(ctx, "abab abc baba", p[2]);
fixedPoints(ctx, "ababab abc bababa", p[3]);
@ -870,6 +1334,7 @@ void drawDualLimitCurve(DrawingContext *ctx)
fixedPoints(ctx, "bababa bca ababab", p[15]);
fixedPoints(ctx, "baba bca abab", p[16]);
fixedPoints(ctx, "ba bca ab", p[17]);
*/
/*
fixedPoints(ctx, "abc", p[0]);
@ -907,6 +1372,7 @@ void drawDualLimitCurve(DrawingContext *ctx)
// drawCovector(ctx, l[i][2]);
}*/
/*
fixedPoints(ctx, "abc", ptmp);
drawCovector(ctx, cross(ptmp[0], ptmp[1]));
fixedPoints(ctx, "bca", ptmp);
@ -919,6 +1385,7 @@ void drawDualLimitCurve(DrawingContext *ctx)
drawCovector(ctx, cross(ptmp[0], ptmp[1]));
fixedPoints(ctx, "acaba", ptmp);
drawCovector(ctx, cross(ptmp[0], ptmp[1]));
*/
cairo_restore(C);
}
@ -974,7 +1441,7 @@ void drawLimitCurve(DrawingContext *ctx)
void drawCoxeterOrbit(DrawingContext *ctx)
{
gsl_matrix *rot = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
gsl_vector *eval = getTempVector(ctx->ws);
gsl_matrix *coxeter_fixedpoints = getTempMatrix(ctx->ws);
gsl_vector *startpoint_coxeterbasis = getTempVector(ctx->ws);
@ -991,7 +1458,7 @@ void drawCoxeterOrbit(DrawingContext *ctx)
int first = 1;
cairo_save(C);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
cairo_set_source_rgb(C, 0, 0, 1);
@ -1002,7 +1469,7 @@ void drawCoxeterOrbit(DrawingContext *ctx)
wordEigenvalues(ctx, "abc", ev);
LOOP(i) LOOP(j) gsl_matrix_set(coxeter_fixedpoints, j, i, cox[0][i].x[j]);
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(elements[0]);
for(int i = 1; i < ctx->n_group_elements; i++)
multiply(gen[ctx->group[i].letter], elements[ctx->group[i].parent->id], elements[i]);
@ -1016,7 +1483,7 @@ void drawCoxeterOrbit(DrawingContext *ctx)
gsl_vector_set(startpoint_drawbasis, 2, 1);
solve(ctx->cob, startpoint_drawbasis, startpoint_globalbasis, ctx->ws);
solve(coxeter_fixedpoints, startpoint_globalbasis, startpoint_coxeterbasis, ctx->ws);
// solve(coxeter_fixedpoints, startpoint_globalbasis, startpoint_coxeterbasis, ctx->ws);
// LOOP(i) start.x[i] = gsl_vector_get(startpoint_coxeterbasis, i);
LOOP(i) start.x[i] = gsl_vector_get(startpoint_globalbasis, i);
@ -1056,7 +1523,7 @@ void drawCoxeterOrbit(DrawingContext *ctx)
// LOOP(i) drawVector(ctx, abcb[i]);
cairo_restore(C);
releaseTempMatrices(ctx->ws, 5 + ctx->n_group_elements);
releaseTempMatrices(ctx->ws, 8 + ctx->n_group_elements);
releaseTempVectors(ctx->ws, 4);
}
@ -1065,7 +1532,7 @@ void drawText(DrawingContext *ctx)
cairo_move_to(ctx->cairo, 15, 30);
cairo_set_source_rgb(ctx->cairo, 0, 0, 0);
char buf[100];
sprintf(buf, "t = exp(%.8f) = %.8f, marking = (%.5f, %.5f)", log(ctx->parameter), ctx->parameter, ctx->marking.x, ctx->marking.y);
sprintf(buf, "t = exp(%.8f) = %.8f, s = exp(%.8f) = %.8f, marking = (%.5f, %.5f)", log(ctx->parameter), ctx->parameter, log(ctx->parameter2), ctx->parameter2, ctx->marking.x, ctx->marking.y);
cairo_show_text(ctx->cairo, buf);
}
@ -1111,6 +1578,7 @@ void draw(DrawingContext *ctx)
if(ctx->show_marking)
{
cairo_set_source_rgb(C, 0, 0, 1);
drawPoint(ctx, ctx->marking);
}
if(ctx->show_coxeter_orbit)
@ -1120,7 +1588,7 @@ void draw(DrawingContext *ctx)
cairo_identity_matrix(C); // text is in screen coordinates
if(ctx->show_text)
drawText(ctx);
drawText(ctx);
cairo_surface_flush(cairo_get_target(C));
}

107
limit_set.c
View File

@ -5,6 +5,7 @@ static int compareAngle(const void *x, const void *y)
return ((double*)x)[2] > ((double*)y)[2] ? 1 : -1;
}
// might need a rewrite
void cartanMatrix(gsl_matrix *cartan, double a1, double a2, double a3, double s)
{
gsl_matrix_set(cartan, 0, 0, 2);
@ -20,11 +21,58 @@ void cartanMatrix(gsl_matrix *cartan, double a1, double a2, double a3, double s)
gsl_matrix_set(cartan, 2, 2, 2);
}
void initializeTriangleGenerators(gsl_matrix **gen, gsl_matrix *cartan)
void initializeTriangleGenerators(gsl_matrix **gen, double a1, double a2, double a3, double s, double t, workspace_t *ws)
{
LOOP(i) gsl_matrix_set_identity(gen[i]);
LOOP(i) LOOP(j) *gsl_matrix_ptr(gen[i], j, j) = -1.0;
LOOP(i) LOOP(j) *gsl_matrix_ptr(gen[i], i, j) += gsl_matrix_get(cartan, i, j);
gsl_matrix *reflection_gen[3];
LOOP(i) {
reflection_gen[i] = gsl_matrix_alloc(3, 3);
gsl_matrix_set_identity(reflection_gen[i]);
}
double rho[3];
rho[0] = sqrt(s*s + 2*s*cos(a1) + 1);
rho[1] = sqrt(s*s + 2*s*cos(a2) + 1);
rho[2] = sqrt(s*s + 2*s*cos(a3) + 1);
gsl_matrix_set(reflection_gen[0], 0, 0, -1.0);
gsl_matrix_set(reflection_gen[0], 0, 1, rho[2]*t);
gsl_matrix_set(reflection_gen[0], 0, 2, rho[1]/t);
gsl_matrix_set(reflection_gen[1], 1, 0, rho[2]/t);
gsl_matrix_set(reflection_gen[1], 1, 1, -1.0);
gsl_matrix_set(reflection_gen[1], 1, 2, rho[0]*t);
gsl_matrix_set(reflection_gen[2], 2, 0, rho[1]*t);
gsl_matrix_set(reflection_gen[2], 2, 1, rho[0]/t);
gsl_matrix_set(reflection_gen[2], 2, 2, -1.0);
LOOP(i) {
gsl_matrix_set_identity(gen[i]);
gsl_matrix_set(gen[i], (i+1)%3, (i+1)%3, s);
gsl_matrix_set(gen[i], (i+2)%3, (i+2)%3, 1/s);
gsl_matrix_set_identity(gen[i+3]);
gsl_matrix_set(gen[i+3], (i+1)%3, (i+1)%3, 1/s);
gsl_matrix_set(gen[i+3], (i+2)%3, (i+2)%3, s);
}
LOOP(i) {
multiply_left(reflection_gen[i], gen[(i+2)%3], ws);
multiply_right(gen[(i+2)%3], reflection_gen[(i+1)%3], ws);
multiply_left(reflection_gen[(i+1)%3], gen[(i+2)%3+3], ws);
multiply_right(gen[(i+2)%3+3], reflection_gen[i], ws);
}
LOOP(i) gsl_matrix_free(reflection_gen[i]);
}
void initializeTriangleGeneratorsCurrent(gsl_matrix **gen, DrawingContext *ctx)
{
double angle[3];
LOOP(i) angle[i] = 2*M_PI*ctx->k[i]/ctx->p[i];
initializeTriangleGenerators(gen, angle[0], angle[1], angle[2], ctx->parameter2, ctx->parameter, ctx->ws);
}
int computeLimitCurve(DrawingContext *ctx)
@ -38,7 +86,7 @@ int computeLimitCurve(DrawingContext *ctx)
gsl_matrix *coxeter = getTempMatrix(ctx->ws);
gsl_matrix *coxeter_fixedpoints = getTempMatrix(ctx->ws);
gsl_matrix *fixedpoints = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
gsl_matrix **elements = getTempMatrices(ctx->ws, ctx->n_group_elements);
groupelement_t *group = ctx->group;
int success = 0;
@ -50,31 +98,43 @@ int computeLimitCurve(DrawingContext *ctx)
// do first in the Fuchsian positive case to get the angles
cartanMatrix(cartan_pos, M_PI/ctx->p[0], M_PI/ctx->p[1], M_PI/ctx->p[2], 1.0);
initializeTriangleGenerators(gen, cartan_pos);
initializeTriangleGenerators(gen, 2*M_PI/ctx->p[0], 2*M_PI/ctx->p[1], 2*M_PI/ctx->p[2], 1.0, 1.0, ctx->ws);
gsl_matrix_set_identity(elements[0]);
for(int i = 1; i < ctx->n_group_elements; i++)
multiply(gen[group[i].letter], elements[group[i].parent->id], elements[i]);
for(int i = 1; i < ctx->n_group_elements; i++) {
if(group[i].length % 2)
continue;
int letter = ROTATION_LETTER(group[i].letter, group[i].parent->letter);
multiply(gen[letter], elements[group[i].parent->parent->id], elements[i]);
}
diagonalize_symmetric_form(cartan_pos, cob_pos, ws);
multiply_many(ws, coxeter_pos, 3, gen[0], gen[1], gen[2]);
multiply_many(ws, coxeter_pos, 3, gen[2], gen[1], gen[0]);
int ev_count_pos = real_eigenvectors(coxeter_pos, coxeter_fixedpoints_pos, ws);
if(ev_count_pos != 3)
goto error_out;
int n = 0;
for(int i = 0; i < ctx->n_group_elements; i++) {
if(group[i].length % 2)
continue;
multiply_many(ws, fixedpoints_pos, 3, cob_pos, elements[i], coxeter_fixedpoints_pos);
ctx->limit_curve[3*i+2] = atan2(
ctx->limit_curve[3*n+2] = atan2(
gsl_matrix_get(fixedpoints_pos, 2, column)/gsl_matrix_get(fixedpoints_pos, 0, column),
gsl_matrix_get(fixedpoints_pos, 1, column)/gsl_matrix_get(fixedpoints_pos, 0, column));
n++;
}
// now do it again to calculate x and y coordinates
initializeTriangleGenerators(gen, ctx->cartan);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(elements[0]);
for(int i = 1; i < ctx->n_group_elements; i++)
multiply(gen[group[i].letter], elements[group[i].parent->id], elements[i]);
for(int i = 1; i < ctx->n_group_elements; i++) {
if(group[i].length % 2)
continue;
int letter = ROTATION_LETTER(group[i].letter, group[i].parent->letter);
multiply(gen[letter], elements[group[i].parent->parent->id], elements[i]);
}
multiply_many(ws, coxeter, 3, gen[0], gen[1], gen[2]);
multiply_many(ws, coxeter, 3, gen[2], gen[1], gen[0]);
int ev_count = real_eigenvectors(coxeter, coxeter_fixedpoints, ws);
if(ev_count == 1)
@ -82,11 +142,17 @@ int computeLimitCurve(DrawingContext *ctx)
if(ev_count == 0)
goto error_out;
ctx->limit_curve_count = 0;
for(int i = 0; i < ctx->n_group_elements; i++) {
if(group[i].length % 2)
continue;
multiply_many(ws, fixedpoints, 3, ctx->cob, elements[i], coxeter_fixedpoints);
x = ctx->limit_curve[3*i ] = gsl_matrix_get(fixedpoints, 0, column)/gsl_matrix_get(fixedpoints, 2, column);
y = ctx->limit_curve[3*i+1] = gsl_matrix_get(fixedpoints, 1, column)/gsl_matrix_get(fixedpoints, 2, column);
x = ctx->limit_curve[3*ctx->limit_curve_count ] = gsl_matrix_get(fixedpoints, 0, column)/gsl_matrix_get(fixedpoints, 2, column);
y = ctx->limit_curve[3*ctx->limit_curve_count+1] = gsl_matrix_get(fixedpoints, 1, column)/gsl_matrix_get(fixedpoints, 2, column);
ctx->limit_curve_count++;
if((x - ctx->marking.x)*(x - ctx->marking.x) + (y - ctx->marking.y)*(y - ctx->marking.y) < 25e-10)
{
@ -95,19 +161,16 @@ int computeLimitCurve(DrawingContext *ctx)
fputc('a' + cur->letter, stdout); // bcbcbca, bacbcacab, bc bca cb
fputc('\n',stdout);
}
// bca abc acb = abc
}
qsort(ctx->limit_curve, ctx->n_group_elements, 3*sizeof(double), compareAngle);
qsort(ctx->limit_curve, ctx->limit_curve_count, 3*sizeof(double), compareAngle);
ctx->limit_curve_count = ctx->n_group_elements;
// ctx->limit_curve_count = ctx->n_group_elements;
success = 1;
error_out:
releaseTempMatrices(ctx->ws, 11+ctx->n_group_elements);
releaseTempMatrices(ctx->ws, 14+ctx->n_group_elements);
return success;
}

2
linalg.h
View File

@ -13,7 +13,7 @@
#define ERROR(condition, msg, ...) if(condition){fprintf(stderr, msg, ##__VA_ARGS__); exit(1);}
#define FCMP(x, y) gsl_fcmp(x, y, 1e-10)
#define MAX_TEMP_MATRICES 1200000
#define MAX_TEMP_MATRICES 600000
#define MAX_TEMP_VECTORS 100
typedef struct _workspace {

742
main.c
View File

@ -11,458 +11,474 @@
#include "linalg.h"
#define TOGGLE(a) do { (a) = !(a); } while(0)
#define SIGN(x) ((x) > 0 ? 1.0 : -1.0)
DrawingContext *screen_context;
// setup everything except cairo and dim, which will be provided by the graphics system
void setupContext(DrawingContext *ctx, int argc, char *argv[])
{
ctx->n_group_elements = NUM_GROUP_ELEMENTS;
ctx->p[0] = atoi(argv[1]);
ctx->p[1] = atoi(argv[2]);
ctx->p[2] = atoi(argv[3]);
ctx->k[0] = atoi(argv[4]);
ctx->k[1] = atoi(argv[5]);
ctx->k[2] = atoi(argv[6]);
if(argc > 7)
ctx->parameter = atof(argv[7]);
else
ctx->parameter = 1.0;
// ctx->parameter = 2.77;
// ctx->parameter = 0.1;
ctx->show_boxes = 0;
ctx->show_boxes2 = 0;
ctx->show_attractors = 0;
ctx->show_reflectors = 0;
ctx->show_rotated_reflectors = 0;
ctx->show_limit= 0;
ctx->show_dual_limit= 0;
ctx->show_text = 1;
ctx->mode = 0;
ctx->use_rotation_basis = 0;
ctx->limit_with_lines = 1;
ctx->use_repelling = 0;
ctx->show_marking = 1;
ctx->marking.x = -0.73679;
ctx->marking.y = -0.01873;
ctx->show_coxeter_orbit = 0;
ctx->n_group_elements = NUM_GROUP_ELEMENTS;
ctx->n_group_elements_combinatorial = NUM_GROUP_ELEMENTS_COMBINATORIAL;
ctx->p[0] = atoi(argv[1]);
ctx->p[1] = atoi(argv[2]);
ctx->p[2] = atoi(argv[3]);
ctx->k[0] = atoi(argv[4]);
ctx->k[1] = atoi(argv[5]);
ctx->k[2] = atoi(argv[6]);
if(argc > 7)
ctx->parameter = atof(argv[7]);
else
ctx->parameter = 1.0;
if(argc > 8)
ctx->parameter2 = atof(argv[8]);
else
ctx->parameter2 = 1.0;
if(argc > 12) {
ctx->movie_filename = argv[9];
ctx->movie_parameter_duration = atof(argv[10]);
ctx->movie_parameter2_duration = atof(argv[11]);
ctx->movie_n_frames = atoi(argv[12]);
} else {
ctx->movie_n_frames = 0;
}
// ctx->parameter = 2.77;
// ctx->parameter = 0.1;
ctx->show_boxes = 0;
ctx->show_boxes2 = 0;
ctx->show_attractors = 0;
ctx->show_reflectors = 0;
ctx->show_rotated_reflectors = 0;
ctx->show_limit = 0;
ctx->show_dual_limit = 0;
ctx->show_text = 1;
ctx->mode = 0;
ctx->use_rotation_basis = 1;
ctx->limit_with_lines = 0;
ctx->use_repelling = 0;
ctx->show_marking = 0;
ctx->marking.x = -0.73679;
ctx->marking.y = -0.01873;
ctx->show_coxeter_orbit = 0;
ctx->limit_curve = malloc(3*ctx->n_group_elements*sizeof(double));
ctx->limit_curve_count = -1;
ctx->limit_curve = malloc(3*ctx->n_group_elements*sizeof(double));
ctx->limit_curve_count = -1;
ctx->group = malloc(ctx->n_group_elements*sizeof(groupelement_t));
generate_triangle_group(ctx->group, ctx->n_group_elements, ctx->p[0], ctx->p[1], ctx->p[2]);
ctx->group = malloc(ctx->n_group_elements_combinatorial*sizeof(groupelement_t));
generate_triangle_group(ctx->group, ctx->n_group_elements_combinatorial, ctx->p[0], ctx->p[1], ctx->p[2]);
// the temporary stuff
ctx->cartan = gsl_matrix_alloc(3, 3);
ctx->cob = gsl_matrix_alloc(3, 3);
ctx->ws = workspace_alloc(3);
// the temporary stuff
ctx->cartan = gsl_matrix_alloc(3, 3);
ctx->cob = gsl_matrix_alloc(3, 3);
ctx->ws = workspace_alloc(3);
}
void destroyContext(DrawingContext *ctx)
{
free(ctx->limit_curve);
free(ctx->group);
free(ctx->limit_curve);
free(ctx->group);
gsl_matrix_free(ctx->cartan);
gsl_matrix_free(ctx->cob);
gsl_matrix_free(ctx->cartan);
gsl_matrix_free(ctx->cob);
workspace_free(ctx->ws);
workspace_free(ctx->ws);
}
void computeRotationMatrix(DrawingContext *ctx, gsl_matrix *result, const char *type)
void computeMatrix(DrawingContext *ctx, gsl_matrix *result, const char *type)
{
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix **gen = getTempMatrices(ctx->ws, 3);
gsl_matrix **gen = getTempMatrices(ctx->ws, 6);
// ERROR(strlen(type) != 2, "Invalid call of computeRotationMatrix()\n");
// ERROR(strlen(type) != 2, "Invalid call of computeRotationMatrix()\n");
initializeTriangleGenerators(gen, ctx->cartan);
gsl_matrix_set_identity(tmp);
for(int i = 0; i < strlen(type); i++)
multiply_right(tmp, gen[type[i]-'a'], ctx->ws);
initializeTriangleGeneratorsCurrent(gen, ctx);
gsl_matrix_set_identity(result);
for(int i = 0; i < strlen(type); i++) {
if(type[i] >= 'a' && type[i] <= 'c')
multiply_right(result, gen[type[i]-'a'], ctx->ws);
else if(type[i] >= 'A' && type[i] <= 'C')
multiply_right(result, gen[type[i]-'A'+3], ctx->ws);
}
rotation_frame(tmp, result, ctx->ws);
releaseTempMatrices(ctx->ws, 6);
}
releaseTempMatrices(ctx->ws, 4);
void computeRotationMatrixFrame(DrawingContext *ctx, gsl_matrix *result, const char *type)
{
gsl_matrix *tmp = getTempMatrix(ctx->ws);
computeMatrix(ctx, tmp, type);
rotation_frame(tmp, result, ctx->ws);
releaseTempMatrices(ctx->ws, 1);
}
void computeBoxTransform(DrawingContext *ctx, char *word1, char *word2, gsl_matrix *result)
{
vector_t p[2][3],i[2];
vector_t std[4] = {
{-1, -1, 1},
{-1, 1, 1},
{1, 1, 1},
{1, -1, 1}
};
vector_t p[2][3],i[2];
vector_t std[4] = {
{-1, -1, 1},
{-1, 1, 1},
{1, 1, 1},
{1, -1, 1}
};
gsl_vector **vertices = getTempVectors(ctx->ws, 4);
gsl_vector **std_vertices = getTempVectors(ctx->ws, 4);
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix *to_frame = getTempMatrix(ctx->ws);
gsl_matrix *to_std_frame = getTempMatrix(ctx->ws);
gsl_vector **vertices = getTempVectors(ctx->ws, 4);
gsl_vector **std_vertices = getTempVectors(ctx->ws, 4);
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix *to_frame = getTempMatrix(ctx->ws);
gsl_matrix *to_std_frame = getTempMatrix(ctx->ws);
fixedPoints(ctx, word1, p[0]);
fixedPoints(ctx, word2, p[1]);
fixedPoints(ctx, word1, p[0]);
fixedPoints(ctx, word2, p[1]);
// intersect attracting line with neutral line of the other element
for(int j = 0; j < 2; j++)
i[j] = cross(cross(p[j%2][0],p[j%2][1]),cross(p[(j+1)%2][0],p[(j+1)%2][2]));
// intersect attracting line with neutral line of the other element
for(int j = 0; j < 2; j++)
i[j] = cross(cross(p[j%2][0],p[j%2][1]),cross(p[(j+1)%2][0],p[(j+1)%2][2]));
// box consists of p[0][0], i[0], p[1][0], i[1]
// box consists of p[0][0], i[0], p[1][0], i[1]
for(int i = 0; i < 4; i++)
vectorToGsl(std[i], std_vertices[i]);
for(int i = 0; i < 4; i++)
vectorToGsl(std[i], std_vertices[i]);
vectorToGsl(p[0][0], vertices[0]);
vectorToGsl(i[0], vertices[1]);
vectorToGsl(p[1][0], vertices[2]);
vectorToGsl(i[1], vertices[3]);
vectorToGsl(p[0][0], vertices[0]);
vectorToGsl(i[0], vertices[1]);
vectorToGsl(p[1][0], vertices[2]);
vectorToGsl(i[1], vertices[3]);
projective_frame(std_vertices, to_std_frame, ctx->ws);
projective_frame(vertices, to_frame, ctx->ws);
invert(to_frame, tmp, ctx->ws);
multiply(to_std_frame, tmp, result);
projective_frame(std_vertices, to_std_frame, ctx->ws);
projective_frame(vertices, to_frame, ctx->ws);
invert(to_frame, tmp, ctx->ws);
multiply(to_std_frame, tmp, result);
/*
LOOP(i) {
LOOP(j) {
printf("%.4f ", gsl_matrix_get(result, i, j));
}
printf("\n");
}*/
/*
LOOP(i) {
LOOP(j) {
printf("%.4f ", gsl_matrix_get(result, i, j));
}
printf("\n");
}*/
releaseTempVectors(ctx->ws, 8);
releaseTempMatrices(ctx->ws, 3);
releaseTempVectors(ctx->ws, 8);
releaseTempMatrices(ctx->ws, 3);
}
void updateMatrices(DrawingContext *ctx)
{
double angle[3];
LOOP(i) angle[i] = M_PI*ctx->k[i]/ctx->p[i];
cartanMatrix(ctx->cartan, angle[0], angle[1], angle[2], ctx->parameter);
double angle[3];
LOOP(i) angle[i] = M_PI*ctx->k[i]/ctx->p[i];
cartanMatrix(ctx->cartan, angle[0], angle[1], angle[2], ctx->parameter);
gsl_matrix *tmp = getTempMatrix(ctx->ws);
gsl_matrix *tmp = getTempMatrix(ctx->ws);
int nmodes = 5;
if(ctx->use_rotation_basis % 5 == 0) {
gsl_matrix_set(tmp, 0, 0, 0.0);
gsl_matrix_set(tmp, 0, 1, sqrt(3.0)/2.0);
gsl_matrix_set(tmp, 0, 2, -sqrt(3.0)/2.0);
gsl_matrix_set(tmp, 1, 0, 1.0);
gsl_matrix_set(tmp, 1, 1, -0.5);
gsl_matrix_set(tmp, 1, 2, -0.5);
gsl_matrix_set(tmp, 2, 0, 1.0);
gsl_matrix_set(tmp, 2, 1, 1.0);
gsl_matrix_set(tmp, 2, 2, 1.0);
gsl_matrix_memcpy(ctx->cob, tmp);
} else if(ctx->use_rotation_basis % 5 == 1) {
gsl_matrix_memcpy(ctx->cob, ctx->cartan); // is this a good choice of basis for any reason?
} else if(ctx->use_rotation_basis % 5 == 2) {
computeRotationMatrix(ctx, tmp, "ba");
invert(tmp, ctx->cob, ctx->ws);
} else if(ctx->use_rotation_basis % 5 == 3) {
computeBoxTransform(ctx, "bca", "abc", ctx->cob);
// computeBoxTransform(ctx, "cab", "bca", ctx->cob);
// computeBoxTransform(ctx, "acb", "cba", ctx->cob);
} else {
cartanMatrix(tmp, M_PI/ctx->p[0], M_PI/ctx->p[1], M_PI/ctx->p[2], 1.0);
diagonalize_symmetric_form(tmp, ctx->cob, ctx->ws);
}
if(ctx->use_rotation_basis % nmodes == 0) {
gsl_matrix_set(tmp, 0, 0, 0.0);
gsl_matrix_set(tmp, 0, 1, sqrt(3.0)/2.0);
gsl_matrix_set(tmp, 0, 2, -sqrt(3.0)/2.0);
gsl_matrix_set(tmp, 1, 0, 1.0);
gsl_matrix_set(tmp, 1, 1, -0.5);
gsl_matrix_set(tmp, 1, 2, -0.5);
gsl_matrix_set(tmp, 2, 0, 1.0);
gsl_matrix_set(tmp, 2, 1, 1.0);
gsl_matrix_set(tmp, 2, 2, 1.0);
gsl_matrix_memcpy(ctx->cob, tmp);
} else if(ctx->use_rotation_basis % nmodes == 1) {
gsl_matrix_set(tmp, 0, 0, 1.0);
gsl_matrix_set(tmp, 0, 1, -1.0);
gsl_matrix_set(tmp, 0, 2, 0.0);
gsl_matrix_set(tmp, 1, 0, 1.0);
gsl_matrix_set(tmp, 1, 1, 1.0);
gsl_matrix_set(tmp, 1, 2, 0.0);
gsl_matrix_set(tmp, 2, 0, 0.0);
gsl_matrix_set(tmp, 2, 1, 0.0);
gsl_matrix_set(tmp, 2, 2, 1.0);
gsl_matrix_memcpy(ctx->cob, ctx->cartan); // is this a good choice of basis for any reason?
multiply_left(tmp, ctx->cob, ctx->ws);
} else if(ctx->use_rotation_basis % nmodes == 2) {
computeRotationMatrixFrame(ctx, tmp, "C");
invert(tmp, ctx->cob, ctx->ws);
} else if(ctx->use_rotation_basis % nmodes == 3) {
computeBoxTransform(ctx, "acb", "cba", ctx->cob);
// computeBoxTransform(ctx, "cab", "bca", ctx->cob);
// computeBoxTransform(ctx, "acb", "cba", ctx->cob);
} else {
cartanMatrix(tmp, M_PI/ctx->p[0], M_PI/ctx->p[1], M_PI/ctx->p[2], 1.0);
diagonalize_symmetric_form(tmp, ctx->cob, ctx->ws);
}
releaseTempMatrices(ctx->ws, 1);
releaseTempMatrices(ctx->ws, 1);
}
void output_info(DrawingContext *ctx)
{
vector_t p[4][3];
point_t pt;
vector_t p[4][3];
point_t pt;
fixedPoints(ctx, "abc", p[0]);
fixedPoints(ctx, "bca", p[1]);
fixedPoints(ctx, "cab", p[2]);
fixedPoints(ctx, "abc", p[0]);
fixedPoints(ctx, "bca", p[1]);
fixedPoints(ctx, "cab", p[2]);
pt = vectorToPoint(ctx, p[0][0]);
printf("(abc)-+ = (%f %f)\n", pt.x, pt.y);
pt = vectorToPoint(ctx, p[1][0]);
printf("(bca)-+ = (%f %f)\n", pt.x, pt.y);
pt = vectorToPoint(ctx, p[0][0]);
printf("(abc)-+ = (%f %f)\n", pt.x, pt.y);
pt = vectorToPoint(ctx, p[1][0]);
printf("(bca)-+ = (%f %f)\n", pt.x, pt.y);
}
void print(DrawingContext *screen)
{
DrawingContext file;
DimensionsInfo dim;
cairo_surface_t *surface;
DrawingContext file;
DimensionsInfo dim;
cairo_surface_t *surface;
char filename[100];
time_t t = time(NULL);
strftime(filename, sizeof(filename), "screenshot_%Y%m%d_%H%M%S.pdf", localtime(&t));
char filename[100];
time_t t = time(NULL);
strftime(filename, sizeof(filename), "screenshot_%Y%m%d_%H%M%S.pdf", localtime(&t));
memcpy(&file, screen, sizeof(file));
memcpy(&file, screen, sizeof(file));
dim.width = screen->dim->width;
dim.height = screen->dim->width / sqrt(2.0);
dim.matrix = screen->dim->matrix;
dim.matrix.y0 += ((double)dim.height - (double)screen->dim->height) / 2.0; // recenter vertically
updateDimensions(&dim);
file.dim = &dim;
dim.width = screen->dim->width;
dim.height = screen->dim->width / sqrt(2.0);
dim.matrix = screen->dim->matrix;
dim.matrix.y0 += ((double)dim.height - (double)screen->dim->height) / 2.0; // recenter vertically
updateDimensions(&dim);
file.dim = &dim;
surface = cairo_pdf_surface_create(filename, (double)dim.width, (double)dim.height);
surface = cairo_pdf_surface_create(filename, (double)dim.width, (double)dim.height);
file.cairo = cairo_create(surface);
draw(&file);
draw(&file);
cairo_destroy(file.cairo);
cairo_surface_destroy(surface);
cairo_destroy(file.cairo);
cairo_surface_destroy(surface);
printf("Wrote sceenshot to file: %s\n", filename);
printf("Wrote sceenshot to file: %s\n", filename);
}
int processEvent(GraphicsInfo *info, XEvent *ev)
{
int state;
unsigned long key;
char filename[100];
int state;
unsigned long key;
char filename[100];
// fprintf(stderr, "Event: %d\n", ev->type);
// fprintf(stderr, "Event: %d\n", ev->type);
switch(ev->type) {
case ButtonPress:
state = ev->xbutton.state & (ShiftMask | LockMask | ControlMask);
switch(ev->type) {
case ButtonPress:
state = ev->xbutton.state & (ShiftMask | LockMask | ControlMask);
if(ev->xbutton.button == 1 && state & ShiftMask) {
screen_context->marking.x = (double)ev->xbutton.x;
screen_context->marking.y = (double)ev->xbutton.y;
printf("mouse button pressed: %f, %f\n", screen_context->marking.x, screen_context->marking.y);
cairo_set_matrix(screen_context->cairo, &screen_context->dim->matrix);
cairo_device_to_user(screen_context->cairo, &screen_context->marking.x, &screen_context->marking.y);
printf("mouse button pressed transformed: %f, %f\n", screen_context->marking.x, screen_context->marking.y);
return STATUS_REDRAW;
}
break;
if(ev->xbutton.button == 1 && state & ShiftMask) {
screen_context->marking.x = (double)ev->xbutton.x;
screen_context->marking.y = (double)ev->xbutton.y;
printf("mouse button pressed: %f, %f\n", screen_context->marking.x, screen_context->marking.y);
cairo_set_matrix(screen_context->cairo, &screen_context->dim->matrix);
cairo_device_to_user(screen_context->cairo, &screen_context->marking.x, &screen_context->marking.y);
printf("mouse button pressed transformed: %f, %f\n", screen_context->marking.x, screen_context->marking.y);
return STATUS_REDRAW;
}
break;
case KeyPress:
state = ev->xkey.state & (ShiftMask | LockMask | ControlMask);
key = XkbKeycodeToKeysym(ev->xkey.display, ev->xkey.keycode, 0, !!(state & ShiftMask));
printf("Key pressed: %ld\n", key);
case KeyPress:
state = ev->xkey.state & (ShiftMask | LockMask | ControlMask);
key = XkbKeycodeToKeysym(ev->xkey.display, ev->xkey.keycode, 0, !!(state & ShiftMask));
printf("Key pressed: %ld\n", key);
switch(key) {
case XK_Down:
screen_context->parameter /= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Up:
screen_context->parameter *= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Left:
screen_context->parameter /= exp(0.00002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Right:
screen_context->parameter *= exp(0.00002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Page_Down:
screen_context->parameter /= exp(0.02);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Page_Up:
screen_context->parameter *= exp(0.02);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case ' ':
screen_context->parameter = 5.57959706;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Return:
// screen_context->parameter = 2.76375163;
screen_context->parameter = 5.29063366;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case 'm':
printf("matrix.xx = %f;\n", info->dim->matrix.xx);
printf("matrix.xy = %f;\n", info->dim->matrix.xy);
printf("matrix.x0 = %f;\n", info->dim->matrix.x0);
printf("matrix.yx = %f;\n", info->dim->matrix.yx);
printf("matrix.yy = %f;\n", info->dim->matrix.yy);
printf("matrix.y0 = %f;\n", info->dim->matrix.y0);
break;
case 'i':
output_info(screen_context);
break;
case 'b':
TOGGLE(screen_context->show_boxes);
break;
case 'B':
TOGGLE(screen_context->show_boxes2);
break;
case 'a':
TOGGLE(screen_context->show_attractors);
break;
case 'r':
TOGGLE(screen_context->show_reflectors);
break;
case 'x':
TOGGLE(screen_context->show_rotated_reflectors);
break;
case 'L':
TOGGLE(screen_context->limit_with_lines);
break;
case 'l':
TOGGLE(screen_context->show_limit);
break;
case 'd':
TOGGLE(screen_context->show_dual_limit);
break;
case 'R':
screen_context->use_rotation_basis++;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case 'p':
print(screen_context);
break;
case 'M':
/*
screen_context->limit_with_lines = 0;
double parameter_start = screen_context->parameter;
for(int i = 0; i <= 1300; i++) {
if(i < 400)
screen_context->parameter = exp(log(parameter_start)+0.002*i);
else if(i < 500)
screen_context->parameter = exp(log(parameter_start)+0.002*400);
else
screen_context->parameter = exp(log(parameter_start)+0.002*(900-i));
updateMatrices(screen_context);
computeLimitCurve(screen_context);
draw(screen_context);
sprintf(filename, "movie3/test%03d.png", i);
cairo_surface_write_to_png(info->buffer_surface, filename);
printf("Finished drawing %s\n", filename);
}
*/
/*
screen_context->limit_with_lines = 0;
double parameter_start = screen_context->parameter;
for(int i = 0; i <= 1300; i++) {
if(i < 400)
screen_context->parameter = exp(0.003*i);
else if(i < 500)
screen_context->parameter = exp(0.003*400);
else
screen_context->parameter = exp(0.003*(900-i));
updateMatrices(screen_context);
computeLimitCurve(screen_context);
draw(screen_context);
sprintf(filename, "movie5/test%03d.png", i);
cairo_surface_write_to_png(info->buffer_surface, filename);
printf("Finished drawing %s\n", filename);
}
*/
screen_context->limit_with_lines = 0;
for(int i = 0; i <= 1000; i++) {
screen_context->parameter = exp(log(3.0)*((double)i/500-1));
updateMatrices(screen_context);
computeLimitCurve(screen_context);
draw(screen_context);
sprintf(filename, "output/movie8/test%04d.png", i);
// cairo_surface_write_to_png(info->buffer_surface, filename);
printf("Finished drawing %s\n", filename);
}
case 'f':
TOGGLE(screen_context->use_repelling);
computeLimitCurve(screen_context);
break;
case 't':
TOGGLE(screen_context->show_text);
break;
case 'c':
TOGGLE(screen_context->show_coxeter_orbit);
break;
case '1': case '2': case '3': case '4': case '5': case '6': case '7': case '8': case '9': case '0':
screen_context->mode = key - '0';
break;
}
switch(key) {
case XK_Down:
if(ev->xkey.state & ShiftMask)
screen_context->parameter /= exp(0.00005);
else
screen_context->parameter /= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Up:
if(ev->xkey.state & ShiftMask)
screen_context->parameter *= exp(0.00005);
else
screen_context->parameter *= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Left:
if(ev->xkey.state & ShiftMask)
screen_context->parameter2 /= exp(0.00005);
else
screen_context->parameter2 /= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Right:
if(ev->xkey.state & ShiftMask)
screen_context->parameter2 *= exp(0.00005);
else
screen_context->parameter2 *= exp(0.002);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Page_Down:
screen_context->parameter /= exp(0.02);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Page_Up:
screen_context->parameter *= exp(0.02);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case ' ':
screen_context->parameter = 5.57959706;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case XK_Return:
// screen_context->parameter = 2.76375163;
screen_context->parameter = 5.29063366;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case 'm':
printf("matrix.xx = %f;\n", info->dim->matrix.xx);
printf("matrix.xy = %f;\n", info->dim->matrix.xy);
printf("matrix.x0 = %f;\n", info->dim->matrix.x0);
printf("matrix.yx = %f;\n", info->dim->matrix.yx);
printf("matrix.yy = %f;\n", info->dim->matrix.yy);
printf("matrix.y0 = %f;\n", info->dim->matrix.y0);
break;
case 'i':
output_info(screen_context);
break;
case 'b':
TOGGLE(screen_context->show_boxes);
break;
case 'B':
TOGGLE(screen_context->show_boxes2);
break;
case 'a':
TOGGLE(screen_context->show_attractors);
break;
case 'r':
TOGGLE(screen_context->show_reflectors);
break;
case 'x':
TOGGLE(screen_context->show_rotated_reflectors);
break;
case 'L':
TOGGLE(screen_context->limit_with_lines);
break;
case 'l':
TOGGLE(screen_context->show_limit);
break;
case 'd':
TOGGLE(screen_context->show_dual_limit);
break;
case 'R':
screen_context->use_rotation_basis++;
updateMatrices(screen_context);
computeLimitCurve(screen_context);
break;
case 'p':
print(screen_context);
break;
case 'M':
screen_context->limit_with_lines = 0;
double parameter_start = screen_context->parameter;
double parameter2_start = screen_context->parameter2;
for(int i = 0; i <= screen_context->movie_n_frames; i++) {
screen_context->parameter = SIGN(parameter_start)*exp(log(fabs(parameter_start)) +
i*screen_context->movie_parameter_duration/screen_context->movie_n_frames);
screen_context->parameter2 = SIGN(parameter2_start)*exp(log(fabs(parameter2_start)) +
i*screen_context->movie_parameter2_duration/screen_context->movie_n_frames);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
draw(screen_context);
sprintf(filename, "output/%s%03d.png", screen_context->movie_filename, i);
cairo_surface_write_to_png(info->buffer_surface, filename);
printf("Finished drawing %s\n", filename);
}
return STATUS_REDRAW;
case 'f':
TOGGLE(screen_context->use_repelling);
computeLimitCurve(screen_context);
break;
case 't':
TOGGLE(screen_context->show_text);
break;
case 'c':
TOGGLE(screen_context->show_coxeter_orbit);
break;
case '1': case '2': case '3': case '4': case '5': case '6': case '7': case '8': case '9': case '0':
screen_context->mode = key - '0';
break;
}
return STATUS_NOTHING;
return STATUS_REDRAW;
}
return STATUS_NOTHING;
}
int main(int argc, char *argv[])
{
GraphicsInfo *info;
GraphicsInfo *info;
screen_context = malloc(sizeof(DrawingContext));
setupContext(screen_context, argc, argv);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
screen_context = malloc(sizeof(DrawingContext));
setupContext(screen_context, argc, argv);
updateMatrices(screen_context);
computeLimitCurve(screen_context);
info = initCairo(0, KeyPressMask, 200, 200, "Triangle group");
if(!info)
return 1;
info = initCairo(0, KeyPressMask, 200, 200, "Triangle group");
if(!info)
return 1;
/*
info->dim->matrix.xx = 274.573171;
info->dim->matrix.xy = 0.000000;
info->dim->matrix.x0 = 583.073462;
info->dim->matrix.yx = 0.000000;
info->dim->matrix.yy = 274.573171;
info->dim->matrix.y0 = 777.225293;
*/
/*
info->dim->matrix.xx = 274.573171;
info->dim->matrix.xy = 0.000000;
info->dim->matrix.x0 = 583.073462;
info->dim->matrix.yx = 0.000000;
info->dim->matrix.yy = 274.573171;
info->dim->matrix.y0 = 777.225293;
*/
info->dim->matrix.xx = 274.573171;
info->dim->matrix.xy = 0.000000;
info->dim->matrix.x0 = 910.073462;
info->dim->matrix.yx = 0.000000;
info->dim->matrix.yy = 274.573171;
info->dim->matrix.y0 = 509.225293;
info->dim->matrix.xx = 274.573171;
info->dim->matrix.xy = 0.000000;
info->dim->matrix.x0 = 910.073462;
info->dim->matrix.yx = 0.000000;
info->dim->matrix.yy = 274.573171;
info->dim->matrix.y0 = 509.225293;
updateDimensions(info->dim);
updateDimensions(info->dim);
screen_context->dim = info->dim;
screen_context->cairo = info->buffer_context;
screen_context->dim = info->dim;
screen_context->cairo = info->buffer_context;
startTimer(info);
startTimer(info);
while(1) {
int result = checkEvents(info, processEvent, NULL);
if(result == STATUS_QUIT)
return 0;
else if(result == STATUS_REDRAW) {
struct timeval current_time;
double start_time, intermediate_time, end_time;
gettimeofday(&current_time, 0);
start_time = current_time.tv_sec + current_time.tv_usec*1e-6;
while(1) {
int result = checkEvents(info, processEvent, NULL);
if(result == STATUS_QUIT)
return 0;
else if(result == STATUS_REDRAW) {
struct timeval current_time;
double start_time, intermediate_time, end_time;
gettimeofday(&current_time, 0);
start_time = current_time.tv_sec + current_time.tv_usec*1e-6;
draw(screen_context);
draw(screen_context);
gettimeofday(&current_time, 0);
intermediate_time = current_time.tv_sec + current_time.tv_usec*1e-6;
gettimeofday(&current_time, 0);
intermediate_time = current_time.tv_sec + current_time.tv_usec*1e-6;
cairo_set_source_surface(info->front_context, info->buffer_surface, 0, 0);
cairo_paint(info->front_context);
cairo_set_source_surface(info->front_context, info->buffer_surface, 0, 0);
cairo_paint(info->front_context);
gettimeofday(&current_time, 0);
end_time = current_time.tv_sec + current_time.tv_usec*1e-6;
printf("drawing finished in %.2f milliseconds, of which %.2f milliseconds were buffer switching\n", (end_time - start_time) * 1000, (end_time - intermediate_time) * 1000);
}
waitUpdateTimer(info);
gettimeofday(&current_time, 0);
end_time = current_time.tv_sec + current_time.tv_usec*1e-6;
printf("drawing finished in %.2f milliseconds, of which %.2f milliseconds were buffer switching\n", (end_time - start_time) * 1000, (end_time - intermediate_time) * 1000);
}
waitUpdateTimer(info);
}
free(screen_context);
destroyCairo(info);
destroyContext(screen_context);
free(screen_context);
destroyCairo(info);
destroyContext(screen_context);
return 0;
return 0;
}

17
main.h
View File

@ -12,6 +12,11 @@
#define LOOP(i) for(int i = 0; i < 3; i++)
#define NUM_GROUP_ELEMENTS 10000
#define NUM_GROUP_ELEMENTS_COMBINATORIAL 100000
// (0,1) -> 2, (1,2) -> 0, (2,0) -> 1
// (1,0) -> 5, (2,1) -> 3, (0,2) -> 4
#define ROTATION_LETTER(x,y) (((y)-(x)+3)%3 == 1 ? ((y)+1)%3 : ((x)+1)%3+3)
typedef struct {
double x[3];
@ -30,7 +35,13 @@ typedef struct {
// a priori parameter
int p[3],k[3];
double parameter;
double parameter2;
char *movie_filename;
double movie_parameter_duration;
double movie_parameter2_duration;
int movie_n_frames;
int n_group_elements;
int n_group_elements_combinatorial;
int show_boxes;
int show_boxes2;
int show_attractors;
@ -66,7 +77,8 @@ typedef enum {
// implemented in limit_set.c
void cartanMatrix(gsl_matrix *cartan, double a1, double a2, double a3, double s);
void initializeTriangleGenerators(gsl_matrix **gen, gsl_matrix *cartan);
void initializeTriangleGenerators(gsl_matrix **gen, double a1, double a2, double a3, double s, double t, workspace_t *ws);
void initializeTriangleGeneratorsCurrent(gsl_matrix **gen, DrawingContext *ctx);
int computeLimitCurve(DrawingContext *ctx);
// implemented in draw.c
@ -96,7 +108,8 @@ void setupContext(DrawingContext *ctx, int argc, char *argv[]);
void destroyContext(DrawingContext *ctx);
void print(DrawingContext *screen);
int processEvent(GraphicsInfo *info, XEvent *ev);
void computeRotationMatrix(DrawingContext *ctx, gsl_matrix *result, const char *type);
void computeMatrix(DrawingContext *ctx, gsl_matrix *result, const char *type);
void computeRotationMatrixFrame(DrawingContext *ctx, gsl_matrix *result, const char *type);
void updateMatrices(DrawingContext *ctx);
static vector_t vectorFromGsl(gsl_vector *v)

1
triangle.h
View File

@ -14,6 +14,7 @@ typedef struct _groupelement {
struct _groupelement *parent;
struct _groupelement *inverse;
int letter;
int visited;
} groupelement_t;
int generate_triangle_group(groupelement_t *group, int size, int k1, int k2, int k3);