fork download
copy 
  1. #include <iostream>
  2. #include <vector>
  3. #include <cmath>
  4. #include <mpi.h>
  5.  
  6. const double PI = 3.14159265358979323846;
  7.  
  8. // Function to initialize temperature
  9. void initialize_temperature(std::vector<double>& T, int M, int rank, int size, double dx) {
  10.     int local_size = M / size;
  11.     for (int i = 0; i < local_size; ++i) {
  12.         double x = (rank * local_size + i) * dx;
  13.         T[i] = 100 * sin(PI * x); // Initial temperature distribution
  14.     }
  15. }
  16.  
  17. // Function to print temperature (gathered at rank 0)
  18. void print_temperature(const std::vector<double>& T, int M, int rank, int size, double dx) {
  19.     if (rank == 0) {
  20.         std::vector<double> full_T(M + 1, 0.0);
  21.         MPI_Gather(T.data(), M / size, MPI_DOUBLE, full_T.data(), M / size, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  22.  
  23.         std::cout << "Temperature:\n";
  24.         for (int i = 0; i <= M; ++i) {
  25.             std::cout << i * dx << " " << full_T[i] << "\n";
  26.         }
  27.         std::cout << "\n";
  28.     } else {
  29.         MPI_Gather(T.data(), M / size, MPI_DOUBLE, nullptr, 0, MPI_DOUBLE, 0, MPI_COMM_WORLD);
  30.     }
  31. }
  32.  
  33. // Function to perform the finite difference simulation
  34. void heat_conduction_simulation(int M, int time_steps, double alpha, double dx, double dt, int rank, int size) {
  35.     int local_size = M / size;
  36.     std::vector<double> T(local_size, 0.0);     // Temperature at current time step
  37.     std::vector<double> T_new(local_size, 0.0); // Temperature at next time step
  38.  
  39.     initialize_temperature(T, M, rank, size, dx);
  40.  
  41.     // Print initial state (gathered at rank 0)
  42.     if (rank == 0) std::cout << "Initial state:\n";
  43.     print_temperature(T, M, rank, size, dx);
  44.  
  45.     // Buffers for boundary exchange
  46.     double left_boundary = 0.0, right_boundary = 0.0;
  47.  
  48.     // Time-stepping loop
  49.     for (int t = 0; t < time_steps; ++t) {
  50.         // Exchange boundary data
  51.         if (rank > 0) {
  52.             MPI_Send(&T[0], 1, MPI_DOUBLE, rank - 1, 0, MPI_COMM_WORLD);
  53.             MPI_Recv(&left_boundary, 1, MPI_DOUBLE, rank - 1, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
  54.         }
  55.         if (rank < size - 1) {
  56.             MPI_Send(&T[local_size - 1], 1, MPI_DOUBLE, rank + 1, 0, MPI_COMM_WORLD);
  57.             MPI_Recv(&right_boundary, 1, MPI_DOUBLE, rank + 1, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
  58.         }
  59.  
  60.         // Update interior points
  61.         for (int i = 1; i < local_size - 1; ++i) {
  62.             T_new[i] = T[i] + alpha * dt / (dx * dx) * (T[i + 1] - 2 * T[i] + T[i - 1]);
  63.         }
  64.  
  65.         // Update boundary points
  66.         if (rank > 0) {
  67.             T_new[0] = T[0] + alpha * dt / (dx * dx) * (T[1] - 2 * T[0] + left_boundary);
  68.         }
  69.         if (rank < size - 1) {
  70.             T_new[local_size - 1] = T[local_size - 1] + alpha * dt / (dx * dx) *
  71.                                      (right_boundary - 2 * T[local_size - 1] + T[local_size - 2]);
  72.         }
  73.  
  74.         // Swap T and T_new
  75.         T.swap(T_new);
  76.  
  77.         // Print intermediate states (optional, at rank 0)
  78.         if (t % (time_steps / 10) == 0 && rank == 0) {
  79.             std::cout << "Time step " << t << ":\n";
  80.             print_temperature(T, M, rank, size, dx);
  81.         }
  82.     }
  83.  
  84.     // Final state (gathered at rank 0)
  85.     if (rank == 0) std::cout << "Final state:\n";
  86.     print_temperature(T, M, rank, size, dx);
  87. }
  88.  
  89. int main(int argc, char** argv) {
  90.     MPI_Init(&argc, &argv);
  91.  
  92.     int rank, size;
  93.     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
  94.     MPI_Comm_size(MPI_COMM_WORLD, &size);
  95.  
  96.     int N = 1;                  // Rod length in meters
  97.     int M = 50;                 // Number of spatial divisions
  98.     double dx = static_cast<double>(N) / M; // Spatial step size
  99.     double alpha = 0.01;        // Thermal diffusivity
  100.     double dt = 0.0001;         // Time step size
  101.     int time_steps = 1000;      // Number of time steps
  102.  
  103.     if (M % size != 0) {
  104.         if (rank == 0) {
  105.             std::cerr << "Error: Number of spatial divisions (M) must be divisible by number of processes.\n";
  106.         }
  107.         MPI_Finalize();
  108.         return 1;
  109.     }
  110.  
  111.     heat_conduction_simulation(M, time_steps, alpha, dx, dt, rank, size);
  112.  
  113.     MPI_Finalize();
  114.     return 0;
  115. }
Success #stdin #stdout #stderr 0.28s 40588KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
Standard output is empty
stderr
copy 
Error: unexpected symbol in "const double"
Execution halted
https://ideone.com/lBeJiV
language:
R (R 3.5.2)
created:
2 days ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!