1 /*
2 * Copyright (c) 2004, Bull S.A.. All rights reserved.
3 * Created by: Sebastien Decugis
4
5 * This program is free software; you can redistribute it and/or modify it
6 * under the terms of version 2 of the GNU General Public License as
7 * published by the Free Software Foundation.
8 *
9 * This program is distributed in the hope that it would be useful, but
10 * WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
12 *
13 * You should have received a copy of the GNU General Public License along
14 * with this program; if not, write the Free Software Foundation, Inc.,
15 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
16 *
17
18 * This file is a scalability test for the pthread_cond_timedwait function.
19 *
20 * It aims to measure the time between end of timeout and actual wakeup
21 * with different factors.
22
23 * The steps are:
24 * -*- Number of threads waiting on the conditionnal variable
25 * -> for an increaing number of threads,
26 * -> create the other threads which will do a pthread_cond_wait on the same cond/mutex
27 * -> When the threads are waiting, create one thread which will measure the time
28 * -> once the timeout has expired and the measure is done, broadcast the condition.
29 * -> do each measure 10 times (with different attributes i.e.)
30 *
31 * -*- other possible influencial parameters
32 * -> To be defined.
33 */
34
35 /* We are testing conformance to IEEE Std 1003.1, 2003 Edition */
36 #define _POSIX_C_SOURCE 200112L
37
38 #ifndef WITHOUT_XOPEN
39 #define _XOPEN_SOURCE 600
40 #endif
41
42 /********************************************************************************************/
43 /****************************** standard includes *****************************************/
44 /********************************************************************************************/
45 #include <pthread.h>
46 #include <stdarg.h>
47 #include <stdio.h>
48 #include <stdlib.h>
49 #include <unistd.h>
50
51 #include <time.h>
52 #include <errno.h>
53 #include <math.h>
54
55 /********************************************************************************************/
56 /****************************** Test framework *****************************************/
57 /********************************************************************************************/
58 #include "testfrmw.h"
59 #include "testfrmw.c"
60 /* This header is responsible for defining the following macros:
61 * UNRESOLVED(ret, descr);
62 * where descr is a description of the error and ret is an int (error code for example)
63 * FAILED(descr);
64 * where descr is a short text saying why the test has failed.
65 * PASSED();
66 * No parameter.
67 *
68 * Both three macros shall terminate the calling process.
69 * The testcase shall not terminate in any other maneer.
70 *
71 * The other file defines the functions
72 * void output_init()
73 * void output(char * string, ...)
74 *
75 * Those may be used to output information.
76 */
77
78 /********************************************************************************************/
79 /********************************** Configuration ******************************************/
80 /********************************************************************************************/
81 #ifndef SCALABILITY_FACTOR
82 #define SCALABILITY_FACTOR 1
83 #endif
84 #ifndef VERBOSE
85 #define VERBOSE 1
86 #endif
87
88 #ifndef WITHOUT_ALTCLK
89 #define USE_ALTCLK /* make tests with MONOTONIC CLOCK if supported */
90 #endif
91
92 #define MES_TIMEOUT (1000000) /* ns, offset for the pthread_cond_timedwait call */
93
94 #ifdef PLOT_OUTPUT
95 #undef VERBOSE
96 #define VERBOSE 0
97 #endif
98
99 // #define USE_CANCEL /* Will cancel the threads instead of signaling the cond */#define
100
101 /********************************************************************************************/
102 /*********************************** Test case *****************************************/
103 /********************************************************************************************/
104
105 long altclk_ok, pshared_ok;
106
107 typedef struct {
108 pthread_cond_t *cnd;
109 pthread_mutex_t *mtx;
110 int *predicate;
111 int *tnum;
112 } test_t;
113
114 struct {
115 int mutex_type;
116 int pshared;
117 clockid_t cid;
118 char *desc;
119 } test_scenar[] = {
120 {
121 PTHREAD_MUTEX_DEFAULT, PTHREAD_PROCESS_PRIVATE, CLOCK_REALTIME,
122 "Default"}
123 , {
124 PTHREAD_MUTEX_DEFAULT, PTHREAD_PROCESS_SHARED, CLOCK_REALTIME,
125 "PShared"}
126 #ifndef WITHOUT_XOPEN
127 , {
128 PTHREAD_MUTEX_NORMAL, PTHREAD_PROCESS_PRIVATE, CLOCK_REALTIME,
129 "Normal"}
130 , {
131 PTHREAD_MUTEX_NORMAL, PTHREAD_PROCESS_SHARED, CLOCK_REALTIME,
132 "Normal+PShared"}
133 , {
134 PTHREAD_MUTEX_ERRORCHECK, PTHREAD_PROCESS_PRIVATE,
135 CLOCK_REALTIME, "Errorcheck"}
136 , {
137 PTHREAD_MUTEX_ERRORCHECK, PTHREAD_PROCESS_SHARED,
138 CLOCK_REALTIME, "Errorcheck+PShared"}
139 , {
140 PTHREAD_MUTEX_RECURSIVE, PTHREAD_PROCESS_PRIVATE,
141 CLOCK_REALTIME, "Recursive"}
142 , {
143 PTHREAD_MUTEX_RECURSIVE, PTHREAD_PROCESS_SHARED, CLOCK_REALTIME,
144 "Recursive+PShared"}
145 #endif
146 #ifdef USE_ALTCLK
147 , {
148 PTHREAD_MUTEX_DEFAULT, PTHREAD_PROCESS_PRIVATE, CLOCK_MONOTONIC,
149 "Monotonic"}
150 , {
151 PTHREAD_MUTEX_DEFAULT, PTHREAD_PROCESS_SHARED, CLOCK_MONOTONIC,
152 "PShared+Monotonic"}
153 #ifndef WITHOUT_XOPEN
154 , {
155 PTHREAD_MUTEX_NORMAL, PTHREAD_PROCESS_PRIVATE, CLOCK_MONOTONIC,
156 "Normal+Monotonic"}
157 , {
158 PTHREAD_MUTEX_NORMAL, PTHREAD_PROCESS_SHARED, CLOCK_MONOTONIC,
159 "Normal+PShared+Monotonic"}
160 , {
161 PTHREAD_MUTEX_ERRORCHECK, PTHREAD_PROCESS_PRIVATE,
162 CLOCK_MONOTONIC, "Errorcheck+Monotonic"}
163 , {
164 PTHREAD_MUTEX_ERRORCHECK, PTHREAD_PROCESS_SHARED,
165 CLOCK_MONOTONIC, "Errorcheck+PShared+Monotonic"}
166 , {
167 PTHREAD_MUTEX_RECURSIVE, PTHREAD_PROCESS_PRIVATE,
168 CLOCK_MONOTONIC, "Recursive+Monotonic"}
169 , {
170 PTHREAD_MUTEX_RECURSIVE, PTHREAD_PROCESS_SHARED,
171 CLOCK_MONOTONIC, "Recursive+PShared+Monotonic"}
172 #endif
173 #endif
174 };
175
176 #define NSCENAR (sizeof(test_scenar) / sizeof(test_scenar[0]))
177
178 pthread_attr_t ta;
179
180 /* The next structure is used to save the tests measures */
181 typedef struct __mes_t {
182 int nthreads;
183 long _data[NSCENAR]; /* As we store µsec values, a long type should be amply enough. */
184 struct __mes_t *next;
185 } mes_t;
186
187 /**** do_measure
188 * This function will do a timedwait on the cond cnd after locking mtx.
189 * Once the timedwait times out, it will read the clock cid then
190 * compute the difference and put it into ts.
191 * This function must be called once test is ready, as the timeout will be very short. */
do_measure(pthread_mutex_t * mtx,pthread_cond_t * cnd,clockid_t cid,struct timespec * ts)192 void do_measure(pthread_mutex_t * mtx, pthread_cond_t * cnd, clockid_t cid,
193 struct timespec *ts)
194 {
195 int ret, rc;
196
197 struct timespec ts_cnd, ts_clk;
198
199 /* lock the mutex */
200 ret = pthread_mutex_lock(mtx);
201 if (ret != 0) {
202 UNRESOLVED(ret, "Unable to lock mutex");
203 }
204
205 /* Prepare the timeout parameter */
206 ret = clock_gettime(cid, &ts_cnd);
207 if (ret != 0) {
208 UNRESOLVED(errno, "Unable to read clock");
209 }
210
211 ts_cnd.tv_nsec += MES_TIMEOUT;
212 while (ts_cnd.tv_nsec >= 1000000000) {
213 ts_cnd.tv_nsec -= 1000000000;
214 ts_cnd.tv_sec++;
215 }
216
217 /* Do the timedwait */
218 do {
219 rc = pthread_cond_timedwait(cnd, mtx, &ts_cnd);
220 /* Re-read the clock as soon as possible */
221 ret = clock_gettime(cid, &ts_clk);
222 if (ret != 0) {
223 UNRESOLVED(errno, "Unable to read clock");
224 }
225 }
226 while (rc == 0);
227 if (rc != ETIMEDOUT) {
228 UNRESOLVED(rc,
229 "Timedwait returned an unexpected error (expected ETIMEDOUT)");
230 }
231
232 /* Compute the difference time */
233 ts->tv_sec = ts_clk.tv_sec - ts_cnd.tv_sec;
234 ts->tv_nsec = ts_clk.tv_nsec - ts_cnd.tv_nsec;
235 if (ts->tv_nsec < 0) {
236 ts->tv_nsec += 1000000000;
237 ts->tv_sec -= 1;
238 }
239
240 if (ts->tv_sec < 0) {
241 FAILED
242 ("The function returned from wait with a timeout before the time has passed\n");
243 }
244
245 /* unlock the mutex */
246 ret = pthread_mutex_unlock(mtx);
247 if (ret != 0) {
248 UNRESOLVED(ret, "Unable to unlock mutex");
249 }
250
251 return;
252 }
253
waiter(void * arg)254 void *waiter(void *arg)
255 {
256 test_t *dt = (test_t *) arg;
257
258 int ret;
259
260 ret = pthread_mutex_lock(dt->mtx);
261 if (ret != 0) {
262 UNRESOLVED(ret, "Mutex lock failed in waiter");
263 }
264 #ifdef USE_CANCEL
265 pthread_cleanup_push((void *)pthread_mutex_unlock, (void *)(dt->mtx));
266 #endif
267
268 /* This thread is ready to wait */
269 *(dt->tnum) += 1;
270
271 do {
272 ret = pthread_cond_wait(dt->cnd, dt->mtx);
273 }
274 while ((ret == 0) && (*(dt->predicate) == 0));
275 if (ret != 0) {
276 UNRESOLVED(ret, "pthread_cond_wait failed in waiter");
277 }
278 #ifdef USE_CANCEL
279 pthread_cleanup_pop(0); /* We could put 1 and avoid the next line, but we would miss the return code */
280 #endif
281
282 ret = pthread_mutex_unlock(dt->mtx);
283 if (ret != 0) {
284 UNRESOLVED(ret, "Mutex unlock failed in waiter");
285 }
286
287 return NULL;
288 }
289
290 /**** do_threads_test
291 * This function is responsible for all the testing with a given # of threads
292 * nthreads is the amount of threads to create.
293 * the return value is:
294 * < 0 if function was not able to create enough threads.
295 * cumulated # of nanoseconds otherwise.
296 */
do_threads_test(int nthreads,mes_t * measure)297 long do_threads_test(int nthreads, mes_t * measure)
298 {
299 int ret;
300
301 int scal, i, j;
302
303 pthread_t *th;
304
305 int s;
306 pthread_mutexattr_t ma;
307 pthread_condattr_t ca;
308
309 pthread_cond_t cnd;
310 pthread_mutex_t mtx;
311 int predicate;
312 int tnum;
313
314 test_t td;
315
316 struct timespec ts, ts_cumul;
317
318 td.mtx = &mtx;
319 td.cnd = &cnd;
320 td.predicate = &predicate;
321 td.tnum = &tnum;
322
323 /* Allocate space for the threads structures */
324 th = (pthread_t *) calloc(nthreads, sizeof(pthread_t));
325 if (th == NULL) {
326 UNRESOLVED(errno, "Not enough memory for thread storage");
327 }
328 #ifdef PLOT_OUTPUT
329 output("%d", nthreads);
330 #endif
331 /* For each test scenario (mutex and cond attributes) */
332 for (s = 0; s < NSCENAR; s++) {
333 /* Initialize the attributes */
334 ret = pthread_mutexattr_init(&ma);
335 if (ret != 0) {
336 UNRESOLVED(ret,
337 "Unable to initialize mutex attribute object");
338 }
339
340 ret = pthread_condattr_init(&ca);
341 if (ret != 0) {
342 UNRESOLVED(ret,
343 "Unable to initialize cond attribute object");
344 }
345
346 /* Set the attributes according to the scenar and the impl capabilities */
347 ret = pthread_mutexattr_settype(&ma, test_scenar[s].mutex_type);
348 if (ret != 0) {
349 UNRESOLVED(ret, "Unable to set mutex type");
350 }
351
352 if (pshared_ok > 0) {
353 ret =
354 pthread_mutexattr_setpshared(&ma,
355 test_scenar[s].
356 pshared);
357 if (ret != 0) {
358 UNRESOLVED(ret,
359 "Unable to set mutex process-shared");
360 }
361
362 ret =
363 pthread_condattr_setpshared(&ca,
364 test_scenar[s].pshared);
365 if (ret != 0) {
366 UNRESOLVED(ret,
367 "Unable to set cond process-shared");
368 }
369 }
370 #ifdef USE_ALTCLK
371 if (altclk_ok > 0) {
372 ret =
373 pthread_condattr_setclock(&ca, test_scenar[s].cid);
374 if (ret != 0) {
375 UNRESOLVED(ret, "Unable to set clock for cond");
376 }
377 }
378 #endif
379
380 ts_cumul.tv_sec = 0;
381 ts_cumul.tv_nsec = 0;
382
383 #if VERBOSE > 1
384 output("Starting case %s\n", test_scenar[s].desc);
385 #endif
386
387 for (scal = 0; scal < 5 * SCALABILITY_FACTOR; scal++) {
388 /* Initialize the mutex, the cond, and other data */
389 ret = pthread_mutex_init(&mtx, &ma);
390 if (ret != 0) {
391 UNRESOLVED(ret, "Mutex init failed");
392 }
393
394 ret = pthread_cond_init(&cnd, &ca);
395 if (ret != 0) {
396 UNRESOLVED(ret, "Cond init failed");
397 }
398
399 predicate = 0;
400 tnum = 0;
401
402 /* Create the waiter threads */
403 for (i = 0; i < nthreads; i++) {
404 ret =
405 pthread_create(&th[i], &ta, waiter,
406 (void *)&td);
407 if (ret != 0) { /* We reached the limits */
408 #if VERBOSE > 1
409 output
410 ("Limit reached with %i threads\n",
411 i);
412 #endif
413 #ifdef USE_CANCEL
414 for (j = i - 1; j >= 0; j--) {
415 ret = pthread_cancel(th[j]);
416 if (ret != 0) {
417 UNRESOLVED(ret,
418 "Unable to cancel a thread");
419 }
420 }
421 #else
422 predicate = 1;
423 ret = pthread_cond_broadcast(&cnd);
424 if (ret != 0) {
425 UNRESOLVED(ret,
426 "Unable to broadcast the condition");
427 }
428 #endif
429 for (j = i - 1; j >= 0; j--) {
430 ret = pthread_join(th[j], NULL);
431 if (ret != 0) {
432 UNRESOLVED(ret,
433 "Unable to join a canceled thread");
434 }
435 }
436 free(th);
437 return -1;
438 }
439 }
440 /* All waiter threads are created now */
441 #if VERBOSE > 5
442 output("%i waiter threads created successfully\n", i);
443 #endif
444
445 ret = pthread_mutex_lock(&mtx);
446 if (ret != 0) {
447 UNRESOLVED(ret, "Unable to lock mutex");
448 }
449
450 /* Wait for all threads be waiting */
451 while (tnum < nthreads) {
452 ret = pthread_mutex_unlock(&mtx);
453 if (ret != 0) {
454 UNRESOLVED(ret, "Mutex unlock failed");
455 }
456
457 sched_yield();
458
459 ret = pthread_mutex_lock(&mtx);
460 if (ret != 0) {
461 UNRESOLVED(ret, "Mutex lock failed");
462 }
463 }
464
465 /* All threads are now waiting - we do the measure */
466
467 ret = pthread_mutex_unlock(&mtx);
468 if (ret != 0) {
469 UNRESOLVED(ret, "Mutex unlock failed");
470 }
471 #if VERBOSE > 5
472 output("%i waiter threads are waiting; start measure\n",
473 tnum);
474 #endif
475
476 do_measure(&mtx, &cnd, test_scenar[s].cid, &ts);
477
478 #if VERBOSE > 5
479 output("Measure for %s returned %d.%09d\n",
480 test_scenar[s].desc, ts.tv_sec, ts.tv_nsec);
481 #endif
482
483 ts_cumul.tv_sec += ts.tv_sec;
484 ts_cumul.tv_nsec += ts.tv_nsec;
485 if (ts_cumul.tv_nsec >= 1000000000) {
486 ts_cumul.tv_nsec -= 1000000000;
487 ts_cumul.tv_sec += 1;
488 }
489
490 /* We can release the threads */
491 predicate = 1;
492 ret = pthread_cond_broadcast(&cnd);
493 if (ret != 0) {
494 UNRESOLVED(ret,
495 "Unable to broadcast the condition");
496 }
497 #if VERBOSE > 5
498 output("Joining the waiters...\n");
499 #endif
500
501 /* We will join every threads */
502 for (i = 0; i < nthreads; i++) {
503 ret = pthread_join(th[i], NULL);
504 if (ret != 0) {
505 UNRESOLVED(ret,
506 "Unable to join a thread");
507 }
508
509 }
510
511 /* Destroy everything */
512 ret = pthread_mutex_destroy(&mtx);
513 if (ret != 0) {
514 UNRESOLVED(ret, "Unable to destroy mutex");
515 }
516
517 ret = pthread_cond_destroy(&cnd);
518 if (ret != 0) {
519 UNRESOLVED(ret, "Unable to destroy cond");
520 }
521 }
522
523 #ifdef PLOT_OUTPUT
524 output(" %d.%09d", ts_cumul.tv_sec, ts_cumul.tv_nsec);
525 #endif
526
527 measure->_data[s] = ts_cumul.tv_sec * 1000000 + (ts_cumul.tv_nsec / 1000); /* We reduce precision */
528
529 /* Destroy the mutex attributes */
530 ret = pthread_mutexattr_destroy(&ma);
531 if (ret != 0) {
532 UNRESOLVED(ret, "Unable to destroy mutex attribute");
533 }
534
535 ret = pthread_condattr_destroy(&ca);
536 if (ret != 0) {
537 UNRESOLVED(ret, "Unable to destroy cond attribute");
538 }
539 }
540
541 /* Free the memory */
542 free(th);
543
544 #if VERBOSE > 2
545 output("%5d threads; %d.%09d s (%i loops)\n", nthreads, ts_cumul.tv_sec,
546 ts_cumul.tv_nsec, scal);
547 #endif
548
549 #ifdef PLOT_OUTPUT
550 output("\n");
551 #endif
552
553 return ts_cumul.tv_sec * 1000000000 + ts_cumul.tv_nsec;
554 }
555
556 /* Forward declaration */
557 int parse_measure(mes_t * measures);
558
559 /* Main
560 */
main(int argc,char * argv[])561 int main(int argc, char *argv[])
562 {
563 int ret, nth;
564 long dur;
565
566 /* Initialize the measure list */
567 mes_t sentinel;
568 mes_t *m_cur, *m_tmp;
569 m_cur = &sentinel;
570 m_cur->next = NULL;
571
572 /* Initialize the output */
573 output_init();
574
575 /* Test machine capabilities */
576 /* -> clockid_t; pshared; ... */
577 altclk_ok = sysconf(_SC_CLOCK_SELECTION);
578 if (altclk_ok > 0)
579 altclk_ok = sysconf(_SC_MONOTONIC_CLOCK);
580 #ifndef USE_ALTCLK
581 if (altclk_ok > 0)
582 output
583 ("Implementation supports the MONOTONIC CLOCK but option is disabled in test.\n");
584 #endif
585
586 pshared_ok = sysconf(_SC_THREAD_PROCESS_SHARED);
587
588 #if VERBOSE > 0
589 output("Test starting\n");
590 output(" Process-shared primitive %s be tested\n",
591 (pshared_ok > 0) ? "will" : "won't");
592 output(" Alternative clock for cond %s be tested\n",
593 (altclk_ok > 0) ? "will" : "won't");
594 #endif
595
596 /* Prepare thread attribute */
597 ret = pthread_attr_init(&ta);
598 if (ret != 0) {
599 UNRESOLVED(ret, "Unable to initialize thread attributes");
600 }
601
602 ret = pthread_attr_setstacksize(&ta, sysconf(_SC_THREAD_STACK_MIN));
603 if (ret != 0) {
604 UNRESOLVED(ret, "Unable to set stack size to minimum value");
605 }
606 #ifdef PLOT_OUTPUT
607 output("# COLUMNS %d #threads", NSCENAR + 1);
608 for (nth = 0; nth < NSCENAR; nth++)
609 output(" %s", test_scenar[nth].desc);
610 output("\n");
611 #endif
612
613 /* Do the testing */
614 nth = 0;
615 do {
616 nth += 100 * SCALABILITY_FACTOR;
617
618 /* Create a new measure item */
619 m_tmp = malloc(sizeof(mes_t));
620 if (m_tmp == NULL) {
621 UNRESOLVED(errno,
622 "Unable to alloc memory for measure saving");
623 }
624 m_tmp->nthreads = nth;
625 m_tmp->next = NULL;
626
627 /* Run the test */
628 dur = do_threads_test(nth, m_tmp);
629
630 /* If test was success, add this measure to the list. Otherwise, free the mem */
631 if (dur >= 0) {
632 m_cur->next = m_tmp;
633 m_cur = m_tmp;
634 } else {
635 free(m_tmp);
636 }
637 }
638 while (dur >= 0);
639
640 /* We will now parse the results to determine if the measure is ~ constant or is growing. */
641
642 ret = parse_measure(&sentinel);
643
644 /* Free the memory from the list */
645 m_cur = sentinel.next;
646 while (m_cur != NULL) {
647 m_tmp = m_cur;
648 m_cur = m_tmp->next;
649 free(m_tmp);
650 }
651
652 if (ret != 0) {
653 FAILED("This function is not scalable");
654 }
655 #if VERBOSE > 0
656 output("The function is scalable\n");
657 #endif
658
659 PASSED;
660 }
661
662 /***
663 * The next function will seek for the better model for each series of measurements.
664 *
665 * The tested models are: -- X = # threads; Y = latency
666 * -> Y = a; -- Error is r1 = avg((Y - Yavg)²);
667 * -> Y = aX + b; -- Error is r2 = avg((Y -aX -b)²);
668 * -- where a = avg ((X - Xavg)(Y - Yavg)) / avg((X - Xavg)²)
669 * -- Note: We will call _q = sum((X - Xavg) * (Y - Yavg));
670 * -- and _d = sum((X - Xavg)²);
671 * -- and b = Yavg - a * Xavg
672 * -> Y = c * X^a;-- Same as previous, but with log(Y) = a log(X) + b; and b = log(c). Error is r3
673 * -> Y = exp(aX + b); -- log(Y) = aX + b. Error is r4
674 *
675 * We compute each error factor (r1, r2, r3, r4) then search which is the smallest (with ponderation).
676 * The function returns 0 when r1 is the best for all cases (latency is constant) and !0 otherwise.
677 */
678
679 struct row {
680 long X; /* the X values -- copied from function argument */
681 long Y[NSCENAR]; /* the Y values -- copied from function argument */
682 long _x; /* Value X - Xavg */
683 long _y[NSCENAR]; /* Value Y - Yavg */
684 double LnX; /* Natural logarithm of X values */
685 double LnY[NSCENAR]; /* Natural logarithm of Y values */
686 double _lnx; /* Value LnX - LnXavg */
687 double _lny[NSCENAR]; /* Value LnY - LnYavg */
688 };
689
parse_measure(mes_t * measures)690 int parse_measure(mes_t * measures)
691 {
692 int ret, i, r;
693
694 mes_t *cur;
695
696 double Xavg, Yavg[NSCENAR];
697 double LnXavg, LnYavg[NSCENAR];
698
699 int N;
700
701 double r1[NSCENAR], r2[NSCENAR], r3[NSCENAR], r4[NSCENAR];
702
703 /* Some more intermediate vars */
704 long double _q[3][NSCENAR];
705 long double _d[3][NSCENAR];
706
707 long double t; /* temp value */
708
709 struct row *Table = NULL;
710
711 /* Initialize the datas */
712 Xavg = 0.0;
713 LnXavg = 0.0;
714 for (i = 0; i < NSCENAR; i++) {
715 Yavg[i] = 0.0;
716 LnYavg[i] = 0.0;
717 r1[i] = 0.0;
718 r2[i] = 0.0;
719 r3[i] = 0.0;
720 r4[i] = 0.0;
721 _q[0][i] = 0.0;
722 _q[1][i] = 0.0;
723 _q[2][i] = 0.0;
724 _d[0][i] = 0.0;
725 _d[1][i] = 0.0;
726 _d[2][i] = 0.0;
727 }
728 N = 0;
729 cur = measures;
730
731 #if VERBOSE > 1
732 output("Data analysis starting\n");
733 #endif
734
735 /* We start with reading the list to find:
736 * -> number of elements, to assign an array
737 * -> average values
738 */
739 while (cur->next != NULL) {
740 cur = cur->next;
741
742 N++;
743
744 Xavg += (double)cur->nthreads;
745 LnXavg += log((double)cur->nthreads);
746
747 for (i = 0; i < NSCENAR; i++) {
748 Yavg[i] += (double)cur->_data[i];
749 LnYavg[i] += log((double)cur->_data[i]);
750 }
751 }
752
753 /* We have the sum; we can divide to obtain the average values */
754 Xavg /= N;
755 LnXavg /= N;
756
757 for (i = 0; i < NSCENAR; i++) {
758 Yavg[i] /= N;
759 LnYavg[i] /= N;
760 }
761
762 #if VERBOSE > 1
763 output(" Found %d rows and %d columns\n", N, NSCENAR);
764 #endif
765
766 /* We will now alloc the array ... */
767 Table = calloc(N, sizeof(struct row));
768 if (Table == NULL) {
769 UNRESOLVED(errno, "Unable to alloc space for results parsing");
770 }
771
772 /* ... and fill it */
773 N = 0;
774 cur = measures;
775
776 while (cur->next != NULL) {
777 cur = cur->next;
778
779 Table[N].X = (long)cur->nthreads;
780 Table[N]._x = Table[N].X - Xavg;
781 Table[N].LnX = log((double)cur->nthreads);
782 Table[N]._lnx = Table[N].LnX - LnXavg;
783 for (i = 0; i < NSCENAR; i++) {
784 Table[N].Y[i] = cur->_data[i];
785 Table[N]._y[i] = Table[N].Y[i] - Yavg[i];
786 Table[N].LnY[i] = log((double)cur->_data[i]);
787 Table[N]._lny[i] = Table[N].LnY[i] - LnYavg[i];
788 }
789
790 N++;
791 }
792
793 /* We won't need the list anymore -- we'll work with the array which should be faster. */
794 #if VERBOSE > 1
795 output(" Data was stored in an array.\n");
796 #endif
797
798 /* We need to read the full array at least twice to compute all the error factors */
799
800 /* In the first pass, we'll compute:
801 * -> r1 for each scenar.
802 * -> "a" factor for linear (0), power (1) and exponential (2) approximations -- with using the _d and _q vars.
803 */
804 #if VERBOSE > 1
805 output("Starting first pass...\n");
806 #endif
807 for (r = 0; r < N; r++) {
808 for (i = 0; i < NSCENAR; i++) {
809 r1[i] +=
810 ((double)Table[r]._y[i] / N) *
811 (double)Table[r]._y[i];
812
813 _q[0][i] += Table[r]._y[i] * Table[r]._x;
814 _d[0][i] += Table[r]._x * Table[r]._x;
815
816 _q[1][i] += Table[r]._lny[i] * Table[r]._lnx;
817 _d[1][i] += Table[r]._lnx * Table[r]._lnx;
818
819 _q[2][i] += Table[r]._lny[i] * Table[r]._x;
820 _d[2][i] += Table[r]._x * Table[r]._x;
821 }
822 }
823
824 /* First pass is terminated; a2 = _q[0]/_d[0]; a3 = _q[1]/_d[1]; a4 = _q[2]/_d[2] */
825
826 /* In the first pass, we'll compute:
827 * -> r2, r3, r4 for each scenar.
828 */
829
830 #if VERBOSE > 1
831 output("Starting second pass...\n");
832 #endif
833 for (r = 0; r < N; r++) {
834 for (i = 0; i < NSCENAR; i++) {
835 /* r2 = avg((y - ax -b)²); t = (y - ax - b) = (y - yavg) - a (x - xavg); */
836 t = (Table[r]._y[i] -
837 ((_q[0][i] * Table[r]._x) / _d[0][i]));
838 r2[i] += t * t / N;
839
840 /* r3 = avg((y - c.x^a) ²);
841 t = y - c * x ^ a
842 = y - log (LnYavg - (_q[1]/_d[1]) * LnXavg) * x ^ (_q[1]/_d[1])
843 */
844 t = (Table[r].Y[i]
845 - (logl(LnYavg[i] - (_q[1][i] / _d[1][i]) * LnXavg)
846 * powl(Table[r].X, (_q[1][i] / _d[1][i]))
847 ));
848 r3[i] += t * t / N;
849
850 /* r4 = avg((y - exp(ax+b))²);
851 t = y - exp(ax+b)
852 = y - exp(_q[2]/_d[2] * x + (LnYavg - (_q[2]/_d[2] * Xavg)));
853 = y - exp(_q[2]/_d[2] * (x - Xavg) + LnYavg);
854 */
855 t = (Table[r].Y[i]
856 - expl((_q[2][i] / _d[2][i]) * Table[r]._x +
857 LnYavg[i]));
858 r4[i] += t * t / N;
859
860 }
861 }
862
863 #if VERBOSE > 1
864 output("All computing terminated.\n");
865 #endif
866 ret = 0;
867 for (i = 0; i < NSCENAR; i++) {
868 #if VERBOSE > 1
869 output("\nScenario: %s\n", test_scenar[i].desc);
870
871 output(" Model: Y = k\n");
872 output(" k = %g\n", Yavg[i]);
873 output(" Divergence %g\n", r1[i]);
874
875 output(" Model: Y = a * X + b\n");
876 output(" a = %Lg\n", _q[0][i] / _d[0][i]);
877 output(" b = %Lg\n",
878 Yavg[i] - ((_q[0][i] / _d[0][i]) * Xavg));
879 output(" Divergence %g\n", r2[i]);
880
881 output(" Model: Y = c * X ^ a\n");
882 output(" a = %Lg\n", _q[1][i] / _d[1][i]);
883 output(" c = %Lg\n",
884 logl(LnYavg[i] - (_q[1][i] / _d[1][i]) * LnXavg));
885 output(" Divergence %g\n", r2[i]);
886
887 output(" Model: Y = exp(a * X + b)\n");
888 output(" a = %Lg\n", _q[2][i] / _d[2][i]);
889 output(" b = %Lg\n",
890 LnYavg[i] - ((_q[2][i] / _d[2][i]) * Xavg));
891 output(" Divergence %g\n", r2[i]);
892 #endif
893 /* Compare r1 to other values, with some ponderations */
894 if ((r1[i] > 1.1 * r2[i]) || (r1[i] > 1.2 * r3[i])
895 || (r1[i] > 1.3 * r4[i]))
896 ret++;
897 #if VERBOSE > 1
898 else
899 output(" Sanction: OK\n");
900 #endif
901 }
902
903 /* We need to free the array */
904 free(Table);
905
906 /* We're done */
907 return ret;
908 }
909