fork download
copy 
  1.  
  2. #include <vector>
  3. // #include <cmath>
  4. #include <fstream>
  5. #include <iostream>
  6. #include <iomanip>
  7. #include <time.h>
  8. #define _USE_MATH_DEFINES
  9. #include <math.h>
  10.  
  11. using namespace std;
  12.  
  13. // Function to compute source term q
  14. double sourceTerm(double x, double y) {
  15.     return -(x * x + y * y);
  16. }
  17.  
  18. // Function to initialize phi grid
  19. void initializePhi(vector<vector<double > >& phi, int n) {
  20.     for (int i = 0; i < n; ++i) {
  21.         for (int j = 0; j < n; ++j) {
  22.             phi[i][j] = 0.0;
  23.         }
  24.     }
  25. }
  26.  
  27.  
  28. // Main function
  29. int main(int argc, char** argv) {
  30.  
  31.     const int N = 401;  // Total points in each direction
  32.     const double deltaX = 0.005, deltaY = 0.005;
  33.     const double tol = 1e-4;
  34.  
  35.     clock_t start=clock();
  36.  
  37.     vector<vector<double > > phi(N, vector<double>(N, 0.0)); // Including halo rows
  38.     vector<vector<double > > phi_new(N, vector<double>(N, 0.0));
  39.     initializePhi(phi, N);
  40.     initializePhi(phi_new, N);
  41.  
  42.     double maxError;
  43.     int iterations = 0;
  44.  
  45.  
  46.     do {
  47.         for (int i=1; i < N - 1; i++)
  48.         {
  49.             double y = i*deltaY;
  50.             phi_new[i][0] = sin(2 * M_PI * y);
  51.         }
  52.  
  53.         for (int i = 1; i < N - 1; ++i) {
  54.             for (int j = 1; j < N - 1; ++j) {
  55.                 if((i+j)%2==0)
  56.                     phi_new[i][j] = 0.25 * (phi[i + 1][j] + phi[i - 1][j] + phi[i][j + 1] + phi[i][j - 1] +
  57.                                         deltaX * deltaY * sourceTerm(i * deltaX, j * deltaY));
  58.             }
  59.         }
  60.  
  61.         for (int i = 1; i < N - 1; ++i) {
  62.             for (int j = 1; j < N - 1; ++j) {
  63.                 if((i+j)%2==1)
  64.                     phi_new[i][j] = 0.25 * (phi_new[i + 1][j] + phi_new[i - 1][j] + phi_new[i][j + 1] + phi_new[i][j - 1] +
  65.                                         deltaX * deltaY * sourceTerm(i * deltaX, j * deltaY));
  66.             }
  67.         }
  68.  
  69.  
  70.         for (int i=1; i < N - 1; i++)
  71.         {
  72.             phi_new[i][N-1] = (4 * phi_new[i][N-2] - phi_new[i][N-3])/3.0;
  73.         }
  74.  
  75.         maxError = 0.0;
  76.         for (int i = 0; i < N; ++i) {
  77.             for (int j = 0; j < N; ++j) {
  78.                 maxError = max(maxError, fabs(phi_new[i][j] - phi[i][j])/phi[i][j]);
  79.                 phi[i][j] = phi_new[i][j];
  80.             }
  81.         }
  82.  
  83.         iterations++;
  84.     } while (maxError > tol);
  85.  
  86.     clock_t end=clock();
  87.  
  88.     double time_spent=(double) (end-start)/CLOCKS_PER_SEC;
  89.  
  90.     cout<<"time taken "<<time_spent<<endl;
  91.     string filename = "que_2d_serial.txt";
  92.     FILE *outputFile;
  93.     outputFile = fopen(filename.c_str(), "w");
  94.  
  95.     fprintf(outputFile, "Convergence reached after %d iterations \n \n", iterations);
  96.     cout << "Convergence reached after " << iterations << " iterations." << endl;
  97.  
  98.     fprintf(outputFile, "Total 201 points \n");
  99.     // cout<<"Total 201 points"<<endl;
  100.     fprintf(outputFile, "[");
  101.     for(int i = 0; i < N; i++)
  102.     {
  103.         fprintf(outputFile, "%f, ", (-1 + i * deltaX));
  104.         // cout<< (i * deltaX) <<" ";
  105.     }
  106.     fprintf(outputFile, "] \n \n");
  107.  
  108.     fprintf(outputFile, "Numerical Solution when x = %f \n", -1 + (N/2)*deltaX);
  109.     // cout<<"Numerical Solution when x = "<< (N/2)*deltaX <<endl;
  110.     fprintf(outputFile, "[");
  111.     for(int j = 0; j < N; j++)
  112.     {
  113.         fprintf(outputFile, "%f, ", phi[N/2][j]);
  114.         // cout<<phi[N/2][j]<<" ";
  115.     }
  116.     fprintf(outputFile, "] \n \n");
  117.  
  118.     fprintf(outputFile, "Numerical Solution when y = %f \n", -1 + (N/2)*deltaY);
  119.     // cout<<"Numerical Solution when x = "<< (N/2)*deltaX <<endl;
  120.     fprintf(outputFile, "[");
  121.     for(int i = 0; i < N; i++)
  122.     {
  123.         fprintf(outputFile, "%f, ", phi[i][N/2]);
  124.         // cout<<phi[i][N/2]<<" ";
  125.     }
  126.     fprintf(outputFile, "] \n \n");
  127.  
  128.     fclose(outputFile);
  129.  
  130.     return 0;
  131. }
  132.  
  133.  
Success #stdin #stdout #stderr 0.29s 40896KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
Standard output is empty
stderr
copy 
Error: unexpected '/' in "/"
Execution halted
https://ideone.com/shYkuU
language:
R (R 3.5.2)
created:
10 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!