--- a/+multiblock/DiffOp.m	Fri Jul 27 10:31:51 2018 -0700
+++ b/+multiblock/DiffOp.m	Fri Jul 27 10:39:12 2018 -0700
@@ -201,33 +201,8 @@
             [blockClosure, blockPenalty] = obj.diffOps{I}.boundary_condition(name, type);
 
             % Expand to matrix for full domain.
-            div = obj.blockmatrixDiv;
-            if ~iscell(blockClosure)
-                temp = blockmatrix.zero(div);
-                temp{I,I} = blockClosure;
-                closure = blockmatrix.toMatrix(temp);
-            else
-                for i = 1:length(blockClosure)
-                    temp = blockmatrix.zero(div);
-                    temp{I,I} = blockClosure{i};
-                    closure{i} = blockmatrix.toMatrix(temp);
-                end
-            end
-
-            if ~iscell(blockPenalty)
-                div{2} = size(blockPenalty, 2); % Penalty is a column vector
-                p = blockmatrix.zero(div);
-                p{I} = blockPenalty;
-                penalty = blockmatrix.toMatrix(p);
-            else
-                % TODO: used by beam equation, should be eliminated. SHould only set one BC per call
-                for i = 1:length(blockPenalty)
-                    div{2} = size(blockPenalty{i}, 2); % Penalty is a column vector
-                    p = blockmatrix.zero(div);
-                    p{I} = blockPenalty{i};
-                    penalty{i} = blockmatrix.toMatrix(p);
-                end
-            end
+            closure = multiblock.local2globalClosure(blockClosure, obj.blockmatrixDiv, I);
+            penalty = multiblock.local2globalPenalty(blockPenalty, obj.blockmatrixDiv, I);
         end
 
         function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+multiblock/local2globalClosure.m	Fri Jul 27 10:39:12 2018 -0700
@@ -0,0 +1,11 @@
+
+% Takes the local closure for ice or water and turns it into a closure for the whole system
+%   local -- The local closure
+%   div   -- block matrix division for the diffOp
+%   I     -- Index of blockmatrix block
+function closure = local2globalClosure(local, div, I)
+    closure_bm = blockmatrix.zero(div);
+    closure_bm{I,I} = local;
+
+    closure = blockmatrix.toMatrix(closure_bm);
+end

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+multiblock/local2globalPenalty.m	Fri Jul 27 10:39:12 2018 -0700
@@ -0,0 +1,11 @@
+% Takes the local penalty for ice or water and turns it into a penalty for the whole system
+%   local -- The local penalty
+%   div   -- block matrix division for the diffOp
+%   I     -- Index of blockmatrix block
+function penalty = local2globalPenalty(local, div, I)
+    penaltyDiv = {div{1}, size(local,2)};
+    penalty_bm = blockmatrix.zero(penaltyDiv);
+    penalty_bm{I,1} = local;
+
+    penalty = blockmatrix.toMatrix(penalty_bm);
+end

--- a/+scheme/bcSetup.m	Fri Jul 27 10:31:51 2018 -0700
+++ b/+scheme/bcSetup.m	Fri Jul 27 10:39:12 2018 -0700
@@ -59,6 +59,17 @@
     S = @S_fun;
 end
 
+% Borde man använda eval on här??
+% Borde man dela upp bcSetup i bcSetupSymbolic(name?) och bcSetupGridData
+% och sen skriva en wrapper som sorterar och wrappar de två andra??
+
+% Borde man ha en separat funktion för closure penalty generering
+% och en separat för att bygga ihop penaltyn med data?
+
+% Erbjuda en separat function for att validera en bc specifikation?
+%  alltid kräva alla fields?
+%  literal struct improvement?
+
 function verifyBcFormat(bcs, diffOp)
     for i = 1:length(bcs)
         assertType(bcs{i}, 'struct');