diff --git a/examples/reference-FTCS_2D_closed.cpp b/examples/reference-FTCS_2D_closed.cpp index 2323531..d755e15 100644 --- a/examples/reference-FTCS_2D_closed.cpp +++ b/examples/reference-FTCS_2D_closed.cpp @@ -5,8 +5,8 @@ using namespace std; int main(int argc, char *argv[]) { - int row = 11; - int col = 11; + int row = 50; + int col = 50; int domain_row = 10; int domain_col = 10; @@ -45,14 +45,11 @@ int main(int argc, char *argv[]) { // Simulation Simulation sim = Simulation(grid, bc, FTCS_APPROACH); sim.setTimestep(0.001); - sim.setIterations(7000); - sim.setOutputCSV(CSV_OUTPUT_ON); - sim.setOutputConsole(CONSOLE_OUTPUT_ON); + sim.setIterations(100); + sim.setOutputCSV(CSV_OUTPUT_OFF); + sim.setOutputConsole(CONSOLE_OUTPUT_OFF); // RUN sim.run(); - - cout << grid.getConcentrations() << endl; - } \ No newline at end of file diff --git a/src/FTCS.cpp b/src/FTCS.cpp index afd7404..6de288e 100644 --- a/src/FTCS.cpp +++ b/src/FTCS.cpp @@ -9,6 +9,7 @@ #include #include #include +#include using namespace std; @@ -276,6 +277,8 @@ void FTCS_2D(Grid &grid, Boundary &bc, double ×tep) { // inner cells // these are independent of the boundary condition type + omp_set_num_threads(10); + #pragma omp parallel for for (int row = 1; row < rowMax-1; row++) { for (int col = 1; col < colMax-1; col++) { concentrations_t1(row, col) = grid.getConcentrations()(row, col) @@ -295,6 +298,7 @@ void FTCS_2D(Grid &grid, Boundary &bc, double ×tep) { // left without corners / looping over rows // hold column constant at index 0 int col = 0; + #pragma omp parallel for for (int row = 1; row < rowMax-1; row++) { concentrations_t1(row, col) = grid.getConcentrations()(row,col) + timestep / (deltaCol*deltaCol) @@ -311,6 +315,7 @@ void FTCS_2D(Grid &grid, Boundary &bc, double ×tep) { // right without corners / looping over rows // hold column constant at max index col = colMax-1; + #pragma omp parallel for for (int row = 1; row < rowMax-1; row++) { concentrations_t1(row,col) = grid.getConcentrations()(row,col) + timestep / (deltaCol*deltaCol) @@ -328,6 +333,7 @@ void FTCS_2D(Grid &grid, Boundary &bc, double ×tep) { // top without corners / looping over columns // hold row constant at index 0 int row = 0; + #pragma omp parallel for for (int col=1; coldomainCol = domainLength; - this->deltaCol = double(this->domainCol)/this->col; + this->deltaCol = double(this->domainCol)/double(this->col); } void Grid::setDomain(int domainRow, int domainCol) { if (dim != 2) { throw_invalid_argument("Grid is not two dimensional, you should probably use the 1D domain setter!"); } - if (domainRow < 1 || domainCol < 1) { + if (domainRow <= 0 || domainCol <= 0) { throw_invalid_argument("Given domain size is not positive!"); } diff --git a/src/Simulation.cpp b/src/Simulation.cpp index d07065f..5fd552c 100644 --- a/src/Simulation.cpp +++ b/src/Simulation.cpp @@ -154,7 +154,7 @@ void Simulation::run() { } if (approach == FTCS_APPROACH) { - + auto begin = std::chrono::high_resolution_clock::now(); for (int i = 0; i < iterations; i++) { if (console_output == CONSOLE_OUTPUT_VERBOSE && i > 0) { printConcentrationsConsole(); @@ -165,6 +165,9 @@ void Simulation::run() { FTCS(grid, bc, timestep); } + auto end = std::chrono::high_resolution_clock::now(); + auto milliseconds = std::chrono::duration_cast(end - begin); + std::cout << milliseconds.count() << endl; } else if (approach == BTCS_APPROACH) { diff --git a/test/testGrid.cpp b/test/testGrid.cpp index 55143e4..4963976 100644 --- a/test/testGrid.cpp +++ b/test/testGrid.cpp @@ -1,5 +1,21 @@ #include #include +#include + +TEST_CASE("1D Grid, too small length") { + int l = 2; + CHECK_THROWS(Grid(l)); +} + +TEST_CASE("2D Grid, too small side") { + int r = 2; + int c = 4; + CHECK_THROWS(Grid(r, c)); + + r = 4; + c = 2; + CHECK_THROWS(Grid(r, c)); +} TEST_CASE("1D Grid") { int l = 12; @@ -22,21 +38,214 @@ TEST_CASE("1D Grid") { CHECK_THROWS(grid.getDeltaRow()); } - SUBCASE("") { + SUBCASE("setting concentrations") { + // correct concentrations matrix + MatrixXd concentrations = MatrixXd::Constant(1, l, 3); + CHECK_NOTHROW(grid.setConcentrations(concentrations)); + // false concentrations matrix + MatrixXd wConcentrations = MatrixXd::Constant(2, l, 4); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + } + + SUBCASE("setting alpha") { + // correct alpha matrix + MatrixXd alpha = MatrixXd::Constant(1, l, 3); + CHECK_NOTHROW(grid.setAlpha(alpha)); + + CHECK_THROWS(grid.setAlpha(alpha, alpha)); + + grid.setAlpha(alpha); + CHECK_EQ(grid.getAlpha(), alpha); + CHECK_THROWS(grid.getAlphaX()); + CHECK_THROWS(grid.getAlphaY()); + + // false alpha matrix + MatrixXd wAlpha = MatrixXd::Constant(3, l, 2); + CHECK_THROWS(grid.setAlpha(wAlpha)); + } + + SUBCASE("setting domain") { + int d = 8; + // set 1D domain + CHECK_NOTHROW(grid.setDomain(d)); + + // set 2D domain + CHECK_THROWS(grid.setDomain(d, d)); + + grid.setDomain(d); + CHECK_EQ(grid.getDeltaCol(), double(d)/double(l)); + CHECK_THROWS(grid.getDeltaRow()); + + // set too small domain + d = 0; + CHECK_THROWS(grid.setDomain(d)); + d = -2; + CHECK_THROWS(grid.setDomain(d)); } } TEST_CASE("2D Grid quadratic") { - int r = 12; - int c = 12; + int rc = 12; + Grid grid(rc, rc); - // TODO + SUBCASE("correct construction") { + CHECK_EQ(grid.getDim(), 2); + CHECK_THROWS(grid.getLength()); + CHECK_EQ(grid.getCol(), rc); + CHECK_EQ(grid.getRow(), rc); + + CHECK_EQ(grid.getConcentrations().rows(), rc); + CHECK_EQ(grid.getConcentrations().cols(), rc); + CHECK_THROWS(grid.getAlpha()); + + CHECK_EQ(grid.getAlphaX().rows(), rc); + CHECK_EQ(grid.getAlphaX().cols(), rc); + CHECK_EQ(grid.getAlphaY().rows(), rc); + CHECK_EQ(grid.getAlphaY().cols(), rc); + CHECK_EQ(grid.getDeltaRow(), 1); + CHECK_EQ(grid.getDeltaCol(), 1); + } + + SUBCASE("setting concentrations") { + // correct concentrations matrix + MatrixXd concentrations = MatrixXd::Constant(rc, rc, 2); + CHECK_NOTHROW(grid.setConcentrations(concentrations)); + + + // false concentrations matrix + MatrixXd wConcentrations = MatrixXd::Constant(rc, rc+3, 1); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + wConcentrations = MatrixXd::Constant(rc+3, rc, 4); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + wConcentrations = MatrixXd::Constant(rc+2, rc+4, 6); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + } + + SUBCASE("setting alphas") { + // correct alpha matrices + MatrixXd alphax = MatrixXd::Constant(rc, rc, 2); + MatrixXd alphay = MatrixXd::Constant(rc, rc, 4); + CHECK_NOTHROW(grid.setAlpha(alphax, alphay)); + + CHECK_THROWS(grid.setAlpha(alphax)); + + grid.setAlpha(alphax, alphay); + CHECK_EQ(grid.getAlphaX(), alphax); + CHECK_EQ(grid.getAlphaY(), alphay); + CHECK_THROWS(grid.getAlpha()); + + // false alpha matrices + alphax = MatrixXd::Constant(rc+3, rc+1, 3); + CHECK_THROWS(grid.setAlpha(alphax, alphay)); + alphay = MatrixXd::Constant(rc+2, rc+1, 3); + CHECK_THROWS(grid.setAlpha(alphax, alphay)); + } + + SUBCASE("setting domain") { + int dr = 8; + int dc = 9; + + // set 1D domain + CHECK_THROWS(grid.setDomain(dr)); + + // set 2D domain + CHECK_NOTHROW(grid.setDomain(dr, dc)); + + grid.setDomain(dr, dc); + CHECK_EQ(grid.getDeltaCol(), double(dc)/double(rc)); + CHECK_EQ(grid.getDeltaRow(), double(dr)/double(rc)); + + // set too small domain + dr = 0; + CHECK_THROWS(grid.setDomain(dr, dc)); + dr = 8; + dc = 0; + CHECK_THROWS(grid.setDomain(dr, dc)); + dr = -2; + CHECK_THROWS(grid.setDomain(dr, dc)); + } } TEST_CASE("2D Grid non-quadratic") { int r = 12; int c = 15; + Grid grid(r, c); - // TODO + SUBCASE("correct construction") { + CHECK_EQ(grid.getDim(), 2); + CHECK_THROWS(grid.getLength()); + CHECK_EQ(grid.getCol(), c); + CHECK_EQ(grid.getRow(), r); + + CHECK_EQ(grid.getConcentrations().rows(), r); + CHECK_EQ(grid.getConcentrations().cols(), c); + CHECK_THROWS(grid.getAlpha()); + + CHECK_EQ(grid.getAlphaX().rows(), r); + CHECK_EQ(grid.getAlphaX().cols(), c); + CHECK_EQ(grid.getAlphaY().rows(), r); + CHECK_EQ(grid.getAlphaY().cols(), c); + CHECK_EQ(grid.getDeltaRow(), 1); + CHECK_EQ(grid.getDeltaCol(), 1); + } + + SUBCASE("setting concentrations") { + // correct concentrations matrix + MatrixXd concentrations = MatrixXd::Constant(r, c, 2); + CHECK_NOTHROW(grid.setConcentrations(concentrations)); + + + // false concentrations matrix + MatrixXd wConcentrations = MatrixXd::Constant(r, c+3, 6); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + wConcentrations = MatrixXd::Constant(r+3, c, 3); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + wConcentrations = MatrixXd::Constant(r+2, c+4, 2); + CHECK_THROWS(grid.setConcentrations(wConcentrations)); + } + + SUBCASE("setting alphas") { + // correct alpha matrices + MatrixXd alphax = MatrixXd::Constant(r, c, 2); + MatrixXd alphay = MatrixXd::Constant(r, c, 4); + CHECK_NOTHROW(grid.setAlpha(alphax, alphay)); + + CHECK_THROWS(grid.setAlpha(alphax)); + + grid.setAlpha(alphax, alphay); + CHECK_EQ(grid.getAlphaX(), alphax); + CHECK_EQ(grid.getAlphaY(), alphay); + CHECK_THROWS(grid.getAlpha()); + + // false alpha matrices + alphax = MatrixXd::Constant(r+3, c+1, 3); + CHECK_THROWS(grid.setAlpha(alphax, alphay)); + alphay = MatrixXd::Constant(r+2, c+1, 5); + CHECK_THROWS(grid.setAlpha(alphax, alphay)); + } + + SUBCASE("setting domain") { + int dr = 8; + int dc = 9; + + // set 1D domain + CHECK_THROWS(grid.setDomain(dr)); + + // set 2D domain + CHECK_NOTHROW(grid.setDomain(dr, dc)); + + grid.setDomain(dr, dc); + CHECK_EQ(grid.getDeltaCol(), double(dc)/double(c)); + CHECK_EQ(grid.getDeltaRow(), double(dr)/double(r)); + + // set too small domain + dr = 0; + CHECK_THROWS(grid.setDomain(dr, dc)); + dr = 8; + dc = -1; + CHECK_THROWS(grid.setDomain(dr, dc)); + dr = -2; + CHECK_THROWS(grid.setDomain(dr, dc)); + } } \ No newline at end of file