From c0c855cc98c587f812b0d670f0b8eed33c789429 Mon Sep 17 00:00:00 2001
From: Kilian Weishaupt <kilian.weishaupt@iws.uni-stuttgart.de>
Date: Fri, 17 May 2019 15:08:13 +0200
Subject: [PATCH] [boxlocalAssembler] Fix Dirichlet BCs for explicit assembly

* only write to matrix entries that actually exist
---
 dumux/assembly/boxlocalassembler.hh | 23 +++++++++++------------
 1 file changed, 11 insertions(+), 12 deletions(-)

diff --git a/dumux/assembly/boxlocalassembler.hh b/dumux/assembly/boxlocalassembler.hh
index 81883e967c..29188da305 100644
--- a/dumux/assembly/boxlocalassembler.hh
+++ b/dumux/assembly/boxlocalassembler.hh
@@ -131,12 +131,12 @@ public:
                                    const auto pvIdx)
         {
             res[scvI.dofIndex()][eqIdx] = this->curElemVolVars()[scvI].priVars()[pvIdx] - dirichletValues[pvIdx];
-            for (const auto& scvJ : scvs(this->fvGeometry()))
-            {
-                jac[scvI.dofIndex()][scvJ.dofIndex()][eqIdx] = 0.0;
-                if (scvI.localDofIndex() == scvJ.localDofIndex())
-                    jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
-            }
+
+            auto& row = jac[scvI.dofIndex()];
+            for (auto col = row.begin(); col != row.end(); ++col)
+                row[col.index()][eqIdx] = 0.0;
+
+            jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
 
             // if a periodic dof has Dirichlet values also apply the same Dirichlet values to the other dof
             if (this->assembler().fvGridGeometry().dofOnPeriodicBoundary(scvI.dofIndex()))
@@ -166,12 +166,11 @@ public:
                                    const auto eqIdx,
                                    const auto pvIdx)
         {
-            for (const auto& scvJ : scvs(this->fvGeometry()))
-            {
-                jac[scvI.dofIndex()][scvJ.dofIndex()][eqIdx] = 0.0;
-                if (scvI.localDofIndex() == scvJ.localDofIndex())
-                    jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
-            }
+            auto& row = jac[scvI.dofIndex()];
+            for (auto col = row.begin(); col != row.end(); ++col)
+                row[col.index()][eqIdx] = 0.0;
+
+            jac[scvI.dofIndex()][scvI.dofIndex()][eqIdx][pvIdx] = 1.0;
         };
 
         this->asImp_().evalDirichletBoundaries(applyDirichlet);
-- 
GitLab