--- a/+scheme/Elastic2dVariableAnisotropic.m	Thu Sep 05 16:41:49 2019 -0700
+++ b/+scheme/Elastic2dVariableAnisotropic.m	Thu Sep 05 17:10:58 2019 -0700
@@ -367,6 +367,12 @@
         %       data                is a function returning the data that should be applied at the boundary.
         %       neighbour_scheme    is an instance of Scheme that should be interfaced to.
         %       neighbour_boundary  is a string specifying which boundary to interface to.
+
+        % For displacement bc:
+        % bc = {comp, 'd', dComps},
+        % where
+        % dComps = vector of components with displacement BC. Default: 1:dim.
+        % In this way, we can specify one BC at a time even though the SATs depend on all BC.
         function [closure, penalty] = boundary_condition(obj, boundary, bc, tuning)
             default_arg('tuning', 1.0);
 
@@ -416,7 +422,7 @@
                 % X = Y + Z.
 
                 % Nonsymmetric part goes on all components to
-                % Yield traction in discrete energy rate
+                % yield traction in discrete energy rate
                 for i = 1:dim
                     Y = T{j,i}';
                     X = e*Y;
@@ -424,7 +430,8 @@
                     penalty = penalty - E{i}*RHOi*Hi*X'*e*H_gamma;
                 end
 
-                % Symmetric part only required on components with displacement BC
+                % Symmetric part only required on components with displacement BC.
+                % (Otherwise it's not symmetric.)
                 for i = dComps
                     Z = sparse(m_tot, m_tot);
                     for l = 1:dim