Index: soft/giet_vm/applications/ocean/main.C
===================================================================
--- soft/giet_vm/applications/ocean/main.C	(revision 589)
+++ soft/giet_vm/applications/ocean/main.C	(revision 598)
@@ -1,2 +1,696 @@
-#line 115 "/Users/alain/soc/giet_vm/applications/ocean/null_macros/c.m4.null.GIET"
-
+/*************************************************************************/
+/*                                                                       */
+/*  Copyright (c) 1994 Stanford University                               */
+/*                                                                       */
+/*  All rights reserved.                                                 */
+/*                                                                       */
+/*  Permission is given to use, copy, and modify this software for any   */
+/*  non-commercial purpose as long as this copyright notice is not       */
+/*  removed.  All other uses, including redistribution in whole or in    */
+/*  part, are forbidden without prior written permission.                */
+/*                                                                       */
+/*  This software is provided with absolutely no warranty and no         */
+/*  support.                                                             */
+/*                                                                       */
+/*************************************************************************/
+
+/*************************************************************************/
+/*                                                                       */
+/*  SPLASH Ocean Code                                                    */
+/*                                                                       */
+/*  This application studies the role of eddy and boundary currents in   */
+/*  influencing large-scale ocean movements.  This implementation uses   */
+/*  dynamically allocated four-dimensional arrays for grid data storage. */
+/*                                                                       */
+/*  Command line options:                                                */
+/*                                                                       */
+/*     -mM : Simulate MxM ocean. M must be (power of 2) +2.              */
+/*     -nN : N = number of threads. N must be power of 2.                */
+/*     -eE : E = error tolerance for iterative relaxation.               */
+/*     -rR : R = distance between grid points in meters.                 */
+/*     -tT : T = timestep in seconds.                                    */
+/*     -s  : Print timing statistics.                                    */
+/*     -o  : Print out relaxation residual values.                       */
+/*     -h  : Print out command line options.                             */
+/*                                                                       */
+/*  Default: OCEAN -m130 -n1 -e1e-7 -r20000.0 -t28800.0                  */
+/*                                                                       */
+/*  NOTE: This code works under both the FORK and SPROC models.          */
+/*                                                                       */
+/*************************************************************************/
+
+MAIN_ENV
+
+#define DEFAULT_M        514
+#define DEFAULT_N        4
+#define DEFAULT_E        1e-7
+#define DEFAULT_T    28800.0
+#define DEFAULT_R    20000.0
+#define UP               0
+#define DOWN             1
+#define LEFT             2
+#define RIGHT            3
+#define UPLEFT           4
+#define UPRIGHT          5
+#define DOWNLEFT         6
+#define DOWNRIGHT        7
+#define PAGE_SIZE     4096
+
+#include <stdio.h>
+#include <math.h>
+#include <stdlib.h>
+
+#include "decs.h"
+
+struct multi_struct *multi;
+struct global_struct *global;
+struct locks_struct *locks;
+struct bars_struct *bars;
+
+struct Global_Private *main_gp;
+double ****main_psi;
+double ****main_psim;
+double ***main_psium;
+double ***main_psilm;
+double ***main_psib;
+double ***main_ga;
+double ***main_gb;
+double ****main_work1;
+double ***main_work2;
+double ***main_work3;
+double ****main_work4;
+double ****main_work5;
+double ***main_work6;
+double ****main_work7;
+double ***main_oldga;
+double ***main_oldgb;
+double ****main_q_multi;
+double ****main_rhs_multi;
+double ****temparray;
+double ***tauz;
+long *main_imx;
+long *main_jmx;
+
+long nprocs = DEFAULT_N;
+const double h1 = 1000.0;
+const double h3 = 4000.0;
+const double h = 5000.0;
+const double lf = -5.12e11;
+double res = DEFAULT_R;
+double dtau = DEFAULT_T;
+const double f0 = 8.3e-5;
+const double beta = 2.0e-11;
+const double gpr = 0.02;
+double ysca;
+long oim;
+long jmm1;
+double tolerance = DEFAULT_E;
+const double pi = 3.141592653589793;
+const double t0 = 0.5e-4;
+const double outday0 = 1.0;
+const double outday1 = 2.0;
+const double outday2 = 2.0;
+const double outday3 = 2.0;
+const double maxwork = 10000.0;
+double factjacob;
+double factlap;
+
+//TODO : rÃ©pliquer Ã§a :
+double *main_lev_res;
+double *main_lev_tol;
+double *main_i_int_coeff;
+double *main_j_int_coeff;
+long *main_xpts_per_proc;
+long *main_ypts_per_proc;
+long main_xprocs;
+long main_yprocs;
+long main_numlev;
+double main_eig2;
+long main_im = DEFAULT_M;
+long main_jm;
+
+long minlevel;
+long do_stats = 1;
+long do_output = 0;
+long *ids_procs;
+
+
+__attribute__ ((constructor)) int main(int argc, char *argv[])
+{
+    long i;
+    long j;
+    long k;
+    long x_part;
+    long y_part;
+    long d_size;
+    long itemp;
+    long jtemp;
+    double procsqrt;
+    long temp = 0;
+    double min_total;
+    double max_total;
+    double avg_total;
+    double avg_wait;
+    double max_wait;
+    double min_wait;
+    double min_multi;
+    double max_multi;
+    double avg_multi;
+    double min_frac;
+    double max_frac;
+    double avg_frac;
+    long imax_wait;
+    long imin_wait;
+    long ch;
+    unsigned long long computeend;
+    unsigned long long start;
+    im = main_im;
+    
+    CLOCK(start);
+
+    while ((ch = getopt(argc, argv, "m:n:e:r:t:soh")) != -1) {
+        switch (ch) {
+        case 'm':
+            im = atoi(optarg);
+            if (log_2(im - 2) == -1) {
+                printerr("Grid must be ((power of 2)+2) in each dimension\n");
+                exit(-1);
+            }
+            break;
+        case 'n':
+            nprocs = atoi(optarg);
+            if (nprocs < 1) {
+                printerr("N must be >= 1\n");
+                exit(-1);
+            }
+            if (log_2(nprocs) == -1) {
+                printerr("N must be a power of 2\n");
+                exit(-1);
+            }
+            break;
+        case 'e':
+            tolerance = atof(optarg);
+            break;
+        case 'r':
+            res = atof(optarg);
+            break;
+        case 't':
+            dtau = atof(optarg);
+            break;
+        case 's':
+            do_stats = !do_stats;
+            break;
+        case 'o':
+            do_output = !do_output;
+            break;
+        case 'h':
+            printf("Usage: ocean <options>\n\n");
+            printf("options:\n");
+            printf("  -mM : Simulate MxM ocean.  M must be (power of 2) + 2 (default = %d).\n", DEFAULT_M);
+            printf("  -nN : N = number of threads. N must be power of 2 (default = %d).\n", DEFAULT_N);
+            printf("  -eE : E = error tolerance for iterative relaxation (default = %f).\n", DEFAULT_E);
+            printf("  -rR : R = distance between grid points in meters (default = %f).\n", DEFAULT_R);
+            printf("  -tT : T = timestep in seconds (default = %f).\n", DEFAULT_T);
+            printf("  -s  : Print timing statistics.\n");
+            printf("  -o  : Print out relaxation residual values.\n");
+            printf("  -h  : Print out command line options.\n\n");
+            exit(0);
+            break;
+        }
+    }
+
+    MAIN_INITENV
+    
+    jm = im;
+
+    printf("\n");
+    printf("Ocean simulation with W-cycle multigrid solver\n");
+    printf("    Processors                         : %1ld\n", nprocs);
+    printf("    Grid size                          : %1ld x %1ld\n", im, jm);
+    printf("    Grid resolution (meters)           : %0.2f\n", res);
+    printf("    Time between relaxations (seconds) : %0.0f\n", dtau);
+    printf("    Error tolerance                    : %0.7g\n", tolerance);
+    printf("\n");
+
+    xprocs = 0;
+    yprocs = 0;
+
+    procsqrt = sqrt((double) nprocs);
+    j = (long) procsqrt;
+
+    while ((xprocs == 0) && (j > 0)) {
+        k = nprocs / j;
+        if (k * j == nprocs) {
+            if (k > j) {
+                xprocs = j;
+                yprocs = k;
+            } else {
+                xprocs = k;
+                yprocs = j;
+            }
+        }
+        j--;
+    }
+
+    if (xprocs == 0) {
+        printerr("Could not find factors for subblocking\n");
+        exit(-1);
+    }
+
+    minlevel = 0;
+    itemp = 1;
+    jtemp = 1;
+    numlev = 0;
+    minlevel = 0;
+
+    while (itemp < (im - 2)) {
+        itemp = itemp * 2;
+        jtemp = jtemp * 2;
+        if ((itemp / yprocs > 1) && (jtemp / xprocs > 1)) {
+            numlev++;
+        }
+    }
+
+    if (numlev == 0) {
+        printerr("Must have at least 2 grid points per processor in each dimension\n");
+        exit(-1);
+    }
+
+    main_imx = (long *) G_MALLOC(numlev * sizeof(long), 0);
+    main_jmx = (long *) G_MALLOC(numlev * sizeof(long), 0);
+    main_lev_res = (double *) G_MALLOC(numlev * sizeof(double), 0);
+    main_lev_tol = (double *) G_MALLOC(numlev * sizeof(double), 0);
+    main_i_int_coeff = (double *) G_MALLOC(numlev * sizeof(double), 0);
+    main_j_int_coeff = (double *) G_MALLOC(numlev * sizeof(double), 0);
+    main_xpts_per_proc = (long *) G_MALLOC(numlev * sizeof(long), 0);
+    main_ypts_per_proc = (long *) G_MALLOC(numlev * sizeof(long), 0);
+    ids_procs = (long *) G_MALLOC(nprocs * sizeof(long), 0);
+    
+    imx = main_imx;
+    jmx = main_jmx;
+    lev_res = main_lev_res;
+    lev_tol = main_lev_tol;
+    i_int_coeff = main_i_int_coeff;
+    j_int_coeff = main_j_int_coeff;
+    xpts_per_proc = main_xpts_per_proc;
+    ypts_per_proc = main_ypts_per_proc;
+
+    for (i = 0; i < nprocs; i++) {
+        ids_procs[i] = i;
+    }
+
+    imx[numlev - 1] = im;
+    jmx[numlev - 1] = jm;
+    lev_res[numlev - 1] = res;
+    lev_tol[numlev - 1] = tolerance;
+
+    for (i = numlev - 2; i >= 0; i--) {
+        imx[i] = ((imx[i + 1] - 2) / 2) + 2;
+        jmx[i] = ((jmx[i + 1] - 2) / 2) + 2;
+        lev_res[i] = lev_res[i + 1] * 2;
+    }
+
+    for (i = 0; i < numlev; i++) {
+        xpts_per_proc[i] = (jmx[i] - 2) / xprocs;
+        ypts_per_proc[i] = (imx[i] - 2) / yprocs;
+    }
+    for (i = numlev - 1; i >= 0; i--) {
+        if ((xpts_per_proc[i] < 2) || (ypts_per_proc[i] < 2)) {
+            minlevel = i + 1;
+            break;
+        }
+    }
+
+    for (i = 0; i < numlev; i++) {
+        temp += imx[i];
+    }
+    temp = 0;
+    j = 0;
+    for (k = 0; k < numlev; k++) {
+        for (i = 0; i < imx[k]; i++) {
+            j++;
+            temp += jmx[k];
+        }
+    }
+
+    d_size = nprocs * sizeof(double ***);
+    main_psi = (double ****) G_MALLOC(d_size, 0);
+    main_psim = (double ****) G_MALLOC(d_size, 0);
+    main_work1 = (double ****) G_MALLOC(d_size, 0);
+    main_work4 = (double ****) G_MALLOC(d_size, 0);
+    main_work5 = (double ****) G_MALLOC(d_size, 0);
+    main_work7 = (double ****) G_MALLOC(d_size, 0);
+    temparray = (double ****) G_MALLOC(d_size, -1);
+
+    psi = main_psi;
+    psim = main_psim;
+    work1 = main_work1;
+    work4 = main_work4;
+    work5 = main_work5;
+    work7 = main_work7;
+
+    d_size = 2 * sizeof(double **);
+    for (i = 0; i < nprocs; i++) {
+        psi[i] = (double ***) G_MALLOC(d_size, i);
+        psim[i] = (double ***) G_MALLOC(d_size, i);
+        work1[i] = (double ***) G_MALLOC(d_size, i);
+        work4[i] = (double ***) G_MALLOC(d_size, i);
+        work5[i] = (double ***) G_MALLOC(d_size, i);
+        work7[i] = (double ***) G_MALLOC(d_size, i);
+        temparray[i] = (double ***) G_MALLOC(d_size, i);
+    }
+
+    d_size = nprocs * sizeof(double **);
+    main_psium = (double ***) G_MALLOC(d_size, 0);
+    main_psilm = (double ***) G_MALLOC(d_size, 0);
+    main_psib = (double ***) G_MALLOC(d_size, 0);
+    main_ga = (double ***) G_MALLOC(d_size, 0);
+    main_gb = (double ***) G_MALLOC(d_size, 0);
+    main_work2 = (double ***) G_MALLOC(d_size, 0);
+    main_work3 = (double ***) G_MALLOC(d_size, 0);
+    main_work6 = (double ***) G_MALLOC(d_size, 0);
+    tauz = (double ***) G_MALLOC(d_size, 0);
+    main_oldga = (double ***) G_MALLOC(d_size, 0);
+    main_oldgb = (double ***) G_MALLOC(d_size, 0);
+
+    psium = main_psium;
+    psilm = main_psilm;
+    psib = main_psib;
+    ga = main_ga;
+    gb = main_gb;
+    work2 = main_work2;
+    work3 = main_work3;
+    work6 = main_work6;
+    oldga = main_oldga;
+    oldgb = main_oldgb;
+
+    main_gp = (struct Global_Private *) G_MALLOC((nprocs + 1) * sizeof(struct Global_Private), -1);
+    gp = main_gp;
+
+    for (i = 0; i < nprocs; i++) {
+        gp[i].pad = (char *) G_MALLOC(PAGE_SIZE * sizeof(char), i);
+        gp[i].rel_num_x = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].rel_num_y = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].eist = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].ejst = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].oist = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].ojst = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].rlist = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].rljst = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].rlien = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].rljen = (long *) G_MALLOC(numlev * sizeof(long), i);
+        gp[i].neighbors = (long *) G_MALLOC(8 * sizeof(long), i);
+        gp[i].rownum = (long *) G_MALLOC(sizeof(long), i);
+        gp[i].colnum = (long *) G_MALLOC(sizeof(long), i);
+        gp[i].lpid = (long *) G_MALLOC(sizeof(long), i);
+        gp[i].multi_time = (double *) G_MALLOC(sizeof(double), i);
+        gp[i].total_time = (double *) G_MALLOC(sizeof(double), i);
+        gp[i].sync_time = (double *) G_MALLOC(sizeof(double), i);
+        gp[i].process_time = (double *) G_MALLOC(sizeof(double), i);
+        gp[i].step_start = (double *) G_MALLOC(sizeof(double), i);
+        gp[i].steps_time = (double *) G_MALLOC(10 * sizeof(double), i);
+        *gp[i].multi_time = 0;
+        *gp[i].total_time = 0;
+        *gp[i].sync_time = 0;
+        *gp[i].process_time = 0;
+        *gp[i].lpid = i;
+    }
+
+    subblock();
+
+    x_part = (jm - 2) / xprocs + 2;
+    y_part = (im - 2) / yprocs + 2;
+
+    d_size = x_part * y_part * sizeof(double) + y_part * sizeof(double *);
+
+    global = (struct global_struct *) G_MALLOC(sizeof(struct global_struct), -1);
+
+    for (i = 0; i < nprocs; i++) {
+        psi[i][0] = (double **) G_MALLOC(d_size, i);
+        psi[i][1] = (double **) G_MALLOC(d_size, i);
+        psim[i][0] = (double **) G_MALLOC(d_size, i);
+        psim[i][1] = (double **) G_MALLOC(d_size, i);
+        psium[i] = (double **) G_MALLOC(d_size, i);
+        psilm[i] = (double **) G_MALLOC(d_size, i);
+        psib[i] = (double **) G_MALLOC(d_size, i);
+        ga[i] = (double **) G_MALLOC(d_size, i);
+        gb[i] = (double **) G_MALLOC(d_size, i);
+        work1[i][0] = (double **) G_MALLOC(d_size, i);
+        work1[i][1] = (double **) G_MALLOC(d_size, i);
+        work2[i] = (double **) G_MALLOC(d_size, i);
+        work3[i] = (double **) G_MALLOC(d_size, i);
+        work4[i][0] = (double **) G_MALLOC(d_size, i);
+        work4[i][1] = (double **) G_MALLOC(d_size, i);
+        work5[i][0] = (double **) G_MALLOC(d_size, i);
+        work5[i][1] = (double **) G_MALLOC(d_size, i);
+        work6[i] = (double **) G_MALLOC(d_size, i);
+        work7[i][0] = (double **) G_MALLOC(d_size, i);
+        work7[i][1] = (double **) G_MALLOC(d_size, i);
+        temparray[i][0] = (double **) G_MALLOC(d_size, i);
+        temparray[i][1] = (double **) G_MALLOC(d_size, i);
+        tauz[i] = (double **) G_MALLOC(d_size, i);
+        oldga[i] = (double **) G_MALLOC(d_size, i);
+        oldgb[i] = (double **) G_MALLOC(d_size, i);
+    }
+
+    oim = im;
+    //f = (double *) G_MALLOC(oim*sizeof(double), 0);
+    multi = (struct multi_struct *) G_MALLOC(sizeof(struct multi_struct), -1);
+
+    d_size = numlev * sizeof(double **);
+    if (numlev % 2 == 1) {      /* To make sure that the actual data
+                                   starts double word aligned, add an extra
+                                   pointer */
+        d_size += sizeof(double **);
+    }
+    for (i = 0; i < numlev; i++) {
+        d_size += ((imx[i] - 2) / yprocs + 2) * ((jmx[i] - 2) / xprocs + 2) * sizeof(double) + ((imx[i] - 2) / yprocs + 2) * sizeof(double *);
+    }
+
+    d_size *= nprocs;
+
+    if (nprocs % 2 == 1) {      /* To make sure that the actual data
+                                   starts double word aligned, add an extra
+                                   pointer */
+        d_size += sizeof(double ***);
+    }
+
+    d_size += nprocs * sizeof(double ***);
+    main_q_multi = (double ****) G_MALLOC(d_size, -1);
+    main_rhs_multi = (double ****) G_MALLOC(d_size, -1);
+    q_multi = main_q_multi;
+    rhs_multi = main_rhs_multi;
+
+
+    locks = (struct locks_struct *) G_MALLOC(sizeof(struct locks_struct), -1);
+    bars = (struct bars_struct *) G_MALLOC(sizeof(struct bars_struct), -1);
+
+    LOCKINIT(locks->idlock)
+    LOCKINIT(locks->psiailock)
+    LOCKINIT(locks->psibilock)
+    LOCKINIT(locks->donelock)
+    LOCKINIT(locks->error_lock)
+    LOCKINIT(locks->bar_lock)
+#if defined(MULTIPLE_BARRIERS)
+    BARINIT(bars->iteration, nprocs)
+    BARINIT(bars->gsudn, nprocs)
+    BARINIT(bars->p_setup, nprocs)
+    BARINIT(bars->p_redph, nprocs)
+    BARINIT(bars->p_soln, nprocs)
+    BARINIT(bars->p_subph, nprocs)
+    BARINIT(bars->sl_prini, nprocs)
+    BARINIT(bars->sl_psini, nprocs)
+    BARINIT(bars->sl_onetime, nprocs)
+    BARINIT(bars->sl_phase_1, nprocs)
+    BARINIT(bars->sl_phase_2, nprocs)
+    BARINIT(bars->sl_phase_3, nprocs)
+    BARINIT(bars->sl_phase_4, nprocs)
+    BARINIT(bars->sl_phase_5, nprocs)
+    BARINIT(bars->sl_phase_6, nprocs)
+    BARINIT(bars->sl_phase_7, nprocs)
+    BARINIT(bars->sl_phase_8, nprocs)
+    BARINIT(bars->sl_phase_9, nprocs)
+    BARINIT(bars->sl_phase_10, nprocs)
+    BARINIT(bars->error_barrier, nprocs)
+#else
+    BARINIT(bars->barrier, nprocs)
+#endif
+    link_all();
+
+    multi->err_multi = 0.0;
+    i_int_coeff[0] = 0.0;
+    j_int_coeff[0] = 0.0;
+
+    for (i = 0; i < numlev; i++) {
+        i_int_coeff[i] = 1.0 / (imx[i] - 1);
+        j_int_coeff[i] = 1.0 / (jmx[i] - 1);
+    }
+
+    /*
+       initialize constants and variables
+
+       id is a global shared variable that has fetch-and-add operations
+       performed on it by processes to obtain their pids.   
+     */
+
+    //global->id = 0;
+    global->trackstart = 0;
+    global->psibi = 0.0;
+
+    factjacob = -1. / (12. * res * res);
+    factlap = 1. / (res * res);
+    eig2 = -h * f0 * f0 / (h1 * h3 * gpr);
+
+    jmm1 = jm - 1;
+    ysca = ((double) jmm1) * res;
+    im = (imx[numlev - 1] - 2) / yprocs + 2;
+    jm = (jmx[numlev - 1] - 2) / xprocs + 2;
+    
+    main_im = im;
+    main_jm = jm;
+    main_numlev = numlev;
+    main_xprocs = xprocs;
+    main_yprocs = yprocs;
+    main_eig2 = eig2;
+
+    if (do_output) {
+        printf("              MULTIGRID OUTPUTS\n");
+    }
+
+    CREATE(slave, nprocs);
+    WAIT_FOR_END(nprocs);
+    CLOCK(computeend);
+
+    printf("\n");
+    printf("                PROCESS STATISTICS\n");
+    printf("                  Total          Multigrid         Multigrid\n");
+    printf(" Proc             Time             Time            Fraction\n");
+    printf("    0   %15.0f    %15.0f        %10.3f\n", (*gp[0].total_time), (*gp[0].multi_time), (*gp[0].multi_time) / (*gp[0].total_time));
+
+    if (do_stats) {
+        double phase_time;
+        min_total = max_total = avg_total = (*gp[0].total_time);
+        min_multi = max_multi = avg_multi = (*gp[0].multi_time);
+        min_frac = max_frac = avg_frac = (*gp[0].multi_time) / (*gp[0].total_time);
+        avg_wait = *gp[0].sync_time;
+        max_wait = *gp[0].sync_time;
+        min_wait = *gp[0].sync_time;
+        imax_wait = 0;
+        imin_wait = 0;
+
+        for (i = 1; i < nprocs; i++) {
+            if ((*gp[i].total_time) > max_total) {
+                max_total = (*gp[i].total_time);
+            }
+            if ((*gp[i].total_time) < min_total) {
+                min_total = (*gp[i].total_time);
+            }
+            if ((*gp[i].multi_time) > max_multi) {
+                max_multi = (*gp[i].multi_time);
+            }
+            if ((*gp[i].multi_time) < min_multi) {
+                min_multi = (*gp[i].multi_time);
+            }
+            if ((*gp[i].multi_time) / (*gp[i].total_time) > max_frac) {
+                max_frac = (*gp[i].multi_time) / (*gp[i].total_time);
+            }
+            if ((*gp[i].multi_time) / (*gp[i].total_time) < min_frac) {
+                min_frac = (*gp[i].multi_time) / (*gp[i].total_time);
+            }
+            avg_total += (*gp[i].total_time);
+            avg_multi += (*gp[i].multi_time);
+            avg_frac += (*gp[i].multi_time) / (*gp[i].total_time);
+            avg_wait += (*gp[i].sync_time);
+            if (max_wait < (*gp[i].sync_time)) {
+                max_wait = (*gp[i].sync_time);
+                imax_wait = i;
+            }
+            if (min_wait > (*gp[i].sync_time)) {
+                min_wait = (*gp[i].sync_time);
+                imin_wait = i;
+            }
+        }
+        avg_total = avg_total / nprocs;
+        avg_multi = avg_multi / nprocs;
+        avg_frac = avg_frac / nprocs;
+        avg_wait = avg_wait / nprocs;
+        for (i = 1; i < nprocs; i++) {
+            printf("  %3ld   %15.0f    %15.0f        %10.3f\n", i, (*gp[i].total_time), (*gp[i].multi_time), (*gp[i].multi_time) / (*gp[i].total_time));
+        }
+        printf("  Avg   %15.0f    %15.0f        %10.3f\n", avg_total, avg_multi, avg_frac);
+        printf("  Min   %15.0f    %15.0f        %10.3f\n", min_total, min_multi, min_frac);
+        printf("  Max   %15.0f    %15.0f        %10.3f\n", max_total, max_multi, max_frac);
+        
+        printf("\n\n                  Sync\n");
+        printf(" Proc      Time        Fraction\n");
+        for (i = 0; i < nprocs; i++) {
+            printf("  %ld        %u      %f\n", i, (unsigned int)*gp[i].sync_time, *gp[i].sync_time / ((long)(*gp[i].total_time)));
+        }
+
+        printf("  Avg   %f   %f\n", avg_wait, (double) avg_wait / (long) (computeend - global->trackstart));
+        printf("  Min   %f   %f\n", min_wait, (double) min_wait / (long) (*gp[imin_wait].total_time));
+        printf("  Max   %f   %f\n", max_wait, (double) max_wait / (long) (*gp[imax_wait].total_time));
+
+        printf("\nPhases Avg :\n\n");
+        for (i = 0; i < 10; i++) {
+            phase_time = 0;
+            for (j = 0; j < nprocs; j++) {
+                phase_time += gp[j].steps_time[i];
+            }
+            phase_time /= (double) nprocs;
+            printf("  %d = %f (fraction %f)\n", i + 1, phase_time, phase_time / (long) (computeend - global->trackstart));
+        }
+    }
+    printf("\n");
+
+    global->starttime = start;
+    printf("                       TIMING INFORMATION\n");
+    printf("[NPROCS]           : %16ld\n", nprocs);
+    printf("[START1]           : %16llu\n", global->starttime);
+    printf("[START2]           : %16llu\n", global->trackstart);
+    printf("[END]              : %16llu\n", computeend);
+    printf("[TOTAL]            : %16llu\n", computeend - global->starttime);    // With init
+    printf("[PARALLEL_COMPUTE] : %16llu\n", computeend - global->trackstart);   // Without init
+    printf("(excludes first timestep)\n");
+    printf("\n");
+
+    MAIN_END
+    
+}
+
+long log_2(long number)
+{
+    long cumulative = 1;
+    long out = 0;
+    long done = 0;
+
+    while ((cumulative < number) && (!done) && (out < 50)) {
+        if (cumulative == number) {
+            done = 1;
+        } else {
+            cumulative = cumulative * 2;
+            out++;
+        }
+    }
+
+    if (cumulative == number) {
+        return (out);
+    } else {
+        return (-1);
+    }
+}
+
+void printerr(char *s)
+{
+    fprintf(stderr, "ERROR: %s\n", s);
+}
+
+
+// Local Variables:
+// tab-width: 4
+// c-basic-offset: 4
+// c-file-offsets:((innamespace . 0)(inline-open . 0))
+// indent-tabs-mode: nil
+// End:
+
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=4:softtabstop=4
