calc_integral_openmp/stage_3.cpp

//#define DEBUG

#include "mpi.h"
#include <cstdio>
#include <cmath>
#include <list>
#include <tuple>

#include "data_format.hpp"
#include "data_output.hpp"
#include "function.hpp"
#include "integral_calc.hpp"

Interval_t get_local_interval(double inf, double sup, size_t interval_count, size_t current_rank, size_t proc_count);

int main(int argc, char* argv[]){
  freopen("input.txt", "r", stdin);

// Локальные параметры
  int proc_count;
  int current_rank;
  MPI_Status mpi_status;

  size_t IO_PROC_RANK = 0;
  size_t INTERVAL_MSG = 0;
  size_t RESULT_MSG = 1;

  size_t test_count = 0;
  double inf = 0;
  double sup = 0;
  double param = 0;
  
  std::list<trip_t> integral_list;
//-------------------------------------------------------------

// MPI - часть
  MPI_Init(&argc, &argv); 
  
  MPI_Comm_size(MPI_COMM_WORLD, &proc_count);
  MPI_Comm_rank(MPI_COMM_WORLD, &current_rank);

  std::cin >> inf >> sup >> param;
  std::cin >> test_count;

  for(int i = 0; i < test_count; ++i){

    size_t interval_count = 0;
    double integral_val = 0;
    double calc_param = 0;
    bool out_flag = false;
    
    int mpi_bufpos = 0;
    char mpi_buf[4 * sizeof(double)+1];

    double start = MPI_Wtime();

//// Чтение и обмен входными данными   
    if(current_rank == IO_PROC_RANK){

      std::cin >> calc_param;
      if(calc_param >= 1){
        interval_count = calc_param;
        out_flag = true;
      }
      else{ 
        interval_count = 1.0 / calc_param;
        out_flag = false;
      }

#ifdef DEBUG
      std::cout << test_count << " " << interval_count << '\n';
#endif // DEBUG
       //
 
      MPI_Pack(&inf, 1, MPI_DOUBLE, mpi_buf, sizeof(mpi_buf), &mpi_bufpos, MPI_COMM_WORLD);
      MPI_Pack(&sup, 1, MPI_DOUBLE, mpi_buf, sizeof(mpi_buf), &mpi_bufpos, MPI_COMM_WORLD);
      MPI_Pack(&param, 1, MPI_DOUBLE, mpi_buf, sizeof(mpi_buf), &mpi_bufpos, MPI_COMM_WORLD);
      MPI_Pack(&interval_count, 1, MPI_UINT32_T, mpi_buf, sizeof(mpi_buf), &mpi_bufpos, MPI_COMM_WORLD);
      
    }
    
    MPI_Bcast(&mpi_buf
            , sizeof(mpi_buf)
            , MPI_PACKED
            , IO_PROC_RANK
            , MPI_COMM_WORLD);    

    if(current_rank != IO_PROC_RANK){
      MPI_Unpack(mpi_buf, sizeof(mpi_buf), &mpi_bufpos, &inf, 1, MPI_DOUBLE, MPI_COMM_WORLD);
      MPI_Unpack(mpi_buf, sizeof(mpi_buf), &mpi_bufpos, &sup, 1, MPI_DOUBLE, MPI_COMM_WORLD);
      MPI_Unpack(mpi_buf, sizeof(mpi_buf), &mpi_bufpos, &param, 1, MPI_DOUBLE, MPI_COMM_WORLD);
      MPI_Unpack(mpi_buf, sizeof(mpi_buf), &mpi_bufpos, &interval_count, 1, MPI_UINT32_T, MPI_COMM_WORLD);
    }
////-------------------------------------------------------------

//// Вычисления   
    Interval_t local_interval = get_local_interval(inf, sup, interval_count, current_rank, proc_count);

#ifdef DEBUG
    std::cout << ">>Local inf: " << local_interval.inf << '\n';
    std::cout << ">>Local sup: " << local_interval.sup << '\n';
    std::cout << ">>Local inc_count: " << local_interval.increase_count << '\n';
#endif // DEBUG
    
    double local_res = integral_calc(param, local_interval);

#ifdef DEBUG
    std::cout << "Local res: " << local_res << '\n';
#endif // DEBUG
////-------------------------------------------------------------

//// Обмен результатами вычислений   

    MPI_Reduce(&local_res
            , &integral_val
            , sizeof(local_res)/sizeof(double)
            , MPI_DOUBLE
            , MPI_SUM
            , IO_PROC_RANK
            , MPI_COMM_WORLD);    

    if(current_rank == IO_PROC_RANK){

#ifdef DEBUG
      std::cout << "Результат: " << integral_val << '\n';
      std::cout << "Образец: " << (primal(sup, param) - primal(inf, param)) << '\n';
#endif // DEBUG
      double exact_val = primal(sup, param) - primal(inf, param);
      double err = fabsl(integral_val - exact_val) / fabsl(integral_val);
      double difftime = MPI_Wtime() - start;

      integral_list.push_back(std::make_tuple(interval_count, err, difftime, out_flag));
    }
////-------------------------------------------------------------
  }
//// Вывод результатов
  if(current_rank == IO_PROC_RANK)
    print_table_2(integral_list);

////-------------------------------------------------------------
  
  MPI_Finalize();
//-------------------------------------------------------------
  return 0;
}

Interval_t get_local_interval(double inf, double sup, size_t interval_count, size_t current_rank, size_t proc_count){
  Interval_t res;
  double increase = (sup - inf) / (double) interval_count;
  size_t base = interval_count / proc_count;
  size_t residue = interval_count % proc_count;

  res.increase_count = base + (residue > current_rank ? 1 : 0);
  res.inf = inf + increase * 
    (base * current_rank + (residue > current_rank ? current_rank : residue)); 
  res.sup = res.inf + increase * res.increase_count;
  
  return res;
}