projects / lttv.git / blob
? search:
summary | shortlog | log | commit | commitdiff | tree
blame | history | raw | HEAD
Fix two memory leaks
[lttv.git] / lttv / lttv / sync / event_analysis_chull.c
1 /* This file is part of the Linux Trace Toolkit viewer
2 * Copyright (C) 2009, 2010 Benjamin Poirier <benjamin.poirier@polymtl.ca>
3 *
4 * This program is free software: you can redistribute it and/or modify it
5 * under the terms of the GNU Lesser General Public License as published by
6 * the Free Software Foundation, either version 2.1 of the License, or (at
7 * your option) any later version.
8 *
9 * This program is distributed in the hope that it will be useful, but WITHOUT
10 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12 * License for more details.
13 *
14 * You should have received a copy of the GNU Lesser General Public License
15 * along with this program. If not, see <http://www.gnu.org/licenses/>.
16 */
17 #define _ISOC99_SOURCE
18
19 #ifdef HAVE_CONFIG_H
20 #include <config.h>
21 #endif
22
23 #include <errno.h>
24 #include <inttypes.h>
25 #include <math.h>
26 #include <float.h>
27 #include <stdlib.h>
28 #include <stdio.h>
29 #include <string.h>
30 #include <unistd.h>
31
32 #include "sync_chain.h"
33
34 #include "event_analysis_chull.h"
35
36
37 typedef enum
38 {
39 LOWER,
40 UPPER
41 } HullType;
42
43 typedef enum
44 {
45 MINIMUM,
46 MAXIMUM
47 } LineType;
48
49 #ifdef HAVE_LIBGLPK
50 struct LPAddRowInfo
51 {
52 glp_prob* lp;
53 int boundType;
54 GArray* iArray, * jArray, * aArray;
55 };
56 #endif
57
58
59 // Functions common to all analysis modules
60 static void initAnalysisCHull(SyncState* const syncState);
61 static void destroyAnalysisCHull(SyncState* const syncState);
62
63 static void analyzeMessageCHull(SyncState* const syncState, Message* const
64 message);
65 static AllFactors* finalizeAnalysisCHull(SyncState* const syncState);
66 static void printAnalysisStatsCHull(SyncState* const syncState);
67 static void writeAnalysisTraceTraceForePlotsCHull(SyncState* const syncState,
68 const unsigned int i, const unsigned int j);
69
70 // Functions specific to this module
71 static void openGraphFiles(SyncState* const syncState);
72 static void closeGraphFiles(SyncState* const syncState);
73 static void writeGraphFiles(SyncState* const syncState);
74 static void gfDumpHullToFile(gpointer data, gpointer userData);
75
76 AllFactors* calculateAllFactors(struct _SyncState* const syncState);
77 void calculateFactorsMiddle(PairFactors* const factors);
78 static Factors* calculateFactorsExact(GQueue* const cu, GQueue* const cl, const
79 LineType lineType) __attribute__((pure));
80 static void calculateFactorsFallback(GQueue* const cr, GQueue* const cs,
81 PairFactors* const result);
82 static void grahamScan(GQueue* const hull, Point* const newPoint, const
83 HullType type);
84 static int jointCmp(const Point* const p1, const Point* const p2, const Point*
85 const p3) __attribute__((pure));
86 static double crossProductK(const Point const* p1, const Point const* p2,
87 const Point const* p3, const Point const* p4) __attribute__((pure));
88 static double slope(const Point* const p1, const Point* const p2)
89 __attribute__((pure));
90 static double intercept(const Point* const p1, const Point* const p2)
91 __attribute__((pure));
92 static double verticalDistance(Point* p1, Point* p2, Point* const point)
93 __attribute__((pure));
94
95 static void gfPointDestroy(gpointer data, gpointer userData);
96
97 // The next group of functions is only needed when computing synchronization
98 // accuracy.
99 #ifdef HAVE_LIBGLPK
100 static AllFactors* finalizeAnalysisCHullLP(SyncState* const syncState);
101 static void writeAnalysisTraceTimeBackPlotsCHull(SyncState* const syncState,
102 const unsigned int i, const unsigned int j);
103 static void writeAnalysisTraceTimeForePlotsCHull(SyncState* const syncState,
104 const unsigned int i, const unsigned int j);
105 static void writeAnalysisTraceTraceBackPlotsCHull(SyncState* const syncState,
106 const unsigned int i, const unsigned int j);
107
108 static glp_prob* lpCreateProblem(GQueue* const lowerHull, GQueue* const
109 upperHull);
110 static void gfLPAddRow(gpointer data, gpointer user_data);
111 static Factors* calculateFactorsLP(glp_prob* const lp, const int direction);
112 static void calculateCompleteFactorsLP(glp_prob* const lp, PairFactors*
113 factors);
114 void timeCorrectionLP(glp_prob* const lp, const PairFactors* const lpFactors,
115 const uint64_t time, CorrectedTime* const correctedTime);
116
117 static void gfAddAbsiscaToArray(gpointer data, gpointer user_data);
118 static gint gcfCompareUint64(gconstpointer a, gconstpointer b);
119 #else
120 static inline AllFactors* finalizeAnalysisCHullLP(SyncState* const syncState)
121 {
122 return NULL;
123 }
124 #endif
125
126
127
128 static AnalysisModule analysisModuleCHull= {
129 .name= "chull",
130 .initAnalysis= &initAnalysisCHull,
131 .destroyAnalysis= &destroyAnalysisCHull,
132 .analyzeMessage= &analyzeMessageCHull,
133 .finalizeAnalysis= &finalizeAnalysisCHull,
134 .printAnalysisStats= &printAnalysisStatsCHull,
135 .graphFunctions= {
136 #ifdef HAVE_LIBGLPK
137 .writeTraceTimeBackPlots= &writeAnalysisTraceTimeBackPlotsCHull,
138 .writeTraceTimeForePlots= &writeAnalysisTraceTimeForePlotsCHull,
139 .writeTraceTraceBackPlots= &writeAnalysisTraceTraceBackPlotsCHull,
140 #endif
141 .writeTraceTraceForePlots= &writeAnalysisTraceTraceForePlotsCHull,
142 }
143 };
144
145
146 /*
147 * Analysis module registering function
148 */
149 void registerAnalysisCHull()
150 {
151 g_queue_push_tail(&analysisModules, &analysisModuleCHull);
152 }
153
154
155 /*
156 * Analysis init function
157 *
158 * This function is called at the beginning of a synchronization run for a set
159 * of traces.
160 *
161 * Allocate some of the analysis specific data structures
162 *
163 * Args:
164 * syncState container for synchronization data.
165 * This function allocates or initializes these analysisData
166 * members:
167 * hullArray
168 * dropped
169 */
170 static void initAnalysisCHull(SyncState* const syncState)
171 {
172 unsigned int i, j;
173 AnalysisDataCHull* analysisData;
174
175 analysisData= malloc(sizeof(AnalysisDataCHull));
176 syncState->analysisData= analysisData;
177
178 analysisData->hullArray= malloc(syncState->traceNb * sizeof(GQueue**));
179 for (i= 0; i < syncState->traceNb; i++)
180 {
181 analysisData->hullArray[i]= malloc(syncState->traceNb * sizeof(GQueue*));
182
183 for (j= 0; j < syncState->traceNb; j++)
184 {
185 analysisData->hullArray[i][j]= g_queue_new();
186 }
187 }
188 #ifdef HAVE_LIBGLPK
189 analysisData->lps= NULL;
190 #endif
191
192 if (syncState->stats)
193 {
194 analysisData->stats= calloc(1, sizeof(AnalysisStatsCHull));
195 }
196
197 if (syncState->graphsStream)
198 {
199 analysisData->graphsData= calloc(1, sizeof(AnalysisGraphsDataCHull));
200 openGraphFiles(syncState);
201 }
202 }
203
204
205 /*
206 * Create and open files used to store convex hull points to genereate
207 * graphs. Allocate and populate array to store file pointers.
208 *
209 * Args:
210 * syncState: container for synchronization data
211 */
212 static void openGraphFiles(SyncState* const syncState)
213 {
214 unsigned int i, j;
215 int retval;
216 char* cwd;
217 char name[31];
218 AnalysisDataCHull* analysisData;
219
220 analysisData= (AnalysisDataCHull*) syncState->analysisData;
221
222 cwd= changeToGraphsDir(syncState->graphsDir);
223
224 analysisData->graphsData->hullPoints= malloc(syncState->traceNb *
225 sizeof(FILE**));
226 for (i= 0; i < syncState->traceNb; i++)
227 {
228 analysisData->graphsData->hullPoints[i]= malloc(syncState->traceNb *
229 sizeof(FILE*));
230 for (j= 0; j < syncState->traceNb; j++)
231 {
232 if (i != j)
233 {
234 retval= snprintf(name, sizeof(name),
235 "analysis_chull-%03u_to_%03u.data", j, i);
236 if (retval > sizeof(name) - 1)
237 {
238 name[sizeof(name) - 1]= '\0';
239 }
240 if ((analysisData->graphsData->hullPoints[i][j]= fopen(name, "w")) ==
241 NULL)
242 {
243 g_error(strerror(errno));
244 }
245 }
246 }
247 }
248
249 retval= chdir(cwd);
250 if (retval == -1)
251 {
252 g_error(strerror(errno));
253 }
254 free(cwd);
255 }
256
257
258 /*
259 * Write hull points to files to generate graphs.
260 *
261 * Args:
262 * syncState: container for synchronization data
263 */
264 static void writeGraphFiles(SyncState* const syncState)
265 {
266 unsigned int i, j;
267 AnalysisDataCHull* analysisData;
268
269 analysisData= (AnalysisDataCHull*) syncState->analysisData;
270
271 for (i= 0; i < syncState->traceNb; i++)
272 {
273 for (j= 0; j < syncState->traceNb; j++)
274 {
275 if (i != j)
276 {
277 g_queue_foreach(analysisData->hullArray[i][j],
278 &gfDumpHullToFile,
279 analysisData->graphsData->hullPoints[i][j]);
280 }
281 }
282 }
283 }
284
285
286 /*
287 * A GFunc for g_queue_foreach. Write a hull point to a file used to generate
288 * graphs
289 *
290 * Args:
291 * data: Point*, point to write to the file
292 * userData: FILE*, file pointer where to write the point
293 */
294 static void gfDumpHullToFile(gpointer data, gpointer userData)
295 {
296 Point* point;
297
298 point= (Point*) data;
299 fprintf((FILE*) userData, "%20" PRIu64 " %20" PRIu64 "\n", point->x, point->y);
300 }
301
302
303 /*
304 * Close files used to store convex hull points to generate graphs.
305 * Deallocate array to store file pointers.
306 *
307 * Args:
308 * syncState: container for synchronization data
309 */
310 static void closeGraphFiles(SyncState* const syncState)
311 {
312 unsigned int i, j;
313 AnalysisDataCHull* analysisData;
314 int retval;
315
316 analysisData= (AnalysisDataCHull*) syncState->analysisData;
317
318 if (analysisData->graphsData->hullPoints == NULL)
319 {
320 return;
321 }
322
323 for (i= 0; i < syncState->traceNb; i++)
324 {
325 for (j= 0; j < syncState->traceNb; j++)
326 {
327 if (i != j)
328 {
329 retval= fclose(analysisData->graphsData->hullPoints[i][j]);
330 if (retval != 0)
331 {
332 g_error(strerror(errno));
333 }
334 }
335 }
336 free(analysisData->graphsData->hullPoints[i]);
337 }
338 free(analysisData->graphsData->hullPoints);
339 analysisData->graphsData->hullPoints= NULL;
340 }
341
342
343 /*
344 * Analysis destroy function
345 *
346 * Free the analysis specific data structures
347 *
348 * Args:
349 * syncState container for synchronization data.
350 * This function deallocates these analysisData members:
351 * hullArray
352 * stDev
353 */
354 static void destroyAnalysisCHull(SyncState* const syncState)
355 {
356 unsigned int i, j;
357 AnalysisDataCHull* analysisData;
358
359 analysisData= (AnalysisDataCHull*) syncState->analysisData;
360
361 if (analysisData == NULL)
362 {
363 return;
364 }
365
366 for (i= 0; i < syncState->traceNb; i++)
367 {
368 for (j= 0; j < syncState->traceNb; j++)
369 {
370 g_queue_foreach(analysisData->hullArray[i][j], gfPointDestroy,
371 NULL);
372 g_queue_free(analysisData->hullArray[i][j]);
373 }
374 free(analysisData->hullArray[i]);
375 }
376 free(analysisData->hullArray);
377
378 #ifdef HAVE_LIBGLPK
379 if (analysisData->lps != NULL)
380 {
381 for (i= 0; i < syncState->traceNb; i++)
382 {
383 unsigned int j;
384
385 for (j= 0; j < i; j++)
386 {
387 glp_delete_prob(analysisData->lps[i][j]);
388 }
389 free(analysisData->lps[i]);
390 }
391 free(analysisData->lps);
392
393 /* Be careful, this invalidates all problem objects which still exist.
394 * Don't keep copies of lps past this point. */
395 glp_free_env();
396 }
397 #endif
398
399 if (syncState->stats)
400 {
401 freeAllFactors(analysisData->stats->allFactors, syncState->traceNb);
402 freeAllFactors(analysisData->stats->geoFactors, syncState->traceNb);
403
404 #ifdef HAVE_LIBGLPK
405 freeAllFactors(analysisData->stats->lpFactors, syncState->traceNb);
406 #endif
407
408 free(analysisData->stats);
409 }
410
411 if (syncState->graphsStream)
412 {
413 AnalysisGraphsDataCHull* graphs= analysisData->graphsData;
414
415 if (graphs->hullPoints != NULL)
416 {
417 closeGraphFiles(syncState);
418 }
419
420 freeAllFactors(graphs->allFactors, syncState->traceNb);
421
422 #ifdef HAVE_LIBGLPK
423 freeAllFactors(graphs->lpFactors, syncState->traceNb);
424 #endif
425
426 free(analysisData->graphsData);
427 }
428
429 free(syncState->analysisData);
430 syncState->analysisData= NULL;
431 }
432
433
434 /*
435 * Perform analysis on an event pair.
436 *
437 * Args:
438 * syncState container for synchronization data
439 * message structure containing the events
440 */
441 static void analyzeMessageCHull(SyncState* const syncState, Message* const message)
442 {
443 AnalysisDataCHull* analysisData;
444 Point* newPoint;
445 HullType hullType;
446 GQueue* hull;
447
448 analysisData= (AnalysisDataCHull*) syncState->analysisData;
449
450 newPoint= malloc(sizeof(Point));
451 if (message->inE->traceNum < message->outE->traceNum)
452 {
453 // CA is inE->traceNum
454 newPoint->x= message->inE->cpuTime;
455 newPoint->y= message->outE->cpuTime;
456 hullType= UPPER;
457 g_debug("Reception point hullArray[%lu][%lu] "
458 "x= inE->time= %" PRIu64 " y= outE->time= %" PRIu64,
459 message->inE->traceNum, message->outE->traceNum, newPoint->x,
460 newPoint->y);
461 }
462 else
463 {
464 // CA is outE->traceNum
465 newPoint->x= message->outE->cpuTime;
466 newPoint->y= message->inE->cpuTime;
467 hullType= LOWER;
468 g_debug("Send point hullArray[%lu][%lu] "
469 "x= inE->time= %" PRIu64 " y= outE->time= %" PRIu64,
470 message->inE->traceNum, message->outE->traceNum, newPoint->x,
471 newPoint->y);
472 }
473
474 hull=
475 analysisData->hullArray[message->inE->traceNum][message->outE->traceNum];
476
477 if (hull->length >= 1 && newPoint->x < ((Point*)
478 g_queue_peek_tail(hull))->x)
479 {
480 if (syncState->stats)
481 {
482 analysisData->stats->dropped++;
483 }
484
485 free(newPoint);
486 }
487 else
488 {
489 grahamScan(hull, newPoint, hullType);
490 }
491 }
492
493
494 /*
495 * Construct one half of a convex hull from abscissa-sorted points
496 *
497 * Args:
498 * hull: the points already in the hull
499 * newPoint: a new point to consider
500 * type: which half of the hull to construct
501 */
502 static void grahamScan(GQueue* const hull, Point* const newPoint, const
503 HullType type)
504 {
505 int inversionFactor;
506
507 g_debug("grahamScan(hull (length: %u), newPoint, %s)", hull->length, type
508 == LOWER ? "LOWER" : "UPPER");
509
510 if (type == LOWER)
511 {
512 inversionFactor= 1;
513 }
514 else
515 {
516 inversionFactor= -1;
517 }
518
519 if (hull->length >= 2)
520 {
521 g_debug("jointCmp(hull[%u], hull[%u], newPoint) * inversionFactor = %d * %d = %d",
522 hull->length - 2,
523 hull->length - 1,
524 jointCmp(g_queue_peek_nth(hull, hull->length - 2),
525 g_queue_peek_tail(hull), newPoint),
526 inversionFactor,
527 jointCmp(g_queue_peek_nth(hull, hull->length - 2),
528 g_queue_peek_tail(hull), newPoint) * inversionFactor);
529 }
530 while (hull->length >= 2 && jointCmp(g_queue_peek_nth(hull, hull->length -
531 2), g_queue_peek_tail(hull), newPoint) * inversionFactor <= 0)
532 {
533 g_debug("Removing hull[%u]", hull->length);
534 free((Point*) g_queue_pop_tail(hull));
535
536 if (hull->length >= 2)
537 {
538 g_debug("jointCmp(hull[%u], hull[%u], newPoint) * inversionFactor = %d * %d = %d",
539 hull->length - 2,
540 hull->length - 1,
541 jointCmp(g_queue_peek_nth(hull, hull->length - 2),
542 g_queue_peek_tail(hull), newPoint),
543 inversionFactor,
544 jointCmp(g_queue_peek_nth(hull, hull->length - 2),
545 g_queue_peek_tail(hull), newPoint) * inversionFactor);
546 }
547 }
548 g_queue_push_tail(hull, newPoint);
549 }
550
551
552 /*
553 * Finalize the factor calculations
554 *
555 * Args:
556 * syncState container for synchronization data.
557 *
558 * Returns:
559 * AllFactors* synchronization factors for each trace pair, the caller is
560 * responsible for freeing the structure
561 */
562 static AllFactors* finalizeAnalysisCHull(SyncState* const syncState)
563 {
564 AnalysisDataCHull* analysisData;
565 AllFactors* geoFactors, * lpFactors;
566
567 analysisData= (AnalysisDataCHull*) syncState->analysisData;
568
569 if (syncState->graphsStream && analysisData->graphsData->hullPoints != NULL)
570 {
571 writeGraphFiles(syncState);
572 closeGraphFiles(syncState);
573 }
574
575 geoFactors= calculateAllFactors(syncState);
576 lpFactors= finalizeAnalysisCHullLP(syncState);
577
578 if (syncState->stats)
579 {
580 geoFactors->refCount++;
581 analysisData->stats->geoFactors= geoFactors;
582
583 if (lpFactors != NULL)
584 {
585 lpFactors->refCount++;
586 analysisData->stats->allFactors= lpFactors;
587 }
588 else
589 {
590 geoFactors->refCount++;
591 analysisData->stats->allFactors= geoFactors;
592 }
593 }
594
595 if (syncState->graphsStream)
596 {
597 if (lpFactors != NULL)
598 {
599 lpFactors->refCount++;
600 analysisData->graphsData->allFactors= lpFactors;
601 }
602 else
603 {
604 geoFactors->refCount++;
605 analysisData->graphsData->allFactors= geoFactors;
606 }
607 }
608
609 if (lpFactors != NULL)
610 {
611 freeAllFactors(geoFactors, syncState->traceNb);
612 return lpFactors;
613 }
614 else
615 {
616 freeAllFactors(lpFactors, syncState->traceNb);
617 return geoFactors;
618 }
619 }
620
621
622 /*
623 * Print statistics related to analysis. Must be called after
624 * finalizeAnalysis.
625 *
626 * Args:
627 * syncState container for synchronization data.
628 */
629 static void printAnalysisStatsCHull(SyncState* const syncState)
630 {
631 AnalysisDataCHull* analysisData;
632 unsigned int i, j;
633
634 if (!syncState->stats)
635 {
636 return;
637 }
638
639 analysisData= (AnalysisDataCHull*) syncState->analysisData;
640
641 printf("Convex hull analysis stats:\n");
642 printf("\tout of order packets dropped from analysis: %u\n",
643 analysisData->stats->dropped);
644
645 printf("\tNumber of points in convex hulls:\n");
646
647 for (i= 0; i < syncState->traceNb; i++)
648 {
649 for (j= i + 1; j < syncState->traceNb; j++)
650 {
651 printf("\t\t%3d - %-3d: lower half-hull %-5u upper half-hull %-5u\n",
652 i, j, analysisData->hullArray[j][i]->length,
653 analysisData->hullArray[i][j]->length);
654 }
655 }
656
657 printf("\tIndividual synchronization factors:\n");
658
659 for (i= 0; i < syncState->traceNb; i++)
660 {
661 for (j= i + 1; j < syncState->traceNb; j++)
662 {
663 PairFactors* factorsCHull;
664
665 factorsCHull= &analysisData->stats->allFactors->pairFactors[j][i];
666 printf("\t\t%3d - %-3d: %s", i, j,
667 approxNames[factorsCHull->type]);
668
669 if (factorsCHull->type == EXACT)
670 {
671 printf(" a0= % 7g a1= 1 %c %7g\n",
672 factorsCHull->approx->offset,
673 factorsCHull->approx->drift < 0. ? '-' : '+',
674 fabs(factorsCHull->approx->drift));
675 }
676 else if (factorsCHull->type == ACCURATE)
677 {
678 printf("\n\t\t a0: % 7g to % 7g (delta= %7g)\n",
679 factorsCHull->max->offset, factorsCHull->min->offset,
680 factorsCHull->min->offset - factorsCHull->max->offset);
681 printf("\t\t a1: 1 %+7g to %+7g (delta= %7g)\n",
682 factorsCHull->min->drift - 1., factorsCHull->max->drift -
683 1., factorsCHull->max->drift - factorsCHull->min->drift);
684 }
685 else if (factorsCHull->type == APPROXIMATE)
686 {
687 printf(" a0= % 7g a1= 1 %c %7g error= %7g\n",
688 factorsCHull->approx->offset, factorsCHull->approx->drift
689 - 1. < 0. ? '-' : '+', fabs(factorsCHull->approx->drift -
690 1.), factorsCHull->accuracy);
691 }
692 else if (factorsCHull->type == INCOMPLETE)
693 {
694 printf("\n");
695
696 if (factorsCHull->min->drift != -INFINITY)
697 {
698 printf("\t\t min: a0: % 7g a1: 1 %c %7g\n",
699 factorsCHull->min->offset, factorsCHull->min->drift -
700 1. < 0 ? '-' : '+', fabs(factorsCHull->min->drift -
701 1.));
702 }
703 if (factorsCHull->max->drift != INFINITY)
704 {
705 printf("\t\t max: a0: % 7g a1: 1 %c %7g\n",
706 factorsCHull->max->offset, factorsCHull->max->drift -
707 1. < 0 ? '-' : '+', fabs(factorsCHull->max->drift -
708 1.));
709 }
710 }
711 else if (factorsCHull->type == FAIL)
712 {
713 printf("\n");
714
715 if (factorsCHull->min != NULL && factorsCHull->min->drift != -INFINITY)
716 {
717 printf("\t\t min: a0: % 7g a1: 1 %c %7g\n",
718 factorsCHull->min->offset, factorsCHull->min->drift -
719 1. < 0 ? '-' : '+', fabs(factorsCHull->min->drift -
720 1.));
721 }
722 if (factorsCHull->max != NULL && factorsCHull->max->drift != INFINITY)
723 {
724 printf("\t\t max: a0: % 7g a1: 1 %c %7g\n",
725 factorsCHull->max->offset, factorsCHull->max->drift -
726 1. < 0 ? '-' : '+', fabs(factorsCHull->max->drift -
727 1.));
728 }
729 }
730 else if (factorsCHull->type == ABSENT)
731 {
732 printf("\n");
733 }
734 else
735 {
736 g_assert_not_reached();
737 }
738 }
739 }
740
741 #ifdef HAVE_LIBGLPK
742 printf("\tFactor comparison:\n"
743 "\t\tTrace pair Factors type Differences (lp - chull)\n"
744 "\t\t a0 a1\n"
745 "\t\t Min Max Min Max\n");
746
747 for (i= 0; i < syncState->traceNb; i++)
748 {
749 for (j= 0; j < i; j++)
750 {
751 PairFactors* geoFactors=
752 &analysisData->stats->geoFactors->pairFactors[i][j];
753 PairFactors* lpFactors=
754 &analysisData->stats->lpFactors->pairFactors[i][j];
755
756 printf("\t\t%3d - %-3d ", i, j);
757 if (lpFactors->type == geoFactors->type)
758 {
759 if (lpFactors->type == ACCURATE)
760 {
761 printf("%-13s %-10.4g %-10.4g %-10.4g %.4g\n",
762 approxNames[lpFactors->type],
763 lpFactors->min->offset - geoFactors->min->offset,
764 lpFactors->max->offset - geoFactors->max->offset,
765 lpFactors->min->drift - geoFactors->min->drift,
766 lpFactors->max->drift - geoFactors->max->drift);
767 }
768 else if (lpFactors->type == ABSENT)
769 {
770 printf("%s\n", approxNames[lpFactors->type]);
771 }
772 }
773 else
774 {
775 printf("Different! %s and %s\n", approxNames[lpFactors->type],
776 approxNames[geoFactors->type]);
777 }
778 }
779 }
780 #endif
781 }
782
783
784 /*
785 * A GFunc for g_queue_foreach()
786 *
787 * Args:
788 * data Point*, point to destroy
789 * user_data NULL
790 */
791 static void gfPointDestroy(gpointer data, gpointer userData)
792 {
793 Point* point;
794
795 point= (Point*) data;
796 free(point);
797 }
798
799
800 /*
801 * Find out if a sequence of three points constitutes a "left turn" or a
802 * "right turn".
803 *
804 * Args:
805 * p1, p2, p3: The three points.
806 *
807 * Returns:
808 * < 0 right turn
809 * 0 colinear (unlikely result since this uses floating point
810 * arithmetic)
811 * > 0 left turn
812 */
813 static int jointCmp(const Point const* p1, const Point const* p2, const
814 Point const* p3)
815 {
816 double result;
817 const double fuzzFactor= 0.;
818
819 result= crossProductK(p1, p2, p1, p3);
820 g_debug("crossProductK(p1= (%" PRIu64 ", %" PRIu64 "), "
821 "p2= (%" PRIu64 ", %" PRIu64 "), p1= (%" PRIu64 ", %" PRIu64 "), "
822 "p3= (%" PRIu64 ", %" PRIu64 "))= %g",
823 p1->x, p1->y, p2->x, p2->y, p1->x, p1->y, p3->x, p3->y, result);
824 if (result < fuzzFactor)
825 {
826 return -1;
827 }
828 else if (result > fuzzFactor)
829 {
830 return 1;
831 }
832 else
833 {
834 return 0;
835 }
836 }
837
838
839 /*
840 * Calculate the k component of the cross product of two vectors.
841 *
842 * Args:
843 * p1, p2: start and end points of the first vector
844 * p3, p4: start and end points of the second vector
845 *
846 * Returns:
847 * the k component of the cross product when considering the two vectors to
848 * be in the i-j plane. The direction (sign) of the result can be useful to
849 * determine the relative orientation of the two vectors.
850 */
851 static double crossProductK(const Point const* p1, const Point const* p2,
852 const Point const* p3, const Point const* p4)
853 {
854 return ((double) p2->x - p1->x) * ((double) p4->y - p3->y) - ((double)
855 p2->y - p1->y) * ((double) p4->x - p3->x);
856 }
857
858
859 /*
860 * Analyze the convex hulls to determine the synchronization factors between
861 * each pair of trace.
862 *
863 * Args:
864 * syncState container for synchronization data.
865 *
866 * Returns:
867 * AllFactors*, see the documentation for the member geoFactors of
868 * AnalysisStatsCHull.
869 */
870 AllFactors* calculateAllFactors(SyncState* const syncState)
871 {
872 unsigned int traceNumA, traceNumB;
873 AllFactors* geoFactors;
874 AnalysisDataCHull* analysisData;
875
876 analysisData= (AnalysisDataCHull*) syncState->analysisData;
877
878 // Allocate geoFactors and calculate min and max
879 geoFactors= createAllFactors(syncState->traceNb);
880 for (traceNumA= 0; traceNumA < syncState->traceNb; traceNumA++)
881 {
882 for (traceNumB= 0; traceNumB < traceNumA; traceNumB++)
883 {
884 unsigned int i;
885 GQueue* cs, * cr;
886 const struct
887 {
888 LineType lineType;
889 size_t factorsOffset;
890 } loopValues[]= {
891 {MINIMUM, offsetof(PairFactors, min)},
892 {MAXIMUM, offsetof(PairFactors, max)}
893 };
894
895 cr= analysisData->hullArray[traceNumB][traceNumA];
896 cs= analysisData->hullArray[traceNumA][traceNumB];
897
898 for (i= 0; i < sizeof(loopValues) / sizeof(*loopValues); i++)
899 {
900 g_debug("geoFactors[%u][%u].%s = calculateFactorsExact(cr= "
901 "hullArray[%u][%u], cs= hullArray[%u][%u], %s)",
902 traceNumA, traceNumB, loopValues[i].factorsOffset ==
903 offsetof(PairFactors, min) ? "min" : "max", traceNumB,
904 traceNumA, traceNumA, traceNumB, loopValues[i].lineType ==
905 MINIMUM ? "MINIMUM" : "MAXIMUM");
906 *((Factors**) ((void*)
907 &geoFactors->pairFactors[traceNumA][traceNumB] +
908 loopValues[i].factorsOffset))=
909 calculateFactorsExact(cr, cs, loopValues[i].lineType);
910 }
911 }
912 }
913
914 // Calculate approx when possible
915 for (traceNumA= 0; traceNumA < syncState->traceNb; traceNumA++)
916 {
917 for (traceNumB= 0; traceNumB < traceNumA; traceNumB++)
918 {
919 PairFactors* factorsCHull;
920
921 factorsCHull= &geoFactors->pairFactors[traceNumA][traceNumB];
922 if (factorsCHull->min == NULL && factorsCHull->max == NULL)
923 {
924 factorsCHull->type= APPROXIMATE;
925 calculateFactorsFallback(analysisData->hullArray[traceNumB][traceNumA],
926 analysisData->hullArray[traceNumA][traceNumB],
927 &geoFactors->pairFactors[traceNumA][traceNumB]);
928 }
929 else if (factorsCHull->min != NULL && factorsCHull->max != NULL)
930 {
931 if (factorsCHull->min->drift != -INFINITY &&
932 factorsCHull->max->drift != INFINITY)
933 {
934 factorsCHull->type= ACCURATE;
935 calculateFactorsMiddle(factorsCHull);
936 }
937 else if (factorsCHull->min->drift != -INFINITY ||
938 factorsCHull->max->drift != INFINITY)
939 {
940 factorsCHull->type= INCOMPLETE;
941 }
942 else
943 {
944 factorsCHull->type= ABSENT;
945 }
946 }
947 else
948 {
949 //g_assert_not_reached();
950 factorsCHull->type= FAIL;
951 }
952 }
953 }
954
955 return geoFactors;
956 }
957
958
959 /* Calculate approximative factors based on minimum and maximum limits. The
960 * best approximation to make is the interior bissector of the angle formed by
961 * the minimum and maximum lines.
962 *
963 * The formulae used come from [Haddad, Yoram: Performance dans les systèmes
964 * répartis: des outils pour les mesures, Université de Paris-Sud, Centre
965 * d'Orsay, September 1988] Section 6.1 p.44
966 *
967 * The reasoning for choosing this estimator comes from [Duda, A., Harrus, G.,
968 * Haddad, Y., and Bernard, G.: Estimating global time in distributed systems,
969 * Proc. 7th Int. Conf. on Distributed Computing Systems, Berlin, volume 18,
970 * 1987] p.303
971 *
972 * Args:
973 * factors: contains the min and max limits, used to store the result
974 */
975 void calculateFactorsMiddle(PairFactors* const factors)
976 {
977 double amin, amax, bmin, bmax, bhat;
978
979 amin= factors->max->offset;
980 amax= factors->min->offset;
981 bmin= factors->min->drift;
982 bmax= factors->max->drift;
983
984 g_assert_cmpfloat(bmax, >=, bmin);
985
986 factors->approx= malloc(sizeof(Factors));
987 bhat= (bmax * bmin - 1. + sqrt(1. + pow(bmax, 2.) * pow(bmin, 2.) +
988 pow(bmax, 2.) + pow(bmin, 2.))) / (bmax + bmin);
989 factors->approx->offset= amax - (amax - amin) / 2. * (pow(bhat, 2.) + 1.)
990 / (1. + bhat * bmax);
991 factors->approx->drift= bhat;
992 factors->accuracy= bmax - bmin;
993 }
994
995
996 /*
997 * Analyze the convex hulls to determine the minimum or maximum
998 * synchronization factors between one pair of trace.
999 *
1000 * This implements and improves upon the algorithm in [Haddad, Yoram:
1001 * Performance dans les systèmes répartis: des outils pour les mesures,
1002 * Université de Paris-Sud, Centre d'Orsay, September 1988] Section 6.2 p.47
1003 *
1004 * Some degenerate cases are possible:
1005 * 1) the result is unbounded. In that case, when searching for the maximum
1006 * factors, result->drift= INFINITY; result->offset= -INFINITY. When
1007 * searching for the minimum factors, it is the opposite. It is not
1008 * possible to improve the situation with this data.
1009 * 2) no line can be above the upper hull and below the lower hull. This is
1010 * because the hulls intersect each other or are reversed. This means that
1011 * an assertion was false. Most probably, the clocks are not linear. It is
1012 * possible to repeat the search with another algorithm that will find a
1013 * "best effort" approximation. See calculateFactorsApprox().
1014 *
1015 * Args:
1016 * cu: the upper half-convex hull, the line must pass above this
1017 * and touch it in one point
1018 * cl: the lower half-convex hull, the line must pass below this
1019 * and touch it in one point
1020 * lineType: search for minimum or maximum factors
1021 *
1022 * Returns:
1023 * If a result is found, a struct Factors is allocated, filed with the
1024 * result and returned
1025 * NULL otherwise, degenerate case 2 is in effect
1026 */
1027 static Factors* calculateFactorsExact(GQueue* const cu, GQueue* const cl, const
1028 LineType lineType)
1029 {
1030 GQueue* c1, * c2;
1031 unsigned int i1, i2;
1032 Point* p1, * p2;
1033 double inversionFactor;
1034 Factors* result;
1035
1036 g_debug("calculateFactorsExact(cu= %p, cl= %p, %s)", cu, cl, lineType ==
1037 MINIMUM ? "MINIMUM" : "MAXIMUM");
1038
1039 if (lineType == MINIMUM)
1040 {
1041 c1= cl;
1042 c2= cu;
1043 inversionFactor= -1.;
1044 }
1045 else
1046 {
1047 c1= cu;
1048 c2= cl;
1049 inversionFactor= 1.;
1050 }
1051
1052 i1= 0;
1053 i2= c2->length - 1;
1054
1055 // Check for degenerate case 1
1056 if (c1->length == 0 || c2->length == 0 || ((Point*) g_queue_peek_nth(c1,
1057 i1))->x >= ((Point*) g_queue_peek_nth(c2, i2))->x)
1058 {
1059 result= malloc(sizeof(Factors));
1060 if (lineType == MINIMUM)
1061 {
1062 result->drift= -INFINITY;
1063 result->offset= INFINITY;
1064 }
1065 else
1066 {
1067 result->drift= INFINITY;
1068 result->offset= -INFINITY;
1069 }
1070
1071 return result;
1072 }
1073
1074 do
1075 {
1076 while
1077 (
1078 (int) i2 - 1 > 0
1079 && crossProductK(
1080 g_queue_peek_nth(c1, i1),
1081 g_queue_peek_nth(c2, i2),
1082 g_queue_peek_nth(c1, i1),
1083 g_queue_peek_nth(c2, i2 - 1)) * inversionFactor < 0.
1084 )
1085 {
1086 if (((Point*) g_queue_peek_nth(c1, i1))->x
1087 < ((Point*) g_queue_peek_nth(c2, i2 - 1))->x)
1088 {
1089 i2--;
1090 }
1091 else
1092 {
1093 // Degenerate case 2
1094 return NULL;
1095 }
1096 }
1097 while
1098 (
1099 i1 + 1 < c1->length - 1
1100 && crossProductK(
1101 g_queue_peek_nth(c1, i1),
1102 g_queue_peek_nth(c2, i2),
1103 g_queue_peek_nth(c1, i1 + 1),
1104 g_queue_peek_nth(c2, i2)) * inversionFactor < 0.
1105 )
1106 {
1107 if (((Point*) g_queue_peek_nth(c1, i1 + 1))->x
1108 < ((Point*) g_queue_peek_nth(c2, i2))->x)
1109 {
1110 i1++;
1111 }
1112 else
1113 {
1114 // Degenerate case 2
1115 return NULL;
1116 }
1117 }
1118 } while
1119 (
1120 (int) i2 - 1 > 0
1121 && crossProductK(
1122 g_queue_peek_nth(c1, i1),
1123 g_queue_peek_nth(c2, i2),
1124 g_queue_peek_nth(c1, i1),
1125 g_queue_peek_nth(c2, i2 - 1)) * inversionFactor < 0.
1126 );
1127
1128 p1= g_queue_peek_nth(c1, i1);
1129 p2= g_queue_peek_nth(c2, i2);
1130
1131 g_debug("Resulting points are: c1[i1]: x= %" PRIu64 " y= %" PRIu64
1132 " c2[i2]: x= %" PRIu64 " y= %" PRIu64 "", p1->x, p1->y, p2->x, p2->y);
1133
1134 result= malloc(sizeof(Factors));
1135 result->drift= slope(p1, p2);
1136 result->offset= intercept(p1, p2);
1137
1138 g_debug("Resulting factors are: drift= %g offset= %g", result->drift,
1139 result->offset);
1140
1141 return result;
1142 }
1143
1144
1145 /*
1146 * Analyze the convex hulls to determine approximate synchronization factors
1147 * between one pair of trace when there is no line that can fit in the
1148 * corridor separating them.
1149 *
1150 * This implements the algorithm in [Ashton, P.: Algorithms for Off-line Clock
1151 * Synchronisation, University of Canterbury, December 1995, 26] Section 4.2.2
1152 * p.7
1153 *
1154 * For each point p1 in cr
1155 * For each point p2 in cs
1156 * errorMin= 0
1157 * Calculate the line paramaters
1158 * For each point p3 in each convex hull
1159 * If p3 is on the wrong side of the line
1160 * error+= distance
1161 * If error < errorMin
1162 * Update results
1163 *
1164 * Args:
1165 * cr: the upper half-convex hull
1166 * cs: the lower half-convex hull
1167 * result: a pointer to the pre-allocated struct where the results
1168 * will be stored
1169 */
1170 static void calculateFactorsFallback(GQueue* const cr, GQueue* const cs,
1171 PairFactors* const result)
1172 {
1173 unsigned int i, j, k;
1174 double errorMin;
1175 Factors* approx;
1176
1177 errorMin= INFINITY;
1178 approx= malloc(sizeof(Factors));
1179
1180 for (i= 0; i < cs->length; i++)
1181 {
1182 for (j= 0; j < cr->length; j++)
1183 {
1184 double error;
1185 Point p1, p2;
1186
1187 error= 0.;
1188
1189 if (((Point*) g_queue_peek_nth(cs, i))->x < ((Point*)g_queue_peek_nth(cr, j))->x)
1190 {
1191 p1= *(Point*)g_queue_peek_nth(cs, i);
1192 p2= *(Point*)g_queue_peek_nth(cr, j);
1193 }
1194 else
1195 {
1196 p1= *(Point*)g_queue_peek_nth(cr, j);
1197 p2= *(Point*)g_queue_peek_nth(cs, i);
1198 }
1199
1200 // The lower hull should be above the point
1201 for (k= 0; k < cs->length; k++)
1202 {
1203 if (jointCmp(&p1, &p2, g_queue_peek_nth(cs, k)) < 0.)
1204 {
1205 error+= verticalDistance(&p1, &p2, g_queue_peek_nth(cs, k));
1206 }
1207 }
1208
1209 // The upper hull should be below the point
1210 for (k= 0; k < cr->length; k++)
1211 {
1212 if (jointCmp(&p1, &p2, g_queue_peek_nth(cr, k)) > 0.)
1213 {
1214 error+= verticalDistance(&p1, &p2, g_queue_peek_nth(cr, k));
1215 }
1216 }
1217
1218 if (error < errorMin)
1219 {
1220 g_debug("Fallback: i= %u j= %u is a better match (error= %g)", i, j, error);
1221 approx->drift= slope(&p1, &p2);
1222 approx->offset= intercept(&p1, &p2);
1223 errorMin= error;
1224 }
1225 }
1226 }
1227
1228 result->approx= approx;
1229 result->accuracy= errorMin;
1230 }
1231
1232
1233 /*
1234 * Calculate the vertical distance between a line and a point
1235 *
1236 * Args:
1237 * p1, p2: Two points defining the line
1238 * point: a point
1239 *
1240 * Return:
1241 * the vertical distance
1242 */
1243 static double verticalDistance(Point* p1, Point* p2, Point* const point)
1244 {
1245 return fabs(slope(p1, p2) * point->x + intercept(p1, p2) - point->y);
1246 }
1247
1248
1249 /*
1250 * Calculate the slope between two points
1251 *
1252 * Args:
1253 * p1, p2 the two points
1254 *
1255 * Returns:
1256 * the slope
1257 */
1258 static double slope(const Point* const p1, const Point* const p2)
1259 {
1260 return ((double) p2->y - p1->y) / (p2->x - p1->x);
1261 }
1262
1263
1264 /* Calculate the y-intercept of a line that passes by two points
1265 *
1266 * Args:
1267 * p1, p2 the two points
1268 *
1269 * Returns:
1270 * the y-intercept
1271 */
1272 static double intercept(const Point* const p1, const Point* const p2)
1273 {
1274 return ((double) p2->y * p1->x - (double) p1->y * p2->x) / ((double) p1->x - p2->x);
1275 }
1276
1277
1278 /*
1279 * Write the analysis-specific graph lines in the gnuplot script.
1280 *
1281 * Args:
1282 * syncState: container for synchronization data
1283 * i: first trace number
1284 * j: second trace number, garanteed to be larger than i
1285 */
1286 void writeAnalysisTraceTraceForePlotsCHull(SyncState* const syncState, const unsigned
1287 int i, const unsigned int j)
1288 {
1289 AnalysisDataCHull* analysisData;
1290 PairFactors* factorsCHull;
1291
1292 analysisData= (AnalysisDataCHull*) syncState->analysisData;
1293
1294 fprintf(syncState->graphsStream,
1295 "\t\"analysis_chull-%1$03d_to_%2$03d.data\" "
1296 "title \"Lower half-hull\" with linespoints "
1297 "linecolor rgb \"#015a01\" linetype 4 pointtype 8 pointsize 0.8, \\\n"
1298 "\t\"analysis_chull-%2$03d_to_%1$03d.data\" "
1299 "title \"Upper half-hull\" with linespoints "
1300 "linecolor rgb \"#003366\" linetype 4 pointtype 10 pointsize 0.8, \\\n",
1301 i, j);
1302
1303 factorsCHull= &analysisData->graphsData->allFactors->pairFactors[j][i];
1304 if (factorsCHull->type == EXACT)
1305 {
1306 fprintf(syncState->graphsStream,
1307 "\t%7g + %7g * x "
1308 "title \"Exact conversion\" with lines "
1309 "linecolor rgb \"black\" linetype 1, \\\n",
1310 factorsCHull->approx->offset, factorsCHull->approx->drift);
1311 }
1312 else if (factorsCHull->type == ACCURATE)
1313 {
1314 fprintf(syncState->graphsStream,
1315 "\t%.2f + %.10f * x "
1316 "title \"Min conversion\" with lines "
1317 "linecolor rgb \"black\" linetype 5, \\\n",
1318 factorsCHull->min->offset, factorsCHull->min->drift);
1319 fprintf(syncState->graphsStream,
1320 "\t%.2f + %.10f * x "
1321 "title \"Max conversion\" with lines "
1322 "linecolor rgb \"black\" linetype 8, \\\n",
1323 factorsCHull->max->offset, factorsCHull->max->drift);
1324 fprintf(syncState->graphsStream,
1325 "\t%.2f + %.10f * x "
1326 "title \"Middle conversion\" with lines "
1327 "linecolor rgb \"black\" linetype 1, \\\n",
1328 factorsCHull->approx->offset, factorsCHull->approx->drift);
1329 }
1330 else if (factorsCHull->type == APPROXIMATE)
1331 {
1332 fprintf(syncState->graphsStream,
1333 "\t%.2f + %.10f * x "
1334 "title \"Fallback conversion\" with lines "
1335 "linecolor rgb \"gray60\" linetype 1, \\\n",
1336 factorsCHull->approx->offset, factorsCHull->approx->drift);
1337 }
1338 else if (factorsCHull->type == INCOMPLETE)
1339 {
1340 if (factorsCHull->min->drift != -INFINITY)
1341 {
1342 fprintf(syncState->graphsStream,
1343 "\t%.2f + %.10f * x "
1344 "title \"Min conversion\" with lines "
1345 "linecolor rgb \"black\" linetype 5, \\\n",
1346 factorsCHull->min->offset, factorsCHull->min->drift);
1347 }
1348
1349 if (factorsCHull->max->drift != INFINITY)
1350 {
1351 fprintf(syncState->graphsStream,
1352 "\t%.2f + %.10f * x "
1353 "title \"Max conversion\" with lines "
1354 "linecolor rgb \"black\" linetype 8, \\\n",
1355 factorsCHull->max->offset, factorsCHull->max->drift);
1356 }
1357 }
1358 else if (factorsCHull->type == FAIL)
1359 {
1360 if (factorsCHull->min != NULL && factorsCHull->min->drift != -INFINITY)
1361 {
1362 fprintf(syncState->graphsStream,
1363 "\t%.2f + %.10f * x "
1364 "title \"Min conversion\" with lines "
1365 "linecolor rgb \"black\" linetype 5, \\\n",
1366 factorsCHull->min->offset, factorsCHull->min->drift);
1367 }
1368
1369 if (factorsCHull->max != NULL && factorsCHull->max->drift != INFINITY)
1370 {
1371 fprintf(syncState->graphsStream,
1372 "\t%.2f + %.10f * x "
1373 "title \"Max conversion\" with lines "
1374 "linecolor rgb \"black\" linetype 8, \\\n",
1375 factorsCHull->max->offset, factorsCHull->max->drift);
1376 }
1377 }
1378 else if (factorsCHull->type == ABSENT)
1379 {
1380 }
1381 else
1382 {
1383 g_assert_not_reached();
1384 }
1385 }
1386
1387
1388 #ifdef HAVE_LIBGLPK
1389 /*
1390 * Create the linear programming problem containing the constraints defined by
1391 * two half-hulls. The objective function and optimization directions are not
1392 * written.
1393 *
1394 * Args:
1395 * syncState: container for synchronization data
1396 * i: first trace number
1397 * j: second trace number, garanteed to be larger than i
1398 * Returns:
1399 * A new glp_prob*, this problem must be freed by the caller with
1400 * glp_delete_prob()
1401 */
1402 static glp_prob* lpCreateProblem(GQueue* const lowerHull, GQueue* const
1403 upperHull)
1404 {
1405 unsigned int it;
1406 const int zero= 0;
1407 const double zeroD= 0.;
1408 glp_prob* lp= glp_create_prob();
1409 unsigned int hullPointNb= g_queue_get_length(lowerHull) +
1410 g_queue_get_length(upperHull);
1411 GArray* iArray= g_array_sized_new(FALSE, FALSE, sizeof(int), hullPointNb +
1412 1);
1413 GArray* jArray= g_array_sized_new(FALSE, FALSE, sizeof(int), hullPointNb +
1414 1);
1415 GArray* aArray= g_array_sized_new(FALSE, FALSE, sizeof(double),
1416 hullPointNb + 1);
1417 struct {
1418 GQueue* hull;
1419 struct LPAddRowInfo rowInfo;
1420 } loopValues[2]= {
1421 {lowerHull, {lp, GLP_UP, iArray, jArray, aArray}},
1422 {upperHull, {lp, GLP_LO, iArray, jArray, aArray}},
1423 };
1424
1425 // Create the LP problem
1426 glp_term_out(GLP_OFF);
1427 if (hullPointNb > 0)
1428 {
1429 glp_add_rows(lp, hullPointNb);
1430 }
1431 glp_add_cols(lp, 2);
1432
1433 glp_set_col_name(lp, 1, "a0");
1434 glp_set_col_bnds(lp, 1, GLP_FR, 0., 0.);
1435 glp_set_col_name(lp, 2, "a1");
1436 glp_set_col_bnds(lp, 2, GLP_LO, 0., 0.);
1437
1438 // Add row constraints
1439 g_array_append_val(iArray, zero);
1440 g_array_append_val(jArray, zero);
1441 g_array_append_val(aArray, zeroD);
1442
1443 for (it= 0; it < sizeof(loopValues) / sizeof(*loopValues); it++)
1444 {
1445 g_queue_foreach(loopValues[it].hull, &gfLPAddRow,
1446 &loopValues[it].rowInfo);
1447 }
1448
1449 g_assert_cmpuint(iArray->len, ==, jArray->len);
1450 g_assert_cmpuint(jArray->len, ==, aArray->len);
1451 g_assert_cmpuint(aArray->len - 1, ==, hullPointNb * 2);
1452
1453 glp_load_matrix(lp, aArray->len - 1, &g_array_index(iArray, int, 0),
1454 &g_array_index(jArray, int, 0), &g_array_index(aArray, double, 0));
1455
1456 glp_scale_prob(lp, GLP_SF_AUTO);
1457
1458 g_array_free(iArray, TRUE);
1459 g_array_free(jArray, TRUE);
1460 g_array_free(aArray, TRUE);
1461
1462 return lp;
1463 }
1464
1465
1466 /*
1467 * A GFunc for g_queue_foreach(). Add constraints and bounds for one row.
1468 *
1469 * Args:
1470 * data Point*, synchronization point for which to add an LP row
1471 * (a constraint)
1472 * user_data LPAddRowInfo*
1473 */
1474 static void gfLPAddRow(gpointer data, gpointer user_data)
1475 {
1476 Point* p= data;
1477 struct LPAddRowInfo* rowInfo= user_data;
1478 int indexes[2];
1479 double constraints[2];
1480
1481 indexes[0]= g_array_index(rowInfo->iArray, int, rowInfo->iArray->len - 1) + 1;
1482 indexes[1]= indexes[0];
1483
1484 if (rowInfo->boundType == GLP_UP)
1485 {
1486 glp_set_row_bnds(rowInfo->lp, indexes[0], GLP_UP, 0., p->y);
1487 }
1488 else if (rowInfo->boundType == GLP_LO)
1489 {
1490 glp_set_row_bnds(rowInfo->lp, indexes[0], GLP_LO, p->y, 0.);
1491 }
1492 else
1493 {
1494 g_assert_not_reached();
1495 }
1496
1497 g_array_append_vals(rowInfo->iArray, indexes, 2);
1498 indexes[0]= 1;
1499 indexes[1]= 2;
1500 g_array_append_vals(rowInfo->jArray, indexes, 2);
1501 constraints[0]= 1.;
1502 constraints[1]= p->x;
1503 g_array_append_vals(rowInfo->aArray, constraints, 2);
1504 }
1505
1506
1507 /*
1508 * Calculate min or max correction factors (as possible) using an LP problem.
1509 *
1510 * Args:
1511 * lp: A linear programming problem with constraints and bounds
1512 * initialized.
1513 * direction: The type of factors desired. Use GLP_MAX for max
1514 * approximation factors (a1, the drift or slope is the
1515 * largest) and GLP_MIN in the other case.
1516 *
1517 * Returns:
1518 * If the calculation was successful, a new Factors struct. Otherwise, NULL.
1519 * The calculation will fail if the hull assumptions are not respected.
1520 */
1521 static Factors* calculateFactorsLP(glp_prob* const lp, const int direction)
1522 {
1523 int retval, status;
1524 Factors* factors;
1525
1526 glp_set_obj_coef(lp, 1, 0.);
1527 glp_set_obj_coef(lp, 2, 1.);
1528
1529 glp_set_obj_dir(lp, direction);
1530 retval= glp_simplex(lp, NULL);
1531 status= glp_get_status(lp);
1532
1533 if (retval == 0 && status == GLP_OPT)
1534 {
1535 factors= malloc(sizeof(Factors));
1536 factors->offset= glp_get_col_prim(lp, 1);
1537 factors->drift= glp_get_col_prim(lp, 2);
1538 }
1539 else
1540 {
1541 factors= NULL;
1542 }
1543
1544 return factors;
1545 }
1546
1547
1548 /*
1549 * Calculate min, max and approx correction factors (as possible) using an LP
1550 * problem.
1551 *
1552 * Args:
1553 * lp A linear programming problem with constraints and bounds
1554 * initialized.
1555 * factors Resulting factors, must be preallocated
1556 */
1557 static void calculateCompleteFactorsLP(glp_prob* const lp, PairFactors* factors)
1558 {
1559 factors->min= calculateFactorsLP(lp, GLP_MIN);
1560 factors->max= calculateFactorsLP(lp, GLP_MAX);
1561
1562 if (factors->min && factors->max)
1563 {
1564 factors->type= ACCURATE;
1565 calculateFactorsMiddle(factors);
1566 }
1567 else if (factors->min || factors->max)
1568 {
1569 factors->type= INCOMPLETE;
1570 }
1571 else
1572 {
1573 factors->type= ABSENT;
1574 }
1575 }
1576
1577
1578 /*
1579 * A GFunc for g_queue_foreach()
1580 *
1581 * Args:
1582 * data Point*, a convex hull point
1583 * user_data GArray*, an array of convex hull point absisca values, as
1584 * uint64_t
1585 */
1586 static void gfAddAbsiscaToArray(gpointer data, gpointer user_data)
1587 {
1588 Point* p= data;
1589 GArray* a= user_data;
1590 uint64_t v= p->x;
1591
1592 g_array_append_val(a, v);
1593 }
1594
1595
1596 /*
1597 * A GCompareFunc for g_array_sort()
1598 *
1599 * Args:
1600 * a, b uint64_t*, absisca values
1601 *
1602 * Returns:
1603 * "returns less than zero for first arg is less than second arg, zero for
1604 * equal, greater zero if first arg is greater than second arg"
1605 * - the great glib documentation
1606 */
1607 static gint gcfCompareUint64(gconstpointer a, gconstpointer b)
1608 {
1609 if (*(uint64_t*) a < *(uint64_t*) b)
1610 {
1611 return -1;
1612 }
1613 else if (*(uint64_t*) a > *(uint64_t*) b)
1614 {
1615 return 1;
1616 }
1617 else
1618 {
1619 return 0;
1620 }
1621 }
1622
1623
1624 /*
1625 * Compute synchronization factors using a linear programming approach.
1626 *
1627 * Args:
1628 * syncState: container for synchronization data
1629 */
1630 static AllFactors* finalizeAnalysisCHullLP(SyncState* const syncState)
1631 {
1632 AnalysisDataCHull* analysisData= syncState->analysisData;
1633 unsigned int i, j;
1634 AllFactors* lpFactorsArray;
1635
1636 lpFactorsArray= createAllFactors(syncState->traceNb);
1637
1638 analysisData->lps= malloc(syncState->traceNb * sizeof(glp_prob**));
1639 for (i= 0; i < syncState->traceNb; i++)
1640 {
1641 analysisData->lps[i]= malloc(i * sizeof(glp_prob*));
1642 }
1643
1644 for (i= 0; i < syncState->traceNb; i++)
1645 {
1646 for (j= 0; j < i; j++)
1647 {
1648 glp_prob* lp;
1649 unsigned int it;
1650 GQueue*** hullArray= analysisData->hullArray;
1651 PairFactors* lpFactors= &lpFactorsArray->pairFactors[i][j];
1652
1653 // Create the LP problem
1654 lp= lpCreateProblem(hullArray[i][j], hullArray[j][i]);
1655 analysisData->lps[i][j]= lp;
1656
1657 // Use the LP problem to find the correction factors for this pair of
1658 // traces
1659 calculateCompleteFactorsLP(lp, lpFactors);
1660
1661 // If possible, compute synchronization accuracy information for
1662 // graphs
1663 if (syncState->graphsStream && lpFactors->type == ACCURATE)
1664 {
1665 int retval;
1666 char* cwd;
1667 char fileName[43];
1668 FILE* fp;
1669 GArray* xValues;
1670
1671 // Open the data file
1672 snprintf(fileName, sizeof(fileName),
1673 "analysis_chull_accuracy-%03u_and_%03u.data", j, i);
1674 fileName[sizeof(fileName) - 1]= '\0';
1675
1676 cwd= changeToGraphsDir(syncState->graphsDir);
1677
1678 if ((fp= fopen(fileName, "w")) == NULL)
1679 {
1680 g_error(strerror(errno));
1681 }
1682 fprintf(fp, "#%-24s %-25s %-25s %-25s\n", "x", "middle", "min", "max");
1683
1684 retval= chdir(cwd);
1685 if (retval == -1)
1686 {
1687 g_error(strerror(errno));
1688 }
1689 free(cwd);
1690
1691 // Build the list of absisca values for the points in the accuracy graph
1692 xValues= g_array_sized_new(FALSE, FALSE, sizeof(uint64_t),
1693 g_queue_get_length(hullArray[i][j]) +
1694 g_queue_get_length(hullArray[j][i]));
1695
1696 g_queue_foreach(hullArray[i][j], &gfAddAbsiscaToArray, xValues);
1697 g_queue_foreach(hullArray[j][i], &gfAddAbsiscaToArray, xValues);
1698
1699 g_array_sort(xValues, &gcfCompareUint64);
1700
1701 /* For each absisca value and each optimisation direction, solve the LP
1702 * and write a line in the data file */
1703 for (it= 0; it < xValues->len; it++)
1704 {
1705 uint64_t time;
1706 CorrectedTime correctedTime;
1707
1708 time= g_array_index(xValues, uint64_t, it);
1709 timeCorrectionLP(lp, lpFactors, time, &correctedTime);
1710 fprintf(fp, "%24" PRIu64 " %24" PRIu64 " %24" PRIu64
1711 "%24" PRIu64 "\n", time, correctedTime.time,
1712 correctedTime.min, correctedTime.max);
1713 }
1714
1715 g_array_free(xValues, TRUE);
1716 fclose(fp);
1717 }
1718 }
1719 }
1720
1721 if (syncState->stats)
1722 {
1723 lpFactorsArray->refCount++;
1724 analysisData->stats->lpFactors= lpFactorsArray;
1725 }
1726
1727 if (syncState->graphsStream)
1728 {
1729 lpFactorsArray->refCount++;
1730 analysisData->graphsData->lpFactors= lpFactorsArray;
1731 }
1732
1733 return lpFactorsArray;
1734 }
1735
1736
1737 /*
1738 * Perform correction on one time value and calculate accuracy bounds.
1739 *
1740 * Args:
1741 * lp: Linear Programming problem containing the coefficients for
1742 * the trace pair between which to perform time correction.
1743 * lpFactors: Correction factors for this trace pair, the factors must be
1744 * of type ACCURATE.
1745 * time: Time value to correct.
1746 * correctedTime: Result of the time correction, preallocated.
1747 */
1748 void timeCorrectionLP(glp_prob* const lp, const PairFactors* const lpFactors,
1749 const uint64_t time, CorrectedTime* const correctedTime)
1750 {
1751 unsigned int it;
1752 const struct
1753 {
1754 int direction;
1755 size_t offset;
1756 } loopValues[]= {
1757 {GLP_MIN, offsetof(CorrectedTime, min)},
1758 {GLP_MAX, offsetof(CorrectedTime, max)}
1759 };
1760
1761 glp_set_obj_coef(lp, 1, 1.);
1762 glp_set_obj_coef(lp, 2, time);
1763
1764 g_assert(lpFactors->type == ACCURATE);
1765
1766 correctedTime->time= lpFactors->approx->offset + lpFactors->approx->drift
1767 * time;
1768
1769 for (it= 0; it < ARRAY_SIZE(loopValues); it++)
1770 {
1771 int status;
1772 int retval;
1773
1774 glp_set_obj_dir(lp, loopValues[it].direction);
1775 retval= glp_simplex(lp, NULL);
1776 status= glp_get_status(lp);
1777
1778 g_assert(retval == 0 && status == GLP_OPT);
1779 *(uint64_t*) ((void*) correctedTime + loopValues[it].offset)=
1780 round(glp_get_obj_val(lp));
1781 }
1782 }
1783
1784
1785 /*
1786 * Write the analysis-specific graph lines in the gnuplot script.
1787 *
1788 * Args:
1789 * syncState: container for synchronization data
1790 * i: first trace number
1791 * j: second trace number, garanteed to be larger than i
1792 */
1793 static void writeAnalysisTraceTimeBackPlotsCHull(SyncState* const syncState,
1794 const unsigned int i, const unsigned int j)
1795 {
1796 if (((AnalysisDataCHull*)
1797 syncState->analysisData)->graphsData->lpFactors->pairFactors[j][i].type
1798 == ACCURATE)
1799 {
1800 fprintf(syncState->graphsStream,
1801 "\t\"analysis_chull_accuracy-%1$03u_and_%2$03u.data\" "
1802 "using 1:(($3 - $2) / clock_freq_%2$u):(($4 - $2) / clock_freq_%2$u) "
1803 "title \"Synchronization accuracy\" "
1804 "with filledcurves linewidth 2 linetype 1 "
1805 "linecolor rgb \"black\" fill solid 0.25 noborder, \\\n", i,
1806 j);
1807 }
1808 }
1809
1810
1811 /*
1812 * Write the analysis-specific graph lines in the gnuplot script.
1813 *
1814 * Args:
1815 * syncState: container for synchronization data
1816 * i: first trace number
1817 * j: second trace number, garanteed to be larger than i
1818 */
1819 static void writeAnalysisTraceTimeForePlotsCHull(SyncState* const syncState,
1820 const unsigned int i, const unsigned int j)
1821 {
1822 if (((AnalysisDataCHull*)
1823 syncState->analysisData)->graphsData->lpFactors->pairFactors[j][i].type
1824 == ACCURATE)
1825 {
1826 fprintf(syncState->graphsStream,
1827 "\t\"analysis_chull_accuracy-%1$03u_and_%2$03u.data\" "
1828 "using 1:(($3 - $2) / clock_freq_%2$u) notitle "
1829 "with lines linewidth 2 linetype 1 "
1830 "linecolor rgb \"gray60\", \\\n"
1831 "\t\"analysis_chull_accuracy-%1$03u_and_%2$03u.data\" "
1832 "using 1:(($4 - $2) / clock_freq_%2$u) notitle "
1833 "with lines linewidth 2 linetype 1 "
1834 "linecolor rgb \"gray60\", \\\n", i, j);
1835 }
1836 }
1837
1838
1839 /*
1840 * Write the analysis-specific graph lines in the gnuplot script.
1841 *
1842 * Args:
1843 * syncState: container for synchronization data
1844 * i: first trace number
1845 * j: second trace number, garanteed to be larger than i
1846 */
1847 static void writeAnalysisTraceTraceBackPlotsCHull(SyncState* const syncState,
1848 const unsigned int i, const unsigned int j)
1849 {
1850 if (((AnalysisDataCHull*)
1851 syncState->analysisData)->graphsData->lpFactors->pairFactors[j][i].type
1852 == ACCURATE)
1853 {
1854 fprintf(syncState->graphsStream,
1855 "\t\"analysis_chull_accuracy-%1$03u_and_%2$03u.data\" "
1856 "using 1:3:4 "
1857 "title \"Synchronization accuracy\" "
1858 "with filledcurves linewidth 2 linetype 1 "
1859 "linecolor rgb \"black\" fill solid 0.25 noborder, \\\n", i, j);
1860 }
1861 }
1862 #endif
LTTng 2.x Viewer (LTTV)
This page took 0.104661 seconds and 4 git commands to generate.