//removes the lowest energy vertical seam from the image
void removeVerticalSeam() {
    double energies[3];
    double min_energy;
    int prev_x;
    int prev_y;
    
    // split up work between processes
    double *my_path_costs;
    double *my_previous_x;
    double *my_previous_y;
    double *temp_path_costs;
    double *temp_previous_x;
    double *temp_previous_y;
    int my_cols = current_width / numprocs;
    int low_cols = my_cols;
    int extra_cols = current_width % numprocs;
    int start;
    int x_offset;
    int recv_cols;

    double left_end_cost, right_end_cost, temp_end_cost;

    if (rank < extra_cols) {
        my_cols++;
        start = rank * my_cols;
    } else {
        start = (extra_cols * (my_cols + 1)) + ((rank - extra_cols) * my_cols);
    }

    //printf("%d %d %d\n", rank, start, my_cols);
    
    my_path_costs = (double *) malloc(my_cols * current_height * sizeof(double));
    my_previous_x = (double *) malloc(my_cols * current_height * sizeof(double));
    my_previous_y = (double *) malloc(my_cols * current_height * sizeof(double));
    //find the lowest cost seam by computing the lowest cost paths to each pixel
    for (int y = 0; y < current_height; y++) {
        //compute the path costs for my columns     
        for (int x = start; x < start + my_cols; x++) {
            //printf("%d %d %d %d %d\n", rank, x, y, (x - start) * current_height + y, my_cols * current_height);
            if (y == 0) {
                path_costs[x * initial_height] = image_energy[x * initial_height];
                my_path_costs[(x - start) * current_height + y] = path_costs[x * initial_height];

                previous_x[x * initial_height] = -1;
                my_previous_x[(x - start) * current_height + y] = previous_x[x * initial_height];

                previous_y[x * initial_height] = -1;
                my_previous_y[(x - start) * current_height + y] = previous_y[x * initial_height];
            } else {
                //the pixel directly above
                energies[1] = path_costs[x * initial_height + y - 1];
                //pixel above to the left
                if (x != 0) {
                    energies[0] = path_costs[(x - 1) * initial_height + y - 1];
                } else {
                    energies[0] = DBL_MAX;
                }
                //pixel above to the right
                if (x != current_width - 1) {
                    energies[2] = path_costs[(x + 1) * initial_height + y - 1];
                } else {
                    energies[2] = DBL_MAX;
                }

                //find the one with the least path cost
                min_energy = energies[0];
                prev_x = x - 1;
                prev_y = y - 1;
                if (energies[1] < min_energy) {
                    min_energy = energies[1];
                    prev_x = x;
                }
                if (energies[2] < min_energy) {
                    min_energy = energies[2];
                    prev_x = x + 1;
                }

                //set the minimum path cost for this pixel
                path_costs[x * initial_height + y] = min_energy + image_energy[x * initial_height + y];
                my_path_costs[(x - start) * current_height + y] = path_costs[x * initial_height + y];

                //set the previous pixel on the minimum path's coordinates for this pixel
                previous_x[x * initial_height + y] = prev_x;
                my_previous_x[(x - start) * current_height + y] = previous_x[x * initial_height + y];

                previous_y[x * initial_height + y] = prev_y;
                my_previous_y[(x - start) * current_height + y] = previous_y[x * initial_height + y];
            }
        }

        //send path cost needed to neighboring processes
        if (numprocs > 1) {
            if (rank != numprocs - 1) {
                //send rightmost cost to following process
                right_end_cost = path_costs[(start + my_cols - 1) * initial_height + y];        
                MPI_Send(&right_end_cost, 1, MPI_DOUBLE, rank + 1, 0, MPI_COMM_WORLD);

                //receive following process's leftmost cost
                MPI_Recv(&temp_end_cost, 1, MPI_DOUBLE, rank + 1, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
//.........这里部分代码省略.........

int main(int argc, char** argv)
{
	MPI_Init(&argc, &argv);
	int size, rank;
	MPI_Comm_size(MPI_COMM_WORLD, &size);
	MPI_Comm_rank(MPI_COMM_WORLD, &rank);

	printf("I am process number - %d\n", rank);

	int rowsA = 9;
	int columnA = 4;
	int columnB = 9;
	int rowsB = columnA;
	TMatrix matrixC = createMatrix(rowsA, columnB);

	if(rank == 0)
	{
		double start, end;
		TMatrix matrixA = createMatrix(rowsA, columnA);
		TMatrix matrixB = createMatrix(rowsB, columnB);

		fillSimpleMatrix(matrixA.data, matrixA.rows, matrixA.columns);
		fillSimpleMatrix(matrixB.data, matrixB.rows, matrixB.columns);
		
		printMatrix(matrixA);
		printf("\n");
		printMatrix(matrixB);
		printf("\n");

		int numberOfSlaves = size - 1;
		start = MPI_Wtime();

		int rowsPerWorker = rowsA / numberOfSlaves;
		int remainingRows = rowsA % numberOfSlaves;
		int offsetRow = 0;
		int messageType = FROM_MASTER;
		
		for(int destination = 1; destination <= numberOfSlaves; destination++)
		{
			int rows = (destination <= remainingRows) ? rowsPerWorker + 1 : rowsPerWorker;
			MPI_Send((void *)&offsetRow, 1, MPI_INT, destination, messageType, MPI_COMM_WORLD);
			MPI_Send((void *)&rows, 1, MPI_INT, destination, messageType, MPI_COMM_WORLD);

			double* temp = matrixA.data + offsetRow;


			MPI_Send((void *)temp, rows * columnA, MPI_DOUBLE, destination, messageType, MPI_COMM_WORLD);
			MPI_Send((void *)matrixB.data, rowsB * columnB, MPI_DOUBLE, destination, messageType, MPI_COMM_WORLD);
			offsetRow += rows;
		}

		messageType = FROM_SLAVE;

		for(int source = 1; source <= numberOfSlaves; source++)
		{
			int rowOffset;
			MPI_Status status;

			MPI_Recv((void*)&rowOffset, 1, MPI_INT, source, messageType, MPI_COMM_WORLD, &status);
			
			int rows;
			MPI_Recv((void*)&rows, 1, MPI_INT, source, messageType, MPI_COMM_WORLD, &status);

			double* temp = (double*)malloc(rows * columnB * sizeof(double));

			MPI_Recv((void*)temp, rows * columnB, MPI_DOUBLE, source, FROM_SLAVE, MPI_COMM_WORLD, &status);

			for(int j = 0; j < rows * columnB; j++)
			{
				matrixC.data[rowOffset * columnB + j] = temp[j];
			}
		}

		printMatrix(matrixC);
		printf("\n");
		end = MPI_Wtime();
		printf("time: %.4f\n", end - start);
	}
	else
	{
		int offsetRow;
		MPI_Status status;
		MPI_Recv((void*)&offsetRow, 1, MPI_INT, 0, FROM_MASTER, MPI_COMM_WORLD, &status);

		int rows;
		MPI_Recv((void*)&rows, 1, MPI_INT, 0, FROM_MASTER, MPI_COMM_WORLD, &status);
		printf("Process - %d, Offset - %d, Rows - %d\n", rank, offsetRow, rows);

		TMatrix aMatrix = createMatrix(rows, columnA);
		MPI_Recv((void*)aMatrix.data, rows * columnA, MPI_DOUBLE, 0, FROM_MASTER, MPI_COMM_WORLD, &status);


		TMatrix bMatrix = createMatrix(rowsB, columnB);
		MPI_Recv((void*)bMatrix.data, rowsB * columnB, MPI_DOUBLE, 0, FROM_MASTER, MPI_COMM_WORLD, &status);

		double* c =(double*)malloc(rows * columnB * sizeof(double));

		for(int i = 0; i < rows; i++)
		{
			for(int j = 0; j < columnB; j++)
//.........这里部分代码省略.........

int main(int argc, char *argv[]) 
{
   //precision to work to
   sscanf(argv[1], "%lf", &PRECISION);

   fill_array();

   long int sum, partial_sum;
   MPI_Status status;
   int ID, root_process, ierr, i, num_procs, an_id, sender, length, flag;

   ierr = MPI_Init(&argc, &argv);
   
   root_process = 0;
   
   /* find out MY process ID, and how many processes were started. */
   
   ierr = MPI_Comm_rank(MPI_COMM_WORLD, &ID);
   ierr = MPI_Comm_size(MPI_COMM_WORLD, &num_procs);

   THREAD_NUM = num_procs-1;

   int INNER_ARRAY_SIZE = n-2;
   int TOTAL_NUMS = (INNER_ARRAY_SIZE)*(INNER_ARRAY_SIZE);
      
   //splits total amount of numbers in matrix into number of threads beng used, and finds the remainder
   int NUMS_PER_THREAD = (int)floor(TOTAL_NUMS / THREAD_NUM);
   int NUMS_REMAINDER = (TOTAL_NUMS % THREAD_NUM);

   if(ID == root_process) {

      // ROOT PROCESS

      //copies original array into a temporary array
      memcpy(temparr, arr, sizeof(arr));

      matrix_work(num_procs, INNER_ARRAY_SIZE, NUMS_PER_THREAD, NUMS_REMAINDER, TOTAL_NUMS);

      while(get_diff()!=0){
         memcpy(arr, temparr, sizeof(arr));
         matrix_work(num_procs, INNER_ARRAY_SIZE, NUMS_PER_THREAD, NUMS_REMAINDER, TOTAL_NUMS);
      }


      flag=1;
      for(an_id = 1; an_id < num_procs; an_id++) {
         ierr = MPI_Send(&flag, 1, MPI_INT, an_id, send_data_tag, MPI_COMM_WORLD);
      }

   }
   else {

      // SLAVE PROCESS

      //while program not done, calculate averages from range in array and send back to root process
      while(get_flag()==0){
         ierr = MPI_Recv(&(arr[0][0]), n*n, MPI_DOUBLE, root_process, send_data_tag, MPI_COMM_WORLD, &status);
      

         //Calculates the start and end numbers in the matrix for the thread to work on
         int START_NUM = 0 + ((ID-1)*NUMS_PER_THREAD);
         int END_NUM = (ID*NUMS_PER_THREAD)-1;

         if (ID == THREAD_NUM){
            END_NUM = END_NUM+ NUMS_REMAINDER;
         }
         
         int range = END_NUM - START_NUM;
         double localnums[range];

         //loops through numbers from start to end that need to be worked on and overwrites the temp array
         int i;
         int count = 0;
         for(i = START_NUM; i <= END_NUM; i++) {
            //get the x & y values for corresponding number in the matrix
            int x = (i%INNER_ARRAY_SIZE)+1;
            int y = ((i-(i%INNER_ARRAY_SIZE))/INNER_ARRAY_SIZE)+1;
            //get four surrounding numbers and average them
            //above
            double a = arr[y-1][x];
            //right
            double b = arr[y][x+1];
            //below
            double c = arr[y+1][x];
            //left
            double d = arr[y][x-1];

            double average = (a+b+c+d)/4;

            localnums[count] = average;
            count++;
         }

         //and finally, send array of new numbers back to the root process
         ierr = MPI_Send(&count, 1, MPI_INT, root_process, return_data_tag, MPI_COMM_WORLD);
         ierr = MPI_Send(&localnums, count, MPI_DOUBLE, root_process, return_data_tag, MPI_COMM_WORLD);
      }
   }

   ierr = MPI_Finalize();
//.........这里部分代码省略.........

double timeStepper::computeDt(int &numReads, int &numWrites)
{
  // Time step control
  array minSpeedTemp,maxSpeedTemp;
  array minSpeed,maxSpeed;
  elemOld->computeMinMaxCharSpeeds(directions::X1,
                                   minSpeedTemp, maxSpeedTemp,
                                   numReads,numWrites
                                  );
  minSpeedTemp = minSpeedTemp/XCoords->dX1;
  maxSpeedTemp = maxSpeedTemp/XCoords->dX1;
  maxSpeed     = af::max(maxSpeedTemp,af::abs(minSpeedTemp));

  if(params::dim>1)
  {
    elemOld->computeMinMaxCharSpeeds(directions::X2,
                                     minSpeedTemp, maxSpeedTemp,
                                     numReads,numWrites
                                    );
    minSpeedTemp = minSpeedTemp/XCoords->dX2;
    maxSpeedTemp = maxSpeedTemp/XCoords->dX2;
    maxSpeed    += af::max(maxSpeedTemp,af::abs(minSpeedTemp));
  }

  if(params::dim>2)
  {
    elemOld->computeMinMaxCharSpeeds(directions::X3,
                                     minSpeedTemp, maxSpeedTemp,
                                     numReads,numWrites);
    minSpeedTemp = minSpeedTemp/XCoords->dX3;
    maxSpeedTemp = maxSpeedTemp/XCoords->dX3;
    maxSpeed    += af::max(maxSpeedTemp,af::abs(minSpeedTemp));
  }
  array maxInvDt_af = af::max(af::max(af::max(maxSpeed,2),1),0);
  double maxInvDt = maxInvDt_af.host<double>()[0];

  /* Use MPI to find minimum over all processors */
  if (world_rank == 0) 
  {
    double temp; 
    for(int i=1;i<world_size;i++)
    {
      MPI_Recv(&temp, 1, MPI_DOUBLE, i, i, PETSC_COMM_WORLD,MPI_STATUS_IGNORE);
      if( maxInvDt < temp)
      {
        maxInvDt = temp;
      }
    }
  }
  else
  {
    MPI_Send(&maxInvDt, 1, MPI_DOUBLE, 0, world_rank, PETSC_COMM_WORLD);
  }
  MPI_Barrier(PETSC_COMM_WORLD);
  MPI_Bcast(&maxInvDt,1,MPI_DOUBLE,0,PETSC_COMM_WORLD);
  MPI_Barrier(PETSC_COMM_WORLD);
  
  double newDt = params::CourantFactor/maxInvDt;
    
  if (newDt > params::maxDtIncrement*dt)
  {
    newDt = params::maxDtIncrement*dt;
  }
  dt = newDt;
}

int MTestTestIntercomm(MPI_Comm comm)
{
    int local_size, remote_size, rank, **bufs, *bufmem, rbuf[2], j;
    int errs = 0, wrank, nsize;
    char commname[MPI_MAX_OBJECT_NAME + 1];
    MPI_Request *reqs;

    MPI_Comm_rank(MPI_COMM_WORLD, &wrank);
    MPI_Comm_size(comm, &local_size);
    MPI_Comm_remote_size(comm, &remote_size);
    MPI_Comm_rank(comm, &rank);
    MPI_Comm_get_name(comm, commname, &nsize);

    MTestPrintfMsg(1, "Testing communication on intercomm '%s', remote_size=%d\n",
                   commname, remote_size);

    reqs = (MPI_Request *) malloc(remote_size * sizeof(MPI_Request));
    if (!reqs) {
        printf("[%d] Unable to allocated %d requests for testing intercomm %s\n",
               wrank, remote_size, commname);
        errs++;
        return errs;
    }
    bufs = (int **) malloc(remote_size * sizeof(int *));
    if (!bufs) {
        printf("[%d] Unable to allocated %d int pointers for testing intercomm %s\n",
               wrank, remote_size, commname);
        errs++;
        return errs;
    }
    bufmem = (int *) malloc(remote_size * 2 * sizeof(int));
    if (!bufmem) {
        printf("[%d] Unable to allocated %d int data for testing intercomm %s\n",
               wrank, 2 * remote_size, commname);
        errs++;
        return errs;
    }

    /* Each process sends a message containing its own rank and the
     * rank of the destination with a nonblocking send.  Because we're using
     * nonblocking sends, we need to use different buffers for each isend */
    /* NOTE: the send buffer access restriction was relaxed in MPI-2.2, although
     * it doesn't really hurt to keep separate buffers for our purposes */
    for (j = 0; j < remote_size; j++) {
        bufs[j] = &bufmem[2 * j];
        bufs[j][0] = rank;
        bufs[j][1] = j;
        MPI_Isend(bufs[j], 2, MPI_INT, j, 0, comm, &reqs[j]);
    }
    MTestPrintfMsg(2, "isends posted, about to recv\n");

    for (j = 0; j < remote_size; j++) {
        MPI_Recv(rbuf, 2, MPI_INT, j, 0, comm, MPI_STATUS_IGNORE);
        if (rbuf[0] != j) {
            printf("[%d] Expected rank %d but saw %d in %s\n", wrank, j, rbuf[0], commname);
            errs++;
        }
        if (rbuf[1] != rank) {
            printf("[%d] Expected target rank %d but saw %d from %d in %s\n",
                   wrank, rank, rbuf[1], j, commname);
            errs++;
        }
    }
    if (errs)
        fflush(stdout);

    MTestPrintfMsg(2, "my recvs completed, about to waitall\n");
    MPI_Waitall(remote_size, reqs, MPI_STATUSES_IGNORE);

    free(reqs);
    free(bufs);
    free(bufmem);

    return errs;
}

void sync_display(void)
{
#ifdef CONF_MPI
    MPI_Status stat;
    int        size = 0;
#endif
    int        rank = 0;

    char name[MAXNAME];
    int  i = 0;

    gethostname(name, MAXNAME);

#ifdef CONF_MPI

    assert_mpi(MPI_Comm_rank(MPI_COMM_WORLD, &rank));
    assert_mpi(MPI_Comm_size(MPI_COMM_WORLD, &size));

    if (rank)
    {
        /* Send the name to the root, recieve the host index. */

        MPI_Send(name, MAXNAME, MPI_BYTE, 0,0, MPI_COMM_WORLD);
        MPI_Recv(&i,   4,       MPI_BYTE, 0,0, MPI_COMM_WORLD, &stat);
    }
    else
    {
        int j, k;

        /* Find a host definition for the root. */

        i = find_display(name);

        /* Recieve a name from each client, send a host definition index. */

        for (j = 1; j < size; ++j)
        {
            MPI_Recv(name, MAXNAME, MPI_BYTE, j,0, MPI_COMM_WORLD, &stat);
            k = find_display(name);
            MPI_Send(&k,   4,       MPI_BYTE, j,0, MPI_COMM_WORLD);
        }
    }
#else

    i = find_display(name);

#endif

    /* If no host definition was found, create a default. */

    if (i == 0)
    {
        i = add_host(DEFAULT_NAME, DEFAULT_X, DEFAULT_Y, DEFAULT_W, DEFAULT_H);
            add_tile(i,            DEFAULT_X, DEFAULT_Y, DEFAULT_W, DEFAULT_H);

        host[i].flags = HOST_FRAMED;
    }

    /* Note the indexed host definition as current. */

    current_host = host[i];

    /* Position the server window, if necessary. */

    if (rank || (current_host.flags & HOST_FRAMED) == 0)
        set_window_pos(current_host.win_x, current_host.win_y);
}

unsigned int master(unsigned int base_dim, unsigned int max_fact, 
		    unsigned int** exponents, mpz_t * As,
		    int comm_size, unsigned int print_fact) {

  unsigned int fact_count = 0;

  MPI_Status status;

  int count;
  int source;
  
  /* Buffer per ricevere gli esponenti */
  unsigned int* buffer_exp;
  /* Buffer per ricevere (A + s) */
  unsigned char buffer_As[BUFFER_DIM];
  init_vector(& buffer_exp, base_dim);

  double t1 = MPI_Wtime();
  double t2;

  int fact_per_rank[comm_size];
  for(int i = 0; i < comm_size; ++i)
    fact_per_rank[i] = 0;

  while(fact_count < max_fact + base_dim) {
    /* Ricevo il vettore di esponenti */
    MPI_Recv(buffer_exp, base_dim, MPI_UNSIGNED,
	     MPI_ANY_SOURCE, ROW_TAG, 
	     MPI_COMM_WORLD, &status);
    source = status.MPI_SOURCE;
    
    for(unsigned int i = 0; i < base_dim; ++i) 
      set_matrix(exponents, fact_count, i, buffer_exp[i]);
    
    /* Ricevo l'mpz contenente (A + s) */
    MPI_Recv(buffer_As, BUFFER_DIM, MPI_UNSIGNED_CHAR, source, 
	     AS_TAG, MPI_COMM_WORLD, &status);
    MPI_Get_count(&status, MPI_UNSIGNED_CHAR, &count);
    mpz_import(As[fact_count], count, 1, 1, 1, 0, buffer_As);
    
    ++fact_count;
    ++fact_per_rank[source];

    if(fact_count % print_fact == 0) {
      t2 = MPI_Wtime() - t1;
      
      printf("#%d/%d in %.6f seconds\n", fact_count, max_fact + base_dim, t2);
    }
  }
  
  /* Spedisco '1' agli slave per indicare la terminazione */
  char stop_signal = '1';
  for(unsigned int i = 1; i < comm_size; ++i)
    MPI_Send(&stop_signal, 1, MPI_CHAR, i, 0, MPI_COMM_WORLD);
  
  printf("#Sending stop_signal\n");

  printf("#Fattorizzazioni per ranks:\n#");
  for(int i = 1; i < comm_size; ++i)
    printf("%d \t", i);
  printf("\n#");
  for(int i = 1; i < comm_size; ++i)
    printf("%d \t", fact_per_rank[i]);
  printf("\n");
 
  return fact_count;
}

    int TransmitProteinsToChildProcesses()
    {
        int numProteins = (int) proteins.size();

        vector< simplethread_handle_t > workerHandles;

        int sourceProcess, batchSize;
        bool IsFinished = false;

        Timer searchTime( true );
        float totalSearchTime = 0.01f;
        float lastUpdate = 0.0f;

        int i = 0;
        int numChildrenFinished = 0;
        while( numChildrenFinished < g_numChildren )
        {
            #ifdef MPI_DEBUG
                cout << g_hostString << " is listening for a child process to offer to search some proteins." << endl;
            #endif
            // Listen for a process requesting proteins. 
            // Extract the number of CPUs available on the process.
            MPI_Recv( &sourceProcess,            1,        MPI_INT,    MPI_ANY_SOURCE,    0xFF, MPI_COMM_WORLD, &st );
            int sourceCPUs = 0;
            MPI_Recv( &sourceCPUs,                1,        MPI_INT,    sourceProcess,    0xFF, MPI_COMM_WORLD, &st );
            
            int pOffset = i;
            // Scale the batchSize with the number of cpus in the requested process.
            batchSize = min( numProteins-i, g_rtConfig->ProteinBatchSize*sourceCPUs );

            stringstream packStream;
            binary_oarchive packArchive( packStream );

            try
            {
                packArchive & pOffset;
                string proteinStream;
                for( int j = i; j < i + batchSize; ++j )
                {
                    proteinStream += ">" + proteins[j].getName() + " " + proteins[j].getDescription() + "\n" 
                                  + proteins[j].getSequence() + "\n";
                }
                packArchive & proteinStream;
            } catch( exception& e )
            {
                cerr << g_hostString << " had an error: " << e.what() << endl;
                exit(1);
            }
            #ifdef MPI_DEBUG
                cout << "Process #" << sourceProcess << " has " << sourceCPUs << " cpus. Sending " << batchSize << " proteins." << endl;
            #endif
            
            if( i < numProteins )
            {
                MPI_Ssend( &batchSize,            1,        MPI_INT,    sourceProcess,    0x99, MPI_COMM_WORLD );

                #ifdef MPI_DEBUG
                    cout << g_hostString << " is sending " << batchSize << " proteins." << endl;
                    Timer sendTime(true);
                #endif

                string pack = packStream.str();
                int len = (int) pack.length();

                MPI_Send( &len,                    1,        MPI_INT,    sourceProcess,    0x00, MPI_COMM_WORLD );
                MPI_Send( (void*) pack.c_str(),    len,    MPI_CHAR,    sourceProcess,    0x01, MPI_COMM_WORLD );

                #ifdef MPI_DEBUG
                    cout << g_hostString << " finished sending " << batchSize << " proteins; " <<
                            sendTime.End() << " seconds elapsed." << endl;
                #endif

                i += batchSize;
            } else
            {
                batchSize = 0;
                MPI_Ssend( &batchSize,    1,    MPI_INT,    sourceProcess,    0x99, MPI_COMM_WORLD );

                #ifdef MPI_DEBUG
                    cout << "Process #" << sourceProcess << " has been informed that all proteins have been searched." << endl;
                #endif

                ++numChildrenFinished;
            }

            totalSearchTime = searchTime.TimeElapsed();
            if( !IsFinished && ( ( totalSearchTime - lastUpdate > g_rtConfig->StatusUpdateFrequency ) || i+1 == numProteins ) )
            {
                if( i+1 == numProteins )
                    IsFinished = true;

                float proteinsPerSec = float(i+1) / totalSearchTime;
                bpt::time_duration estimatedTimeRemaining(0, 0, round((numProteins - i) / proteinsPerSec));

                cout << "Searched " << i << " of " << numProteins << " proteins; "
                     << round(proteinsPerSec) << " per second, "
                     << format_date_time("%H:%M:%S", bpt::time_duration(0, 0, round(totalSearchTime))) << " elapsed, "
                     << format_date_time("%H:%M:%S", estimatedTimeRemaining) << " remaining." << endl;

                lastUpdate = totalSearchTime;
//.........这里部分代码省略.........

    int ReceiveResultsFromChildProcesses(bool firstBatch = false)
    {
        int numSpectra;
        int sourceProcess;

        Timer ResultsTime( true );
        float totalResultsTime = 0.01f;
        float lastUpdate = 0.0f;

        for( int p=0; p < g_numChildren; ++p )
        {
            MPI_Recv( &sourceProcess,        1,        MPI_INT,    MPI_ANY_SOURCE,    0xEE, MPI_COMM_WORLD, &st );

            #ifdef MPI_DEBUG
                cout << g_hostString << " is receiving search results." << endl;
                Timer receiveTime(true);
            #endif

            string pack;
            int len;

            MPI_Recv( &len,                    1,        MPI_INT,    sourceProcess,    0x00, MPI_COMM_WORLD, &st );
            pack.resize( len );
            MPI_Recv( (void*) pack.data(),    len,    MPI_CHAR,    sourceProcess,    0x01, MPI_COMM_WORLD, &st );

            stringstream compressedStream( pack );
            stringstream packStream;
            boost::iostreams::filtering_ostream decompressorStream;
            decompressorStream.push( boost::iostreams::zlib_decompressor() );
            decompressorStream.push( packStream );
            boost::iostreams::copy( compressedStream, decompressorStream );
            decompressorStream.reset();

            binary_iarchive packArchive( packStream );

            try
            {
                SearchStatistics childSearchStats;
                packArchive & numSpectra;
                packArchive & childSearchStats;
                if(firstBatch)
                {
                    searchStatistics = searchStatistics + childSearchStats;
                }
                else 
                {
                    searchStatistics.numPeptidesGenerated += childSearchStats.numPeptidesGenerated;
                    searchStatistics.numVariantsGenerated += childSearchStats.numVariantsGenerated;
                    searchStatistics.numComparisonsDone += childSearchStats.numComparisonsDone;
                    searchStatistics.numPeptidesSkipped += childSearchStats.numPeptidesSkipped;
                }

                //cout << g_hostString << " is unpacking results for " << numSpectra << " spectra." << endl;
                for( SpectraList::iterator sItr = spectra.begin(); sItr != spectra.end(); ++sItr )
                {
                    Spectrum* childSpectrum = new Spectrum;
                    Spectrum* rootSpectrum = *sItr;
                    packArchive & *childSpectrum;
                    rootSpectrum->numTargetComparisons += childSpectrum->numTargetComparisons;
                    rootSpectrum->numDecoyComparisons += childSpectrum->numDecoyComparisons;
                    rootSpectrum->processingTime += childSpectrum->processingTime;

                    rootSpectrum->resultsByCharge.resize(childSpectrum->resultsByCharge.size());
                    for (size_t z=0; z < childSpectrum->resultsByCharge.size(); ++z)
                    {
                        Spectrum::SearchResultSetType& rootResults = rootSpectrum->resultsByCharge[z];
                        Spectrum::SearchResultSetType& childResults = childSpectrum->resultsByCharge[z];

                        BOOST_FOREACH(const Spectrum::SearchResultPtr& result, childResults)
                            rootResults.add( result );

                        if (childResults.bestFullySpecificTarget().get()) rootResults.add(childResults.bestFullySpecificTarget());
                        if (childResults.bestFullySpecificDecoy().get()) rootResults.add(childResults.bestFullySpecificDecoy());
                        if (childResults.bestSemiSpecificTarget().get()) rootResults.add(childResults.bestSemiSpecificTarget());
                        if (childResults.bestSemiSpecificDecoy().get()) rootResults.add(childResults.bestSemiSpecificDecoy());
                        if (childResults.bestNonSpecificTarget().get()) rootResults.add(childResults.bestNonSpecificTarget());
                        if (childResults.bestNonSpecificDecoy().get()) rootResults.add(childResults.bestNonSpecificDecoy());
                    }

                    for(flat_map<int,int>::iterator itr = childSpectrum->mvhScoreDistribution.begin(); itr != childSpectrum->mvhScoreDistribution.end(); ++itr)
                        rootSpectrum->mvhScoreDistribution[(*itr).first] += (*itr).second;
                    for(flat_map<int,int>::iterator itr = childSpectrum->mzFidelityDistribution.begin(); itr != childSpectrum->mzFidelityDistribution.end(); ++itr)
                        rootSpectrum->mzFidelityDistribution[(*itr).first] += (*itr).second;
                    rootSpectrum->scoreHistogram += childSpectrum->scoreHistogram;
                    delete childSpectrum;
                }
                //cout << g_hostString << " is finished unpacking results." << endl;
            } catch( exception& e )
            {
                cerr << g_hostString << " had an error: " << e.what() << endl;
                exit(1);
            }

            #ifdef MPI_DEBUG
                cout << g_hostString << " finished receiving " << numSpectra << " search results; " <<
                        receiveTime.End() << " seconds elapsed.";
            #endif

            totalResultsTime = ResultsTime.TimeElapsed();
            if( ( totalResultsTime - lastUpdate > g_rtConfig->StatusUpdateFrequency ) || p+1 == g_numChildren )
//.........这里部分代码省略.........

void two_d_partitioning(MPI_Comm *comm_new, float *A, int local_rank, int num_procs) {
  MPI_Status status;
  int k, i, j, startingRow, endingRow, numRows, startingColumn, endingColumn, numColumns;
  int n_startingRow, n_startingColumn, n_local_coords[2];
  //long double determinant;
  double start, end, dt;
  int p = (int) sqrt(num_procs);
  int dis, left_rank, right_rank, up_rank, down_rank;
  MPI_Request req;
  
  numRows = n / p;
  numColumns = numRows;

  startingRow = local_coords[1] * numRows;
  endingRow = startingRow + numRows;

  startingColumn = local_coords[0] * numRows;
  endingColumn = startingColumn + numColumns;

  
  start = MPI_Wtime();
  
  for( k = 0; k < n; k++ ) {
    float Akk[1];
    int local_k = k % numRows;
    
    // Send A(k,k) to the right
    start = MPI_Wtime();
    if( k >= startingColumn && k < endingColumn && k >= startingRow && k < endingRow ) {
      send_to(comm_new, 0, A, 1, local_k, local_k, numRows);
      Akk[0] = A[local_k * numRows + local_k];
    } else if( k < startingColumn && k >= startingRow && k < endingRow ) {
      receive_from_left(comm_new, 0, Akk, 1, 0, 0, numRows, k);
    }
    end = MPI_Wtime();
    dt = end - start;
    comm_time += dt;

    // Now calculate the row
    start = MPI_Wtime();
    if( k >= startingColumn && k < endingColumn && k >= startingRow && k < endingRow ) {
      for( j = local_k + 1; j < numColumns; j++ ) {
        A[local_k * numRows + j] /= Akk[0];
      }
    } else if( k >= startingRow && k < endingRow && k < startingColumn ) {
      for( j = 0; j < numColumns; j++ ) {
        A[local_k * numRows + j] /= Akk[0];
      }
    }
    end = MPI_Wtime();
    dt = end - start;
    proc_time += dt;

    // Now calculate the box
    int m, bOutside = 1; 
    float top_row[numRows]; 

    start = MPI_Wtime();
    // k is West of this Partition
    if( k >= startingRow && k < endingRow & k < startingColumn ) {
      send_to(comm_new, 1, A, numColumns, local_k, 0, numRows);
      
      for( m = 0; m < numColumns; m++ ) {
        top_row[m] = A[local_k * numRows + m];
      }
      bOutside = -1;
    } 
    // k is in this BOX
    else if( k >= startingRow && k < endingRow && k >= startingColumn && k < endingColumn ) {
      int size = numColumns - (local_k + 1);
      if( size != 0 ) {
        send_to(comm_new, 1, A, size, local_k, local_k + 1, numRows);

        for( m = 0; m < size; m++ ) {
          top_row[m] = A[local_k * numRows + local_k + 1 + m];
        }
        bOutside = -1;
      }
    } // k is NW of this box 
    else if( k < startingRow && k < startingColumn ) {
      int sender_row = k / numRows;
      int sender_column = k / numColumns;
      int sender_rank = local_coords[0] * sqrt(num_procs) + sender_row;
      
      MPI_Recv(top_row, numColumns, MPI_FLOAT, sender_rank, 0, *comm_new, &status);
      
      bOutside = -1;
    } 
    // k is N of this box
    else if( k < startingRow && k >= startingColumn && k < endingColumn ) {
      int sender_row = k / numRows;
      int sender_column = k / numColumns;
      int sender_rank = sender_column * sqrt(num_procs) + sender_row;
      int size = numColumns - (local_k + 1);
      
      if( size != 0 ) { 
        //top_row = (float *)malloc(sizeof(float) * numberToReceive);
        //printf("%d Waiting to receive from:%d\n", local_rank, sender_rank);
        MPI_Recv(top_row, size, MPI_FLOAT, sender_rank, 0, *comm_new, &status);
        
//.........这里部分代码省略.........

    int TransmitUnpreparedSpectraToChildProcesses()
    {
        int numSpectra = (int) spectra.size();

        int sourceProcess, batchSize;
        bool IsFinished = false;

        Timer PrepareTime( true );
        float totalPrepareTime = 0.01f;
        float lastUpdate = 0.0f;

        int i = 0;
        int numChildrenFinished = 0;
        while( numChildrenFinished < g_numChildren )
        {
            stringstream packStream;
            binary_oarchive packArchive( packStream );
            // For every batch, listen for a worker process that is ready to receive it

            #ifdef MPI_DEBUG
                cout << g_hostString << " is listening for a child process to offer to prepare some spectra." << endl;
            #endif

            if( i < numSpectra )
            {
                batchSize = min( numSpectra-i, g_rtConfig->SpectraBatchSize );

                try
                {
                    packArchive & batchSize;

                    SpectraList::iterator sItr = spectra.begin();
                    advance( sItr, i );
                    for( int j = i; j < i + batchSize; ++j, ++sItr )
                    {
                        packArchive & **sItr;
                    }
                } catch( exception& e )
                {
                    cerr << g_hostString << " had an error: " << e.what() << endl;
                    exit(1);
                }

                i += batchSize;
            } else
            {
                batchSize = 0;
                packArchive & batchSize;

                #ifdef MPI_DEBUG
                    cout << "Process #" << sourceProcess << " has been informed that preparation is complete." << endl;
                #endif

                ++numChildrenFinished;
            }

            MPI_Recv( &sourceProcess,        1,        MPI_INT,    MPI_ANY_SOURCE,    0xFF, MPI_COMM_WORLD, &st );

            stringstream compressedStream;
            boost::iostreams::filtering_ostream compressorStream;
            compressorStream.push( boost::iostreams::zlib_compressor() );
            compressorStream.push( compressedStream );
            boost::iostreams::copy( packStream, compressorStream );
            compressorStream.reset();

            string pack = compressedStream.str();
            int len = (int) pack.length();

            MPI_Send( &len,                    1,        MPI_INT,    sourceProcess,    0x00, MPI_COMM_WORLD );
            MPI_Send( (void*) pack.c_str(),    len,    MPI_CHAR,    sourceProcess,    0x01, MPI_COMM_WORLD );

            totalPrepareTime = PrepareTime.TimeElapsed();
            if( !IsFinished && ( ( totalPrepareTime - lastUpdate > g_rtConfig->StatusUpdateFrequency ) || i == numSpectra ) )
            {
                if( i == numSpectra )
                    IsFinished = true;

                float spectraPerSec = float(i) / totalPrepareTime;
                float estimatedTimeRemaining = float(numSpectra-i) / spectraPerSec;
                cout << "Prepared " << i << " of " << numSpectra << " spectra; " << spectraPerSec <<
                        " per second, " << estimatedTimeRemaining << " seconds remaining." << endl;
                lastUpdate = totalPrepareTime;
            }
        }

        return 0;
    }

int main(int argc, char **argv)
{
	int proc_id, n_procs, i, envios;
	int buffer0[SIZE0], buffer1[SIZE1], buffer2[SIZE2], buffer3[SIZE3], buffer4[SIZE4];
	int buffer5[SIZE5], buffer6[SIZE6], buffer7[SIZE7], buffer8[SIZE8], buffer9[SIZE9];
	MPI_Status estado;

	/* Inicio del entorno MPI */
	MPI_Init (&argc, &argv);

	/* Obtener rango y tamaño del comm_world */
	MPI_Comm_rank(MPI_COMM_WORLD, &proc_id);
	MPI_Comm_size(MPI_COMM_WORLD, &n_procs);

	if (n_procs != 2)
	{
		printf("## Este programa utiliza 2 procesos\n");
		exit(EXIT_FAILURE);
	}

	/* Sincronizar a todos los procesos */
	MPI_Barrier(MPI_COMM_WORLD);

	/* Unicamente el maestro imprime en pantalla los datos iniciales */ 
	if (proc_id == MASTER)
	{
		system("clear");
		printf("## Cálculo de RTT / Envío Maestro-Esclavo\n\n");
		printf("## Total de procesadores: %d\n", n_procs);

		obtener_info_sist();
		/* Leer de stdin la cantidad de envíos a hacer */
		do
		{
			printf("%s", MENU_ENVIOS);
			scanf("%d", &envios);
		} while ((envios != 1) && (envios != 2) && (envios != 3));

		switch (envios)
		{
			case 1:
				envios = 100;
				break;
			case 2:
				envios = 1000;
				break;
			case 3:
				envios = 10000;
				break;
			default:
				break;
		}

		/* Enviar dato a los procesos */
		MPI_Send(&envios, sizeof(envios), MPI_INT, ESCLAVO, 0, MPI_COMM_WORLD);

		printf("## Comienzo del envío con: %d envios;\n", envios);
		printf("## Tamaño de datos:\n");
		printf("## \tSIZE0 = %d\n", SIZE0);
		printf("## \tSIZE1 = %d\n", SIZE1);
		printf("## \tSIZE2 = %d\n", SIZE2);
		printf("## \tSIZE3 = %d\n", SIZE3);
		printf("## \tSIZE4 = %d\n", SIZE4);
		printf("## \tSIZE5 = %d\n", SIZE5);
		printf("## \tSIZE6 = %d\n", SIZE6);
		printf("## \tSIZE7 = %d\n", SIZE7);
		printf("## \tSIZE8 = %d\n", SIZE8);
		printf("## \tSIZE9 = %d\n", SIZE9);
		printf("\n");
	}

	MPI_Barrier(MPI_COMM_WORLD);
	if (proc_id == MASTER)
	{
		printf("## | tam (b) | envios | t_total (seg) | t/envio (s)  |     kB/s   |\n");
		printf("## |---------|--------|---------------|--------------|------------|\n");
		master_func(SIZE0, envios, estado);
		master_func(SIZE1, envios, estado);
		master_func(SIZE2, envios, estado);
		master_func(SIZE3, envios, estado);
		master_func(SIZE4, envios, estado);
		master_func(SIZE5, envios, estado);
		master_func(SIZE6, envios, estado);
		master_func(SIZE7, envios, estado);
		master_func(SIZE8, envios, estado);
		master_func(SIZE9, envios, estado);
		printf("## |_________|________|_______________|______________|____________|\n");
	}
	else
	{
		/* Primero, recibir del maestro el dato de envios */
		MPI_Recv(&envios, sizeof(envios), MPI_INT, MASTER, 0, MPI_COMM_WORLD, &estado);

		/* Ahora realizar el recibo y envío con esa cantidad */
		for(i = 0; i < envios; ++i)
		{
			MPI_Recv(buffer0, SIZE0, MPI_INT, MASTER, 0, MPI_COMM_WORLD, &estado);
			buffer0[0] += 1;
			MPI_Send(buffer0, SIZE0, MPI_INT, MASTER, 0, MPI_COMM_WORLD);
		}
//.........这里部分代码省略.........

//***************************************************************************************************************
int Bellerophon::getChimeras() {
	try {
		
		//create breaking points
		vector<int> midpoints;   midpoints.resize(iters, window);
		for (int i = 1; i < iters; i++) {  midpoints[i] = midpoints[i-1] + increment;  }
	
	#ifdef USE_MPI
		int pid, numSeqsPerProcessor; 
	
		MPI_Comm_rank(MPI_COMM_WORLD, &pid); //find out who we are
		MPI_Comm_size(MPI_COMM_WORLD, &processors); 
		
		numSeqsPerProcessor = iters / processors;
		
		//each process hits this only once
		unsigned long long startPos = pid * numSeqsPerProcessor;
		if(pid == processors - 1){
				numSeqsPerProcessor = iters - pid * numSeqsPerProcessor;
		}
		lines.push_back(linePair(startPos, numSeqsPerProcessor));
		
		//fill pref with scores
		driverChimeras(midpoints, lines[0]);
		
		if (m->control_pressed) { return 0; }
				
		//each process must send its parts back to pid 0
		if (pid == 0) {
			
			//receive results 
			for (int j = 1; j < processors; j++) {
				
				vector<string>  MPIBestSend; 
				for (int i = 0; i < numSeqs; i++) {
				
					if (m->control_pressed) { return 0; }

					MPI_Status status;
					//receive string
					int length;
					MPI_Recv(&length, 1, MPI_INT, j, 2001, MPI_COMM_WORLD, &status);
					
					char* buf = new char[length];
					MPI_Recv(&buf[0], length, MPI_CHAR, j, 2001, MPI_COMM_WORLD, &status);
					
					string temp = buf;
					if (temp.length() > length) { temp = temp.substr(0, length); }
					delete buf;

					MPIBestSend.push_back(temp);
				}
				
				fillPref(j, MPIBestSend);
				
				if (m->control_pressed) { return 0; }
			}

		}else {
			//takes best window for each sequence and turns Preference to string that can be parsed by pid 0.
			//played with this a bit, but it may be better to try user-defined datatypes with set string lengths??
			vector<string> MPIBestSend = getBestWindow(lines[0]);
			pref.clear();
				
			//send your result to parent
			for (int i = 0; i < numSeqs; i++) {
				
				if (m->control_pressed) { return 0; }
				
				int bestLength = MPIBestSend[i].length();
				char* buf = new char[bestLength];
				memcpy(buf, MPIBestSend[i].c_str(), bestLength);
				
				MPI_Send(&bestLength, 1, MPI_INT, 0, 2001, MPI_COMM_WORLD);
				MPI_Send(buf, bestLength, MPI_CHAR, 0, 2001, MPI_COMM_WORLD);
				delete buf;
			}
			
			MPIBestSend.clear();
		}
		MPI_Barrier(MPI_COMM_WORLD); //make everyone wait - just in case
	#else
	
		//divide breakpoints between processors
		#if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
			if(processors == 1){ 
				lines.push_back(linePair(0, iters));	
				
				//fill pref with scores
				driverChimeras(midpoints, lines[0]);
	
			}else{
			
				int numSeqsPerProcessor = iters / processors;
				
				for (int i = 0; i < processors; i++) {
					unsigned long long startPos = i * numSeqsPerProcessor;
					if(i == processors - 1){
						numSeqsPerProcessor = iters - i * numSeqsPerProcessor;
//.........这里部分代码省略.........

//removes the lowest energy vertical seam from the image
void removeHorizontalSeam() {
    double energies[3];
    double min_energy;
    int prev_x;
    int prev_y;
    
    // split up work between processes
    double *my_path_costs;
    double *my_previous_x;
    double *my_previous_y;
    double *temp_path_costs;
    double *temp_previous_x;
    double *temp_previous_y;
    int my_rows = current_height / numprocs;
    int low_rows = my_rows;
    int extra_rows = current_height % numprocs;
    int start;
    int y_offset;
    int recv_rows;

    double top_end_cost, bottom_end_cost, temp_end_cost;

    if (rank < extra_rows) {
        my_rows++;
        start = rank * my_rows;
    } else {
        start = (extra_rows * (my_rows + 1)) + ((rank - extra_rows) * my_rows);
    }
    
    my_path_costs = (double *) malloc(my_rows * current_width * sizeof(double));
    my_previous_x = (double *) malloc(my_rows * current_width * sizeof(double));
    my_previous_y = (double *) malloc(my_rows * current_width * sizeof(double));
    //find the lowest cost seam by computin