Mercurial > repos > public > sbplib

diff +scheme/Hypsyst2d.m @ 997:78db023a7fe3 feature/getBoundaryOp

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add getBoundaryOperator method to all 2d schemes, except obsolete Wave2dCurve and ElastiCurve, which needs a big makeover.
author Martin Almquist <malmquist@stanford.edu>
date Sat, 12 Jan 2019 11:57:50 -0800
parents 706d1c2b4199
children 8d73fcdb07a5
line wrap: on
line diff
--- a/+scheme/Hypsyst2d.m	Tue Jan 08 11:51:24 2019 +0100
+++ b/+scheme/Hypsyst2d.m	Sat Jan 12 11:57:50 2019 -0800
@@ -186,10 +186,10 @@
             params = obj.params;
             x = obj.x;
             y = obj.y;
+            e_ = obj.getBoundaryOperator('e', boundary);
 
             switch boundary
                 case {'w','W','west'}
-                    e_ = obj.e_w;
                     mat = obj.A;
                     boundPos = 'l';
                     Hi = obj.Hxi;
@@ -197,7 +197,6 @@
                     L = obj.evaluateCoefficientMatrix(L,x(1),y);
                     side = max(length(y));
                 case {'e','E','east'}
-                    e_ = obj.e_e;
                     mat = obj.A;
                     boundPos = 'r';
                     Hi = obj.Hxi;
@@ -205,7 +204,6 @@
                     L = obj.evaluateCoefficientMatrix(L,x(end),y);
                     side = max(length(y));
                 case {'s','S','south'}
-                    e_ = obj.e_s;
                     mat = obj.B;
                     boundPos = 'l';
                     Hi = obj.Hyi;
@@ -213,7 +211,6 @@
                     L = obj.evaluateCoefficientMatrix(L,x,y(1));
                     side = max(length(x));
                 case {'n','N','north'}
-                    e_ = obj.e_n;
                     mat = obj.B;
                     boundPos = 'r';
                     Hi = obj.Hyi;
@@ -297,5 +294,56 @@
             signVec = [sum(poseig),sum(zeroeig),sum(negeig)];
         end
 
+        % Returns the boundary operator op for the boundary specified by the string boundary.
+        % op        -- string or a cell array of strings
+        % boundary  -- string
+        function varargout = getBoundaryOperator(obj, op, boundary)
+
+            if ~iscell(op)
+                op = {op};
+            end
+
+            for i = 1:numel(op)
+                switch op{i}
+                case 'e'
+                    switch boundary
+                    case 'w'
+                        e = obj.e_w;
+                    case 'e'
+                        e = obj.e_e;
+                    case 's'
+                        e = obj.e_s;
+                    case 'n'
+                        e = obj.e_n;
+                    otherwise
+                        error('No such boundary: boundary = %s',boundary);
+                    end
+                    varargout{i} = e;
+                end
+            end
+        end
+
+        % Returns square boundary quadrature matrix, of dimension
+        % corresponding to the number of boundary points
+        %
+        % boundary -- string
+        function H_b = getBoundaryQuadrature(obj, boundary)
+
+            e = obj.getBoundaryOperator('e', boundary);
+
+            switch boundary
+                case 'w'
+                    H_b = inv(e'*obj.Hyi*e);
+                case 'e'
+                    H_b = inv(e'*obj.Hyi*e);
+                case 's'
+                    H_b = inv(e'*obj.Hxi*e);
+                case 'n'
+                    H_b = inv(e'*obj.Hxi*e);
+                otherwise
+                    error('No such boundary: boundary = %s',boundary);
+            end
+        end
+
     end
 end
\ No newline at end of file