Fully working BasicSolver

BasicSolver/3D-seq/Makefile

@@ -40,6 +40,18 @@ $(BUILD_DIR)/%.s:  %.c
 
 .PHONY: clean distclean tags info asm format
 
+vis:
+	$(info ===>  GENERATE VISUALIZATION)
+	@gnuplot -e "filename='residual.dat'" ./residual.plot
+
+vis_clean:
+	$(info ===>  CLEAN VISUALIZATION)
+	@rm -f *.dat
+	@rm -f *.vtk
+	@rm -f *.png
+
+clean: vis_clean
+
 clean:
 	$(info ===>  CLEAN)
 	@rm -rf $(BUILD_DIR)

BasicSolver/3D-seq/canal.par

@@ -38,7 +38,7 @@ kmax          50	   # number of interior cells in z-direction
 # Time Data:
 # ---------
 
-te       100.0   # final time
+te       60.0   # final time
 dt       0.02    # time stepsize
 tau      0.5     # safety factor for time stepsize control (<0 constant delt)
 

BasicSolver/3D-seq/dcavity.par

@@ -38,7 +38,7 @@ kmax       128		# number of interior cells in z-direction
 # Time Data:
 # ---------
 
-te       2.0		# final time
+te       10.0		# final time
 dt       0.02	    # time stepsize
 tau      0.5		# safety factor for time stepsize control (<0 constant delt)
 
@@ -46,7 +46,7 @@ tau      0.5		# safety factor for time stepsize control (<0 constant delt)
 # ---------
 
 levels        3         # Multigrid levels
-presmooth     10         # Pre-smoothning iterations
+presmooth     20         # Pre-smoothning iterations
 postsmooth    5         # Post-smoothning iterations
 
 # Pressure Iteration Data:

BasicSolver/3D-seq/residual.plot

@@ -0,0 +1,9 @@
+set terminal png size 1800,768 enhanced font ,12
+set output 'residual.png'
+set datafile separator whitespace
+set xlabel "Timestep"
+set ylabel "Residual"
+
+set logscale y 2
+
+plot 'residual.dat' using 1:2 title "Residual"

BasicSolver/3D-seq/src/main.c

@@ -73,6 +73,9 @@ int main(int argc, char** argv)
     Solver s;
     initParameter(&p);
 
+    FILE* fp;
+    fp = initResidualWriter();
+
     if (argc != 2) {
         printf("Usage: %s <configFile>\n", argv[0]);
         exit(EXIT_SUCCESS);
@@ -90,6 +93,7 @@ int main(int argc, char** argv)
     double te  = d.te;
     double t   = 0.0;
     int nt     = 0;
+    double res = 0.0;
 
     timeStart = getTimeStamp();
     while (t <= te) {
@@ -99,8 +103,11 @@ int main(int argc, char** argv)
         computeFG(&d);
         computeRHS(&d);
         if (nt % 100 == 0) normalizePressure(&d);
-        solve(&s, d.p, d.rhs);
+        res = solve(&s, d.p, d.rhs);
         adaptUV(&d);
+
+        writeResidual(fp, t, res);
+
         t += d.dt;
         nt++;
 
@@ -123,7 +130,8 @@ int main(int argc, char** argv)
     ug = allocate(64, bytesize);
     vg = allocate(64, bytesize);
     wg = allocate(64, bytesize);
-
+    
+    fclose(fp);
     createBulkArrays(&d, pg, ug, vg, wg);
     VtkOptions opts = { .grid = d.grid };
     vtkOpen(&opts, d.problem);

BasicSolver/3D-seq/src/progress.c

@@ -49,3 +49,22 @@ void stopProgress()
     printf("\n");
     fflush(stdout);
 }
+
+FILE* initResidualWriter()
+{
+    FILE* fp;
+    fp = fopen("residual.dat", "w");
+
+    if (fp == NULL) {
+        printf("Error!\n");
+        exit(EXIT_FAILURE);
+    }
+
+    return fp;
+
+}
+
+void writeResidual(FILE* fp, double ts, double res)
+{
+    fprintf(fp, "%f, %f\n", ts, res);
+}

BasicSolver/3D-seq/src/progress.h

@@ -10,5 +10,6 @@
 extern void initProgress(double);
 extern void printProgress(double);
 extern void stopProgress(void);
-
+extern FILE* initResidualWriter(void);
+extern void writeResidual(FILE*, double, double);
 #endif

BasicSolver/3D-seq/src/solver-mg.c

@@ -292,11 +292,13 @@ void initSolver(Solver* s, Discretization* d, Parameter* p)
     }
 }
 
-void solve(Solver* s, double* p, double* rhs)
+double solve(Solver* s, double* p, double* rhs)
 {
     double res = multiGrid(s, p, rhs, 0, s->grid->imax, s->grid->jmax, s->grid->kmax);
 
 #ifdef VERBOSE
     printf("Residuum: %.6f\n", res);
 #endif
+
+return res;
 }

BasicSolver/3D-seq/src/solver-rb.c

@@ -15,7 +15,7 @@ void initSolver(Solver* s, Discretization* d, Parameter* p)
     s->omega   = p->omg;
 }
 
-void solve(Solver* s, double* p, double* rhs)
+double solve(Solver* s, double* p, double* rhs)
 {
     int imax      = s->grid->imax;
     int jmax      = s->grid->jmax;
@@ -96,4 +96,6 @@ void solve(Solver* s, double* p, double* rhs)
 #ifdef VERBOSE
     printf("Solver took %d iterations to reach %f\n", it, sqrt(res));
 #endif
+
+return res;
 }

BasicSolver/3D-seq/src/solver.h

@@ -22,6 +22,6 @@ typedef struct {
 } Solver;
 
 extern void initSolver(Solver*, Discretization*, Parameter*);
-extern void solve(Solver*, double*, double*);
+extern double solve(Solver*, double*, double*);
 
 #endif