diff --git a/EnhancedSolver/2D-seq/src/solver-sor.c b/EnhancedSolver/2D-seq/src/solver-sor.c index f817b5f..165e103 100644 --- a/EnhancedSolver/2D-seq/src/solver-sor.c +++ b/EnhancedSolver/2D-seq/src/solver-sor.c @@ -4,6 +4,7 @@ * Use of this source code is governed by a MIT style * license that can be found in the LICENSE file. */ +#include "grid.h" #include "solver.h" #include "util.h" @@ -13,19 +14,32 @@ void initSolver(Solver* s, Discretization* d, Parameter* p) s->itermax = p->itermax; s->eps = p->eps; s->omega = p->omg; + + Grid* g = s->grid; + int imax = s->grid->imax; + int jmax = s->grid->jmax; + + s->totalFluidCells = 0; + for (int j = 0; j < jmax + 2; j++) { + for (int i = 0; i < imax + 2; i++) { + if (gridIsFluid(g, i, j)) { + s->totalFluidCells++; + } + } + } } -void solveSOR(Solver* solver, double* p, double* rhs) +void solveSOR(Solver* s, double* p, double* rhs) { - int imax = solver->grid->imax; - int jmax = solver->grid->jmax; - double eps = solver->eps; - int itermax = solver->itermax; - double dx2 = solver->grid->dx * solver->grid->dx; - double dy2 = solver->grid->dy * solver->grid->dy; + int imax = s->grid->imax; + int jmax = s->grid->jmax; + double eps = s->eps; + int itermax = s->itermax; + double dx2 = s->grid->dx * s->grid->dx; + double dy2 = s->grid->dy * s->grid->dy; double idx2 = 1.0 / dx2; double idy2 = 1.0 / dy2; - double factor = solver->omega * 0.5 * (dx2 * dy2) / (dx2 + dy2); + double factor = s->omega * 0.5 * (dx2 * dy2) / (dx2 + dy2); double epssq = eps * eps; int it = 0; double res = 1.0; @@ -55,7 +69,7 @@ void solveSOR(Solver* solver, double* p, double* rhs) P(imax + 1, j) = P(imax, j); } - res = res / (double)(imax * jmax); + res = res / (double)(s->totalFluidCells); #ifdef DEBUG printf("%d Residuum: %e\n", it, res); #endif @@ -80,6 +94,7 @@ void solve(Solver* solver, double* p, double* rhs) double factor = solver->omega * 0.5 * (dx2 * dy2) / (dx2 + dy2); double epssq = eps * eps; int it = 0; + Grid* g = solver->grid; double res = 1.0; int pass, jsw, isw; @@ -92,13 +107,15 @@ void solve(Solver* solver, double* p, double* rhs) for (int j = 1; j < jmax + 1; j++) { for (int i = isw; i < imax + 1; i += 2) { + if (gridIsFluid(g, i, j)) { + double r = RHS(i, j) - + ((P(i + 1, j) - 2.0 * P(i, j) + P(i - 1, j)) * idx2 + + (P(i, j + 1) - 2.0 * P(i, j) + P(i, j - 1)) * + idy2); - double r = RHS(i, j) - - ((P(i + 1, j) - 2.0 * P(i, j) + P(i - 1, j)) * idx2 + - (P(i, j + 1) - 2.0 * P(i, j) + P(i, j - 1)) * idy2); - - P(i, j) -= (factor * r); - res += (r * r); + P(i, j) -= (factor * r); + res += (r * r); + } } isw = 3 - isw; } diff --git a/EnhancedSolver/2D-seq/src/solver.h b/EnhancedSolver/2D-seq/src/solver.h index d7ce513..8fff905 100644 --- a/EnhancedSolver/2D-seq/src/solver.h +++ b/EnhancedSolver/2D-seq/src/solver.h @@ -1,5 +1,5 @@ /* - * Copyright (C) NHR@FAU, University Erlangen-Nuremberg. + * Copyright (C) NHR@FAU, University Erlangen-Nuremberg. * All rights reserved. This file is part of nusif-solver. * Use of this source code is governed by a MIT style * license that can be found in the LICENSE file. @@ -17,6 +17,7 @@ typedef struct { double eps, omega, rho; int itermax; int levels; + int totalFluidCells; double **r, **e; } Solver;