Home | History | Annotate | Download | only in pthread_create 
      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  * This scalability sample aims to test the following assertions:
     18  *  -> If pthread_create fails because of a lack of a resource, or
     19  	PTHREAD_THREADS_MAX threads already exist, EAGAIN shall be returned.
     20  *  -> It also tests that the thread creation time does not depend on the # of threads
     21  *     already created.
     22 
     23  * The steps are:
     24  * -> Create threads until failure
     25 
     26  */
     27 
     28  /* We are testing conformance to IEEE Std 1003.1, 2003 Edition */
     29 #define _POSIX_C_SOURCE 200112L
     30 
     31  /* Some routines are part of the XSI Extensions */
     32 #ifndef WITHOUT_XOPEN
     33 #define _XOPEN_SOURCE	600
     34 #endif
     35 /********************************************************************************************/
     36 /****************************** standard includes *****************************************/
     37 /********************************************************************************************/
     38 #include <pthread.h>
     39 #include <stdarg.h>
     40 #include <stdio.h>
     41 #include <stdlib.h>
     42 #include <string.h>
     43 #include <unistd.h>
     44 
     45 #include <sched.h>
     46 #include <semaphore.h>
     47 #include <errno.h>
     48 #include <assert.h>
     49 #include <sys/wait.h>
     50 #include <math.h>
     51 
     52 /********************************************************************************************/
     53 /******************************   Test framework   *****************************************/
     54 /********************************************************************************************/
     55 #include "testfrmw.h"
     56 #include "testfrmw.c"
     57  /* This header is responsible for defining the following macros:
     58   * UNRESOLVED(ret, descr);
     59   *    where descr is a description of the error and ret is an int (error code for example)
     60   * FAILED(descr);
     61   *    where descr is a short text saying why the test has failed.
     62   * PASSED();
     63   *    No parameter.
     64   *
     65   * Both three macros shall terminate the calling process.
     66   * The testcase shall not terminate in any other maneer.
     67   *
     68   * The other file defines the functions
     69   * void output_init()
     70   * void output(char * string, ...)
     71   *
     72   * Those may be used to output information.
     73   */
     74 
     75 /********************************************************************************************/
     76 /********************************** Configuration ******************************************/
     77 /********************************************************************************************/
     78 #ifndef SCALABILITY_FACTOR
     79 #define SCALABILITY_FACTOR 1
     80 #endif
     81 #ifndef VERBOSE
     82 #define VERBOSE 1
     83 #endif
     84 
     85 #define RESOLUTION (5 * SCALABILITY_FACTOR)
     86 
     87 #ifdef PLOT_OUTPUT
     88 #undef VERBOSE
     89 #define VERBOSE 0
     90 #endif
     91 
     92 /********************************************************************************************/
     93 /***********************************    Test cases  *****************************************/
     94 /********************************************************************************************/
     95 
     96 #include "threads_scenarii.c"
     97 
     98 /* This file will define the following objects:
     99  * scenarii: array of struct __scenario type.
    100  * NSCENAR : macro giving the total # of scenarii
    101  * scenar_init(): function to call before use the scenarii array.
    102  * scenar_fini(): function to call after end of use of the scenarii array.
    103  */
    104 
    105 /********************************************************************************************/
    106 /***********************************    Real Test   *****************************************/
    107 /********************************************************************************************/
    108 
    109 /* The next structure is used to save the tests measures */
    110 typedef struct __mes_t {
    111 	int nthreads;
    112 	long _data[NSCENAR];	/* As we store sec values, a long type should be amply enough. */
    113 	struct __mes_t *next;
    114 } mes_t;
    115 
    116 /* Forward declaration */
    117 int parse_measure(mes_t * measures);
    118 
    119 pthread_mutex_t m_synchro = PTHREAD_MUTEX_INITIALIZER;
    120 
    121 void *threaded(void *arg)
    122 {
    123 	int ret = 0;
    124 
    125 	/* Signal we're done */
    126 	do {
    127 		ret = sem_post(&scenarii[sc].sem);
    128 	}
    129 	while ((ret == -1) && (errno == EINTR));
    130 	if (ret == -1) {
    131 		UNRESOLVED(errno, "Failed to wait for the semaphore");
    132 	}
    133 
    134 	/* Wait for all threads being created */
    135 	ret = pthread_mutex_lock(&m_synchro);
    136 	if (ret != 0) {
    137 		UNRESOLVED(ret, "Mutex lock failed");
    138 	}
    139 	(*(int *)arg) += 1;
    140 	ret = pthread_mutex_unlock(&m_synchro);
    141 	if (ret != 0) {
    142 		UNRESOLVED(ret, "Mutex unlock failed");
    143 	}
    144 
    145 	return arg;
    146 }
    147 
    148 int main(int argc, char *argv[])
    149 {
    150 	int ret = 0;
    151 	pthread_t child;
    152 	pthread_t *th;
    153 
    154 	int nthreads, ctl, i, tmp;
    155 
    156 	struct timespec ts_ref, ts_fin;
    157 
    158 	mes_t sentinel;
    159 	mes_t *m_cur, *m_tmp;
    160 
    161 	long PTHREAD_THREADS_MAX = sysconf(_SC_THREAD_THREADS_MAX);
    162 	long my_max = 1000 * SCALABILITY_FACTOR;
    163 
    164 	/* Initialize the measure list */
    165 	m_cur = &sentinel;
    166 	m_cur->next = NULL;
    167 
    168 	/* Initialize output routine */
    169 	output_init();
    170 
    171 	if (PTHREAD_THREADS_MAX > 0)
    172 		my_max = PTHREAD_THREADS_MAX;
    173 
    174 	th = (pthread_t *) calloc(1 + my_max, sizeof(pthread_t));
    175 	if (th == NULL) {
    176 		UNRESOLVED(errno, "Not enough memory for thread storage");
    177 	}
    178 
    179 	/* Initialize thread attribute objects */
    180 	scenar_init();
    181 
    182 #ifdef PLOT_OUTPUT
    183 	printf("# COLUMNS %d #threads", NSCENAR + 1);
    184 	for (sc = 0; sc < NSCENAR; sc++)
    185 		printf(" %i", sc);
    186 	printf("\n");
    187 #endif
    188 
    189 	for (sc = 0; sc < NSCENAR; sc++) {
    190 		if (scenarii[sc].bottom == NULL) {	/* skip the alternate stacks as we could create only 1 */
    191 #if VERBOSE > 0
    192 			output("-----\n");
    193 			output("Starting test with scenario (%i): %s\n", sc,
    194 			       scenarii[sc].descr);
    195 #endif
    196 
    197 			/* Block every (about to be) created threads */
    198 			ret = pthread_mutex_lock(&m_synchro);
    199 			if (ret != 0) {
    200 				UNRESOLVED(ret, "Mutex lock failed");
    201 			}
    202 
    203 			ctl = 0;
    204 			nthreads = 0;
    205 			m_cur = &sentinel;
    206 
    207 			/* Create 1 thread for testing purpose */
    208 			ret =
    209 			    pthread_create(&child, &scenarii[sc].ta, threaded,
    210 					   &ctl);
    211 			switch (scenarii[sc].result) {
    212 			case 0:	/* Operation was expected to succeed */
    213 				if (ret != 0) {
    214 					UNRESOLVED(ret,
    215 						   "Failed to create this thread");
    216 				}
    217 				break;
    218 
    219 			case 1:	/* Operation was expected to fail */
    220 				if (ret == 0) {
    221 					UNRESOLVED(-1,
    222 						   "An error was expected but the thread creation succeeded");
    223 				}
    224 				break;
    225 
    226 			case 2:	/* We did not know the expected result */
    227 			default:
    228 #if VERBOSE > 0
    229 				if (ret == 0) {
    230 					output
    231 					    ("Thread has been created successfully for this scenario\n");
    232 				} else {
    233 					output
    234 					    ("Thread creation failed with the error: %s\n",
    235 					     strerror(ret));
    236 				}
    237 #endif
    238 				;
    239 			}
    240 			if (ret == 0) {	/* The new thread is running */
    241 
    242 				while (1) {	/* we will break */
    243 					/* read clock */
    244 					ret =
    245 					    clock_gettime(CLOCK_REALTIME,
    246 							  &ts_ref);
    247 					if (ret != 0) {
    248 						UNRESOLVED(errno,
    249 							   "Unable to read clock");
    250 					}
    251 
    252 					/* create a new thread */
    253 					ret =
    254 					    pthread_create(&th[nthreads],
    255 							   &scenarii[sc].ta,
    256 							   threaded, &ctl);
    257 
    258 					/* stop here if we've got EAGAIN */
    259 					if (ret == EAGAIN)
    260 						break;
    261 
    262 // temporary hack
    263 					if (ret == ENOMEM)
    264 						break;
    265 					nthreads++;
    266 
    267 					/* FAILED if error is != EAGAIN or nthreads > PTHREAD_THREADS_MAX */
    268 					if (ret != 0) {
    269 						output
    270 						    ("pthread_create returned: %i (%s)\n",
    271 						     ret, strerror(ret));
    272 						FAILED
    273 						    ("pthread_create did not return EAGAIN on a lack of resource");
    274 					}
    275 					if (nthreads > my_max) {
    276 						if (PTHREAD_THREADS_MAX > 0) {
    277 							FAILED
    278 							    ("We were able to create more than PTHREAD_THREADS_MAX threads");
    279 						} else {
    280 							break;
    281 						}
    282 					}
    283 
    284 					/* wait for the semaphore */
    285 					do {
    286 						ret =
    287 						    sem_wait(&scenarii[sc].sem);
    288 					}
    289 					while ((ret == -1) && (errno == EINTR));
    290 					if (ret == -1) {
    291 						UNRESOLVED(errno,
    292 							   "Failed to wait for the semaphore");
    293 					}
    294 
    295 					/* read clock */
    296 					ret =
    297 					    clock_gettime(CLOCK_REALTIME,
    298 							  &ts_fin);
    299 					if (ret != 0) {
    300 						UNRESOLVED(errno,
    301 							   "Unable to read clock");
    302 					}
    303 
    304 					/* add to the measure list if nthreads % resolution == 0 */
    305 					if ((nthreads % RESOLUTION) == 0) {
    306 						if (m_cur->next == NULL) {
    307 							/* Create an empty new element */
    308 							m_tmp = malloc(sizeof(mes_t));
    309 							if (m_tmp == NULL) {
    310 								UNRESOLVED
    311 								    (errno,
    312 								     "Unable to alloc memory for measure saving");
    313 							}
    314 							m_tmp->nthreads =
    315 							    nthreads;
    316 							m_tmp->next = NULL;
    317 							for (tmp = 0;
    318 							     tmp < NSCENAR;
    319 							     tmp++)
    320 								m_tmp->
    321 								    _data[tmp] =
    322 								    0;
    323 							m_cur->next = m_tmp;
    324 						}
    325 
    326 						/* Add this measure to the next element */
    327 						m_cur = m_cur->next;
    328 						m_cur->_data[sc] =
    329 						    ((ts_fin.tv_sec -
    330 						      ts_ref.tv_sec) *
    331 						     1000000) +
    332 						    ((ts_fin.tv_nsec -
    333 						      ts_ref.tv_nsec) / 1000);
    334 
    335 #if VERBOSE > 5
    336 						output
    337 						    ("Added the following measure: sc=%i, n=%i, v=%li\n",
    338 						     sc, nthreads,
    339 						     m_cur->_data[sc]);
    340 #endif
    341 					}
    342 				}
    343 #if VERBOSE > 3
    344 				output
    345 				    ("Could not create anymore thread. Current count is %i\n",
    346 				     nthreads);
    347 #endif
    348 
    349 				/* Unblock every created threads */
    350 				ret = pthread_mutex_unlock(&m_synchro);
    351 				if (ret != 0) {
    352 					UNRESOLVED(ret, "Mutex unlock failed");
    353 				}
    354 
    355 				if (scenarii[sc].detached == 0) {
    356 #if VERBOSE > 3
    357 					output("Joining the threads\n");
    358 #endif
    359 					for (i = 0; i < nthreads; i++) {
    360 						ret = pthread_join(th[i], NULL);
    361 						if (ret != 0) {
    362 							UNRESOLVED(ret,
    363 								   "Unable to join a thread");
    364 						}
    365 					}
    366 
    367 					ret = pthread_join(child, NULL);
    368 					if (ret != 0) {
    369 						UNRESOLVED(ret,
    370 							   "Unalbe to join a thread");
    371 					}
    372 
    373 				}
    374 #if VERBOSE > 3
    375 				output
    376 				    ("Waiting for threads (almost) termination\n");
    377 #endif
    378 				do {
    379 					ret = pthread_mutex_lock(&m_synchro);
    380 					if (ret != 0) {
    381 						UNRESOLVED(ret,
    382 							   "Mutex lock failed");
    383 					}
    384 
    385 					tmp = ctl;
    386 
    387 					ret = pthread_mutex_unlock(&m_synchro);
    388 					if (ret != 0) {
    389 						UNRESOLVED(ret,
    390 							   "Mutex unlock failed");
    391 					}
    392 				} while (tmp != nthreads + 1);
    393 
    394 			}
    395 			/* The thread was created */
    396 		}
    397 	}			/* next scenario */
    398 
    399 	/* Free some memory before result parsing */
    400 	free(th);
    401 
    402 	/* Compute the results */
    403 	ret = parse_measure(&sentinel);
    404 
    405 	/* Free the resources and output the results */
    406 
    407 #if VERBOSE > 5
    408 	printf("Dump : \n");
    409 	printf("%8.8s", "nth");
    410 	for (i = 0; i < NSCENAR; i++)
    411 		printf("|   %2.2i   ", i);
    412 	printf("\n");
    413 #endif
    414 	while (sentinel.next != NULL) {
    415 		m_cur = sentinel.next;
    416 #if (VERBOSE > 5) || defined(PLOT_OUTPUT)
    417 		printf("%8.8i", m_cur->nthreads);
    418 		for (i = 0; i < NSCENAR; i++)
    419 			printf(" %1.1li.%6.6li", m_cur->_data[i] / 1000000,
    420 			       m_cur->_data[i] % 1000000);
    421 		printf("\n");
    422 #endif
    423 		sentinel.next = m_cur->next;
    424 		free(m_cur);
    425 	}
    426 
    427 	scenar_fini();
    428 
    429 #if VERBOSE > 0
    430 	output("-----\n");
    431 	output("All test data destroyed\n");
    432 	output("Test PASSED\n");
    433 #endif
    434 
    435 	PASSED;
    436 }
    437 
    438 /***
    439  * The next function will seek for the better model for each series of measurements.
    440  *
    441  * The tested models are: -- X = # threads; Y = latency
    442  * -> Y = a;      -- Error is r1 = avg((Y - Yavg));
    443  * -> Y = aX + b; -- Error is r2 = avg((Y -aX -b));
    444  *                -- where a = avg ((X - Xavg)(Y - Yavg)) / avg((X - Xavg))
    445  *                --         Note: We will call _q = sum((X - Xavg) * (Y - Yavg));
    446  *                --                       and  _d = sum((X - Xavg));
    447  *                -- and   b = Yavg - a * Xavg
    448  * -> Y = c * X^a;-- Same as previous, but with log(Y) = a log(X) + b; and b = log(c). Error is r3
    449  * -> Y = exp(aX + b); -- log(Y) = aX + b. Error is r4
    450  *
    451  * We compute each error factor (r1, r2, r3, r4) then search which is the smallest (with ponderation).
    452  * The function returns 0 when r1 is the best for all cases (latency is constant) and !0 otherwise.
    453  */
    454 
    455 struct row {
    456 	long X;			/* the X values -- copied from function argument */
    457 	long Y[NSCENAR];	/* the Y values -- copied from function argument */
    458 	long _x[NSCENAR];	/* Value X - Xavg */
    459 	long _y[NSCENAR];	/* Value Y - Yavg */
    460 	double LnX;		/* Natural logarithm of X values */
    461 	double LnY[NSCENAR];	/* Natural logarithm of Y values */
    462 	double _lnx[NSCENAR];	/* Value LnX - LnXavg */
    463 	double _lny[NSCENAR];	/* Value LnY - LnYavg */
    464 };
    465 
    466 int parse_measure(mes_t * measures)
    467 {
    468 	int ret, i, r;
    469 
    470 	mes_t *cur;
    471 
    472 	double Xavg[NSCENAR], Yavg[NSCENAR];
    473 	double LnXavg[NSCENAR], LnYavg[NSCENAR];
    474 
    475 	int N;
    476 
    477 	double r1[NSCENAR], r2[NSCENAR], r3[NSCENAR], r4[NSCENAR];
    478 
    479 	/* Some more intermediate vars */
    480 	long double _q[3][NSCENAR];
    481 	long double _d[3][NSCENAR];
    482 
    483 	long double t;		/* temp value */
    484 
    485 	struct row *Table = NULL;
    486 
    487 	/* This array contains the last element of each serie */
    488 	int array_max[NSCENAR];
    489 
    490 	/* Initialize the datas */
    491 	for (i = 0; i < NSCENAR; i++) {
    492 		array_max[i] = -1;	/* means no data */
    493 		Xavg[i] = 0.0;
    494 		LnXavg[i] = 0.0;
    495 		Yavg[i] = 0.0;
    496 		LnYavg[i] = 0.0;
    497 		r1[i] = 0.0;
    498 		r2[i] = 0.0;
    499 		r3[i] = 0.0;
    500 		r4[i] = 0.0;
    501 		_q[0][i] = 0.0;
    502 		_q[1][i] = 0.0;
    503 		_q[2][i] = 0.0;
    504 		_d[0][i] = 0.0;
    505 		_d[1][i] = 0.0;
    506 		_d[2][i] = 0.0;
    507 	}
    508 	N = 0;
    509 	cur = measures;
    510 
    511 #if VERBOSE > 1
    512 	output("Data analysis starting\n");
    513 #endif
    514 
    515 	/* We start with reading the list to find:
    516 	 * -> number of elements, to assign an array.
    517 	 * -> average values
    518 	 */
    519 	while (cur->next != NULL) {
    520 		cur = cur->next;
    521 
    522 		N++;
    523 
    524 		for (i = 0; i < NSCENAR; i++) {
    525 			if (cur->_data[i] != 0) {
    526 				array_max[i] = N;
    527 				Xavg[i] += (double)cur->nthreads;
    528 				LnXavg[i] += log((double)cur->nthreads);
    529 				Yavg[i] += (double)cur->_data[i];
    530 				LnYavg[i] += log((double)cur->_data[i]);
    531 			}
    532 		}
    533 	}
    534 
    535 	/* We have the sum; we can divide to obtain the average values */
    536 	for (i = 0; i < NSCENAR; i++) {
    537 		if (array_max[i] != -1) {
    538 			Xavg[i] /= array_max[i];
    539 			LnXavg[i] /= array_max[i];
    540 			Yavg[i] /= array_max[i];
    541 			LnYavg[i] /= array_max[i];
    542 		}
    543 	}
    544 
    545 #if VERBOSE > 1
    546 	output(" Found %d rows and %d columns\n", N, NSCENAR);
    547 #endif
    548 
    549 	/* We will now alloc the array ... */
    550 	Table = calloc(N, sizeof(struct row));
    551 	if (Table == NULL) {
    552 		UNRESOLVED(errno, "Unable to alloc space for results parsing");
    553 	}
    554 
    555 	/* ... and fill it */
    556 	N = 0;
    557 	cur = measures;
    558 
    559 	while (cur->next != NULL) {
    560 		cur = cur->next;
    561 
    562 		Table[N].X = (long)cur->nthreads;
    563 		Table[N].LnX = log((double)cur->nthreads);
    564 		for (i = 0; i < NSCENAR; i++) {
    565 			if (array_max[i] > N) {
    566 				Table[N]._x[i] = Table[N].X - Xavg[i];
    567 				Table[N]._lnx[i] = Table[N].LnX - LnXavg[i];
    568 				Table[N].Y[i] = cur->_data[i];
    569 				Table[N]._y[i] = Table[N].Y[i] - Yavg[i];
    570 				Table[N].LnY[i] = log((double)cur->_data[i]);
    571 				Table[N]._lny[i] = Table[N].LnY[i] - LnYavg[i];
    572 			}
    573 		}
    574 
    575 		N++;
    576 	}
    577 
    578 	/* We won't need the list anymore -- we'll work with the array which should be faster. */
    579 #if VERBOSE > 1
    580 	output(" Data was stored in an array.\n");
    581 #endif
    582 
    583 	/* We need to read the full array at least twice to compute all the error factors */
    584 
    585 	/* In the first pass, we'll compute:
    586 	 * -> r1 for each scenar.
    587 	 * -> "a" factor for linear (0), power (1) and exponential (2) approximations -- with using the _d and _q vars.
    588 	 */
    589 #if VERBOSE > 1
    590 	output("Starting first pass...\n");
    591 #endif
    592 	for (i = 0; i < NSCENAR; i++) {
    593 		for (r = 0; r < array_max[i]; r++) {
    594 			r1[i] +=
    595 			    ((double)Table[r]._y[i] / array_max[i]) *
    596 			    (double)Table[r]._y[i];
    597 
    598 			_q[0][i] += Table[r]._y[i] * Table[r]._x[i];
    599 			_d[0][i] += Table[r]._x[i] * Table[r]._x[i];
    600 
    601 			_q[1][i] += Table[r]._lny[i] * Table[r]._lnx[i];
    602 			_d[1][i] += Table[r]._lnx[i] * Table[r]._lnx[i];
    603 
    604 			_q[2][i] += Table[r]._lny[i] * Table[r]._x[i];
    605 			_d[2][i] += Table[r]._x[i] * Table[r]._x[i];
    606 		}
    607 	}
    608 
    609 	/* First pass is terminated; a2 = _q[0]/_d[0]; a3 = _q[1]/_d[1]; a4 = _q[2]/_d[2] */
    610 
    611 	/* In the first pass, we'll compute:
    612 	 * -> r2, r3, r4 for each scenar.
    613 	 */
    614 
    615 #if VERBOSE > 1
    616 	output("Starting second pass...\n");
    617 #endif
    618 	for (i = 0; i < NSCENAR; i++) {
    619 		for (r = 0; r < array_max[i]; r++) {
    620 			/* r2 = avg((y - ax -b));  t = (y - ax - b) = (y - yavg) - a (x - xavg); */
    621 			t = (Table[r]._y[i] -
    622 			     ((_q[0][i] * Table[r]._x[i]) / _d[0][i]));
    623 			r2[i] += t * t / array_max[i];
    624 
    625 			/* r3 = avg((y - c.x^a) );
    626 			   t = y - c * x ^ a
    627 			   = y - log (LnYavg - (_q[1]/_d[1]) * LnXavg) * x ^ (_q[1]/_d[1])
    628 			 */
    629 			t = (Table[r].Y[i]
    630 			     -
    631 			     (logl
    632 			      (LnYavg[i] - (_q[1][i] / _d[1][i]) * LnXavg[i])
    633 			      * powl(Table[r].X, (_q[1][i] / _d[1][i]))
    634 			     ));
    635 			r3[i] += t * t / array_max[i];
    636 
    637 			/* r4 = avg((y - exp(ax+b)));
    638 			   t = y - exp(ax+b)
    639 			   = y - exp(_q[2]/_d[2] * x + (LnYavg - (_q[2]/_d[2] * Xavg)));
    640 			   = y - exp(_q[2]/_d[2] * (x - Xavg) + LnYavg);
    641 			 */
    642 			t = (Table[r].Y[i]
    643 			     - expl((_q[2][i] / _d[2][i]) * Table[r]._x[i] +
    644 				    LnYavg[i]));
    645 			r4[i] += t * t / array_max[i];
    646 
    647 		}
    648 	}
    649 
    650 #if VERBOSE > 1
    651 	output("All computing terminated.\n");
    652 #endif
    653 	ret = 0;
    654 	for (i = 0; i < NSCENAR; i++) {
    655 #if VERBOSE > 1
    656 		output("\nScenario: %s\n", scenarii[i].descr);
    657 
    658 		output(" # of data: %i\n", array_max[i]);
    659 
    660 		output("  Model: Y = k\n");
    661 		output("       k = %g\n", Yavg[i]);
    662 		output("    Divergence %g\n", r1[i]);
    663 
    664 		output("  Model: Y = a * X + b\n");
    665 		output("       a = %Lg\n", _q[0][i] / _d[0][i]);
    666 		output("       b = %Lg\n",
    667 		       Yavg[i] - ((_q[0][i] / _d[0][i]) * Xavg[i]));
    668 		output("    Divergence %g\n", r2[i]);
    669 
    670 		output("  Model: Y = c * X ^ a\n");
    671 		output("       a = %Lg\n", _q[1][i] / _d[1][i]);
    672 		output("       c = %Lg\n",
    673 		       logl(LnYavg[i] - (_q[1][i] / _d[1][i]) * LnXavg[i]));
    674 		output("    Divergence %g\n", r2[i]);
    675 
    676 		output("  Model: Y = exp(a * X + b)\n");
    677 		output("       a = %Lg\n", _q[2][i] / _d[2][i]);
    678 		output("       b = %Lg\n",
    679 		       LnYavg[i] - ((_q[2][i] / _d[2][i]) * Xavg[i]));
    680 		output("    Divergence %g\n", r2[i]);
    681 #endif
    682 
    683 		if (array_max[i] != -1) {
    684 			/* Compare r1 to other values, with some ponderations */
    685 			if ((r1[i] > 1.1 * r2[i]) || (r1[i] > 1.2 * r3[i])
    686 			    || (r1[i] > 1.3 * r4[i]))
    687 				ret++;
    688 #if VERBOSE > 1
    689 			else
    690 				output(" Sanction: OK\n");
    691 #endif
    692 		}
    693 	}
    694 
    695 	/* We need to free the array */
    696 	free(Table);
    697 
    698 	/* We're done */
    699 	return ret;
    700 }
    701