From 6a4ded9ca79474c95a13d7c77fcbdee294c559da Mon Sep 17 00:00:00 2001
From: mejjete <danil.malapura@gmail.com>
Date: Tue, 12 Mar 2024 21:34:40 +0100
Subject: [PATCH] fixed Euler's step and number of iteration for dx = 5

---
 sequential.cpp | 68 ++++++++++++++++++--------------------------------
 1 file changed, 24 insertions(+), 44 deletions(-)

diff --git a/sequential.cpp b/sequential.cpp
index 611cb56..e9a0b98 100644
--- a/sequential.cpp
+++ b/sequential.cpp
@@ -1,39 +1,29 @@
 #include <iostream>
 #include <vector>
-#include <iomanip>
 #include <math.h>
+#include <iomanip>
 
 using namespace std;
 
 int main() 
 {
-    // Length
-    int l = 1925;
-    // Spatial step
-    int dx = 40;
-    // Temporal step
-    int dt = 2500;
-    // Number of spatial steps
-    int kmax = l / dx;
-
-    // Initializing parameters
-    vector<double> list(dt);
-
-    // Arrays to store data
-    double Q[kmax][dt];
-    double P[kmax][dt];
-
-    int i;
+    int l = 1925;             // Edge length
+    int dx = 40;               // Spatial step
+    int dt = 100000;          // Euler's step
+    int kmax = l / dx;        // Number of spatial steps
 
     // Parameters
     double alpha = 0.093;
     double beta = 0.004;
     double gamma = 728.6;
-    double h = 0.05;
-    double sum_of_h = 0;
+    double h = 0.007;
+
+    // Arrays to store data
+    double Q[kmax][dt];
+    double P[kmax][dt];
 
     // Loop over time steps
-    for (i = 0; i < dt; i++) 
+    for (int i = 0; i < dt; i++) 
     {
         // Loop over spatial steps
         for (int k = 0; k < kmax; k++) 
@@ -42,7 +32,7 @@ int main()
             if (k == kmax - 1) 
             {
                 Q[k][i] = 0;
-                P[k][i] = -202.0; // P = -202 when Q = 10
+                P[k][i] = -202.0;   // P = -202 when Q = 10
             } else 
             {
                 Q[k][i] = 0;
@@ -52,31 +42,20 @@ int main()
     }
 
     // Main calculation loop
-    for (i = 2; i < dt - 2; i++) 
+    for (int i = 2; i < dt; i++) 
     {
-        // Add current time step to the list
-        sum_of_h = sum_of_h + h;
-        list.push_back(sum_of_h);
-
-		std::cout << "Q at step " << i;
-		
         // Loop over spatial steps
         for (int k = 1; k < kmax; k++) 
         {
-            // Boundary condition for the last spatial step
-            if (k == kmax - 1) 
-            {
-                Q[k][i + 1] = Q[k][i] + h * (alpha * (P[k - 1][i] - P[k][i]) - beta * Q[k][i] * abs(Q[k][i]));
-            } else 
-            {
-                // Calculate Q and P for non-boundary spatial steps
-                Q[k][i + 1] = Q[k][i] + h * (alpha * (P[k - 1][i] - P[k][i]) - beta * Q[k][i] * abs(Q[k][i]));
+            Q[k][i + 1] = Q[k][i] + h * (alpha * (P[k - 1][i] - P[k][i]) - beta * Q[k][i] * abs(Q[k][i]));
+
+            // DO NOT calculate pressure for the last element such as it is given by default
+            if (k != kmax - 1) 
                 P[k][i + 1] = P[k][i - 2] + h * (gamma * (Q[k][i] - Q[k + 1][i]));
-            }
         }
     }
     
-    // Print out the first 10 values of the timestemp
+    // Print out the first 10 of Q for each spatial step
     std::cout << "FIRST 10 VALUES\n";
     for(int i = 0; i < kmax; i++)
     {
@@ -89,18 +68,19 @@ int main()
         std::cout << std::endl;
     }
 
-    // Print out the last 10 values of the timestemp
-    std::cout << "LAST 10 VALUES\n";
-    for(int i = 0; i < kmax; i++)
+    // Print out the first 10 of Q for each spatial step
+    std::cout << "LAST 10 VALUES:\n";
+    for(int i = 1; i < kmax; i++)
     {
         std::cout.width(5);
         std::cout.flags(std::ios::left);
         std::cout << i << ": ";
- 
+
         for(int k = 0; k < 10; k++)
-            std::cout << std::setw(10) << Q[i][2489 + k];
+            std::cout << std::setw(10) << Q[i][dt - 10 + k];
         std::cout << std::endl;
     }
 
+
     return 0;
 }