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