fix weston launch error
[platform/upstream/weston.git] / tests / matrix-test.c
1 /*
2  * Copyright © 2012 Collabora, Ltd.
3  *
4  * Permission is hereby granted, free of charge, to any person obtaining
5  * a copy of this software and associated documentation files (the
6  * "Software"), to deal in the Software without restriction, including
7  * without limitation the rights to use, copy, modify, merge, publish,
8  * distribute, sublicense, and/or sell copies of the Software, and to
9  * permit persons to whom the Software is furnished to do so, subject to
10  * the following conditions:
11  *
12  * The above copyright notice and this permission notice (including the
13  * next paragraph) shall be included in all copies or substantial
14  * portions of the Software.
15  *
16  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
17  * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
18  * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
19  * NONINFRINGEMENT.  IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS
20  * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
21  * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
22  * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
23  * SOFTWARE.
24  */
25
26 #include "config.h"
27
28 #include <stdio.h>
29 #include <stdlib.h>
30 #include <math.h>
31 #include <unistd.h>
32 #include <signal.h>
33 #include <time.h>
34
35 #include "shared/matrix.h"
36
37 struct inverse_matrix {
38         double LU[16];          /* column-major */
39         unsigned perm[4];       /* permutation */
40 };
41
42 static struct timespec begin_time;
43
44 static void
45 reset_timer(void)
46 {
47         clock_gettime(CLOCK_MONOTONIC, &begin_time);
48 }
49
50 static double
51 read_timer(void)
52 {
53         struct timespec t;
54
55         clock_gettime(CLOCK_MONOTONIC, &t);
56         return (double)(t.tv_sec - begin_time.tv_sec) +
57                1e-9 * (t.tv_nsec - begin_time.tv_nsec);
58 }
59
60 static double
61 det3x3(const float *c0, const float *c1, const float *c2)
62 {
63         return (double)
64                 c0[0] * c1[1] * c2[2] +
65                 c1[0] * c2[1] * c0[2] +
66                 c2[0] * c0[1] * c1[2] -
67                 c0[2] * c1[1] * c2[0] -
68                 c1[2] * c2[1] * c0[0] -
69                 c2[2] * c0[1] * c1[0];
70 }
71
72 static double
73 determinant(const struct weston_matrix *m)
74 {
75         double det = 0;
76 #if 1
77         /* develop on last row */
78         det -= m->d[3 + 0 * 4] * det3x3(&m->d[4], &m->d[8], &m->d[12]);
79         det += m->d[3 + 1 * 4] * det3x3(&m->d[0], &m->d[8], &m->d[12]);
80         det -= m->d[3 + 2 * 4] * det3x3(&m->d[0], &m->d[4], &m->d[12]);
81         det += m->d[3 + 3 * 4] * det3x3(&m->d[0], &m->d[4], &m->d[8]);
82 #else
83         /* develop on first row */
84         det += m->d[0 + 0 * 4] * det3x3(&m->d[5], &m->d[9], &m->d[13]);
85         det -= m->d[0 + 1 * 4] * det3x3(&m->d[1], &m->d[9], &m->d[13]);
86         det += m->d[0 + 2 * 4] * det3x3(&m->d[1], &m->d[5], &m->d[13]);
87         det -= m->d[0 + 3 * 4] * det3x3(&m->d[1], &m->d[5], &m->d[9]);
88 #endif
89         return det;
90 }
91
92 static void
93 print_permutation_matrix(const struct inverse_matrix *m)
94 {
95         const unsigned *p = m->perm;
96         const char *row[4] = {
97                 "1 0 0 0\n",
98                 "0 1 0 0\n",
99                 "0 0 1 0\n",
100                 "0 0 0 1\n"
101         };
102
103         printf("  P =\n%s%s%s%s", row[p[0]], row[p[1]], row[p[2]], row[p[3]]);
104 }
105
106 static void
107 print_LU_decomposition(const struct inverse_matrix *m)
108 {
109         unsigned r, c;
110
111         printf("                            L                      "
112                "                               U\n");
113         for (r = 0; r < 4; ++r) {
114                 double v;
115
116                 for (c = 0; c < 4; ++c) {
117                         if (c < r)
118                                 v = m->LU[r + c * 4];
119                         else if (c == r)
120                                 v = 1.0;
121                         else
122                                 v = 0.0;
123                         printf(" %12.6f", v);
124                 }
125
126                 printf(" | ");
127
128                 for (c = 0; c < 4; ++c) {
129                         if (c >= r)
130                                 v = m->LU[r + c * 4];
131                         else
132                                 v = 0.0;
133                         printf(" %12.6f", v);
134                 }
135                 printf("\n");
136         }
137 }
138
139 static void
140 print_inverse_data_matrix(const struct inverse_matrix *m)
141 {
142         unsigned r, c;
143
144         for (r = 0; r < 4; ++r) {
145                 for (c = 0; c < 4; ++c)
146                         printf(" %12.6f", m->LU[r + c * 4]);
147                 printf("\n");
148         }
149
150         printf("permutation: ");
151         for (r = 0; r < 4; ++r)
152                 printf(" %u", m->perm[r]);
153         printf("\n");
154 }
155
156 static void
157 print_matrix(const struct weston_matrix *m)
158 {
159         unsigned r, c;
160
161         for (r = 0; r < 4; ++r) {
162                 for (c = 0; c < 4; ++c)
163                         printf(" %14.6e", m->d[r + c * 4]);
164                 printf("\n");
165         }
166 }
167
168 static double
169 frand(void)
170 {
171         double r = random();
172         return r / (double)(RAND_MAX / 2) - 1.0f;
173 }
174
175 static void
176 randomize_matrix(struct weston_matrix *m)
177 {
178         unsigned i;
179         for (i = 0; i < 16; ++i)
180 #if 1
181                 m->d[i] = frand() * exp(10.0 * frand());
182 #else
183                 m->d[i] = frand();
184 #endif
185 }
186
187 /* Take a matrix, compute inverse, multiply together
188  * and subtract the identity matrix to get the error matrix.
189  * Return the largest absolute value from the error matrix.
190  */
191 static double
192 test_inverse(struct weston_matrix *m)
193 {
194         unsigned i;
195         struct inverse_matrix q;
196         double errsup = 0.0;
197
198         if (matrix_invert(q.LU, q.perm, m) != 0)
199                 return INFINITY;
200
201         for (i = 0; i < 4; ++i)
202                 inverse_transform(q.LU, q.perm, &m->d[i * 4]);
203
204         m->d[0] -= 1.0f;
205         m->d[5] -= 1.0f;
206         m->d[10] -= 1.0f;
207         m->d[15] -= 1.0f;
208
209         for (i = 0; i < 16; ++i) {
210                 double err = fabs(m->d[i]);
211                 if (err > errsup)
212                         errsup = err;
213         }
214
215         return errsup;
216 }
217
218 enum {
219         TEST_OK,
220         TEST_NOT_INVERTIBLE_OK,
221         TEST_FAIL,
222         TEST_COUNT
223 };
224
225 static int
226 test(void)
227 {
228         struct weston_matrix m;
229         double det, errsup;
230
231         randomize_matrix(&m);
232         det = determinant(&m);
233
234         errsup = test_inverse(&m);
235         if (errsup < 1e-6)
236                 return TEST_OK;
237
238         if (fabs(det) < 1e-5 && isinf(errsup))
239                 return TEST_NOT_INVERTIBLE_OK;
240
241         printf("test fail, det: %g, error sup: %g\n", det, errsup);
242
243         return TEST_FAIL;
244 }
245
246 static int running;
247 static void
248 stopme(int n)
249 {
250         running = 0;
251 }
252
253 static void
254 test_loop_precision(void)
255 {
256         int counts[TEST_COUNT] = { 0 };
257
258         printf("\nRunning a test loop for 10 seconds...\n");
259         running = 1;
260         alarm(10);
261         while (running) {
262                 counts[test()]++;
263         }
264
265         printf("tests: %d ok, %d not invertible but ok, %d failed.\n"
266                "Total: %d iterations.\n",
267                counts[TEST_OK], counts[TEST_NOT_INVERTIBLE_OK],
268                counts[TEST_FAIL],
269                counts[TEST_OK] + counts[TEST_NOT_INVERTIBLE_OK] +
270                counts[TEST_FAIL]);
271 }
272
273 static void __attribute__((noinline))
274 test_loop_speed_matrixvector(void)
275 {
276         struct weston_matrix m;
277         struct weston_vector v = { { 0.5, 0.5, 0.5, 1.0 } };
278         unsigned long count = 0;
279         double t;
280
281         printf("\nRunning 3 s test on weston_matrix_transform()...\n");
282
283         weston_matrix_init(&m);
284
285         running = 1;
286         alarm(3);
287         reset_timer();
288         while (running) {
289                 weston_matrix_transform(&m, &v);
290                 count++;
291         }
292         t = read_timer();
293
294         printf("%lu iterations in %f seconds, avg. %.1f ns/iter.\n",
295                count, t, 1e9 * t / count);
296 }
297
298 static void __attribute__((noinline))
299 test_loop_speed_inversetransform(void)
300 {
301         struct weston_matrix m;
302         struct inverse_matrix inv;
303         struct weston_vector v = { { 0.5, 0.5, 0.5, 1.0 } };
304         unsigned long count = 0;
305         double t;
306
307         printf("\nRunning 3 s test on inverse_transform()...\n");
308
309         weston_matrix_init(&m);
310         matrix_invert(inv.LU, inv.perm, &m);
311
312         running = 1;
313         alarm(3);
314         reset_timer();
315         while (running) {
316                 inverse_transform(inv.LU, inv.perm, v.f);
317                 count++;
318         }
319         t = read_timer();
320
321         printf("%lu iterations in %f seconds, avg. %.1f ns/iter.\n",
322                count, t, 1e9 * t / count);
323 }
324
325 static void __attribute__((noinline))
326 test_loop_speed_invert(void)
327 {
328         struct weston_matrix m;
329         struct inverse_matrix inv;
330         unsigned long count = 0;
331         double t;
332
333         printf("\nRunning 3 s test on matrix_invert()...\n");
334
335         weston_matrix_init(&m);
336
337         running = 1;
338         alarm(3);
339         reset_timer();
340         while (running) {
341                 matrix_invert(inv.LU, inv.perm, &m);
342                 count++;
343         }
344         t = read_timer();
345
346         printf("%lu iterations in %f seconds, avg. %.1f ns/iter.\n",
347                count, t, 1e9 * t / count);
348 }
349
350 static void __attribute__((noinline))
351 test_loop_speed_invert_explicit(void)
352 {
353         struct weston_matrix m;
354         unsigned long count = 0;
355         double t;
356
357         printf("\nRunning 3 s test on weston_matrix_invert()...\n");
358
359         weston_matrix_init(&m);
360
361         running = 1;
362         alarm(3);
363         reset_timer();
364         while (running) {
365                 weston_matrix_invert(&m, &m);
366                 count++;
367         }
368         t = read_timer();
369
370         printf("%lu iterations in %f seconds, avg. %.1f ns/iter.\n",
371                count, t, 1e9 * t / count);
372 }
373
374 int main(void)
375 {
376         struct sigaction ding;
377         struct weston_matrix M;
378         struct inverse_matrix Q;
379         int ret;
380         double errsup;
381         double det;
382
383         ding.sa_handler = stopme;
384         sigemptyset(&ding.sa_mask);
385         ding.sa_flags = 0;
386         sigaction(SIGALRM, &ding, NULL);
387
388         srandom(13);
389
390         M.d[0] = 3.0;   M.d[4] = 17.0;  M.d[8] = 10.0;  M.d[12] = 0.0;
391         M.d[1] = 2.0;   M.d[5] = 4.0;   M.d[9] = -2.0;  M.d[13] = 0.0;
392         M.d[2] = 6.0;   M.d[6] = 18.0;  M.d[10] = -12;  M.d[14] = 0.0;
393         M.d[3] = 0.0;   M.d[7] = 0.0;   M.d[11] = 0.0;  M.d[15] = 1.0;
394
395         ret = matrix_invert(Q.LU, Q.perm, &M);
396         printf("ret = %d\n", ret);
397         printf("det = %g\n\n", determinant(&M));
398
399         if (ret != 0)
400                 return 1;
401
402         print_inverse_data_matrix(&Q);
403         printf("P * A = L * U\n");
404         print_permutation_matrix(&Q);
405         print_LU_decomposition(&Q);
406
407
408         printf("a random matrix:\n");
409         randomize_matrix(&M);
410         det = determinant(&M);
411         print_matrix(&M);
412         errsup = test_inverse(&M);
413         printf("\nThe matrix multiplied by its inverse, error:\n");
414         print_matrix(&M);
415         printf("max abs error: %g, original determinant %g\n", errsup, det);
416
417         test_loop_precision();
418         test_loop_speed_matrixvector();
419         test_loop_speed_inversetransform();
420         test_loop_speed_invert();
421         test_loop_speed_invert_explicit();
422
423         return 0;
424 }