From 61d6cfc5cb355f3fbaf66b8fcf65e5ba1e9f9330 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Max=20L=C3=BCbke?= <mluebke@gfz-potsdam.de>
Date: Tue, 19 Dec 2023 16:04:57 +0100
Subject: [PATCH] fix: remove factor 2 in 'middle' coefficient of explicit part
 (closed)

test: add test case with homogeneous field and same values as const
boundary condition
---
 include/tug/Core/BTCS.hpp |  2 +-
 test/testSimulation.cpp   | 32 ++++++++++++++++++++++++++++++++
 2 files changed, 33 insertions(+), 1 deletion(-)

diff --git a/include/tug/Core/BTCS.hpp b/include/tug/Core/BTCS.hpp
index b844ba4..4ee01ed 100644
--- a/include/tug/Core/BTCS.hpp
+++ b/include/tug/Core/BTCS.hpp
@@ -140,7 +140,7 @@ constexpr T calcExplicitConcentrationsBoundaryConstant(T conc_center, T conc_bc,
                                                        T alpha_neighbor, T sy) {
   return sy * calcAlphaIntercell(alpha_center, alpha_neighbor) * conc_center +
          (1 - sy * (calcAlphaIntercell(alpha_center, alpha_center) +
-                    2 * alpha_center)) *
+                    alpha_center)) *
              conc_center +
          sy * alpha_center * conc_bc;
 }
diff --git a/test/testSimulation.cpp b/test/testSimulation.cpp
index ca6888d..2bd831e 100644
--- a/test/testSimulation.cpp
+++ b/test/testSimulation.cpp
@@ -1,5 +1,6 @@
 #include "TestUtils.hpp"
 
+#include <Eigen/src/Core/Matrix.h>
 #include <doctest/doctest.h>
 #include <stdio.h>
 #include <string>
@@ -119,3 +120,34 @@ TEST_CASE("Simulation environment") {
     CHECK_THROWS(sim.setTimestep(-0.3));
   }
 }
+
+TEST_CASE("Closed Boundaries - no change expected") {
+
+  constexpr std::uint32_t nrows = 5;
+  constexpr std::uint32_t ncols = 5;
+
+  tug::Grid64 grid(nrows, ncols);
+
+  auto concentrations = Eigen::MatrixXd::Constant(nrows, ncols, 1.0);
+  auto alphax = Eigen::MatrixXd::Constant(nrows, ncols, 1E-5);
+  auto alphay = Eigen::MatrixXd::Constant(nrows, ncols, 1E-5);
+
+  grid.setConcentrations(concentrations);
+  grid.setAlpha(alphax, alphay);
+
+  tug::Boundary bc(grid);
+  bc.setBoundarySideConstant(tug::BC_SIDE_LEFT, 1.0);
+  bc.setBoundarySideConstant(tug::BC_SIDE_RIGHT, 1.0);
+  bc.setBoundarySideConstant(tug::BC_SIDE_TOP, 1.0);
+  bc.setBoundarySideConstant(tug::BC_SIDE_BOTTOM, 1.0);
+
+  tug::Simulation<double> sim(grid, bc);
+  sim.setTimestep(1);
+  sim.setIterations(1);
+
+  MatrixXd input_values(concentrations);
+  sim.run();
+
+  CHECK(checkSimilarityV2(input_values, grid.getConcentrations(), 1E-12) ==
+        true);
+}