Mercurial > repos > public > sbplib

changeset 191:7c1d3fc33f90 feature/grids

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Added methods for returning boundary names and boundary coordinates from Cartesian and Curvilinear grids.
author Jonatan Werpers <jonatan@werpers.com>
date Mon, 04 Apr 2016 18:23:50 +0200
parents 2c2ba1f3bbe3
children 3cedd5a596bb
files +grid/Cartesian.m +grid/CartesianTest.m +grid/Curvilinear.m +grid/CurvilinearTest.m
diffstat 4 files changed, 265 insertions(+), 11 deletions(-) [+]
line wrap: on
line diff
--- a/+grid/Cartesian.m	Mon Mar 21 16:33:49 2016 +0100
+++ b/+grid/Cartesian.m	Mon Apr 04 18:23:50 2016 +0200
@@ -125,13 +125,65 @@
 
         % Return the names of all boundaries in this grid.
         function bs = getBoundaryNames(obj)
-            bs = [];
+            switch obj.D()
+                case 1
+                    bs = {'l', 'r'};
+                case 2
+                    bs = {'w', 'e', 's', 'n'};
+                case 3
+                    bs = {'w', 'e', 's', 'n', 'd', 'u'};
+                otherwise
+                    error('not implemented');
+            end
         end
 
         % Return coordinates for the given boundary
-        function b = getBoundary(obj, name)
-            b = [];
+        function X = getBoundary(obj, name)
+            % In what dimension is the boundary?
+            switch name
+                case {'l', 'r', 'w', 'e'}
+                    D = 1;
+                case {'s', 'n'}
+                    D = 2;
+                case {'d', 'u'}
+                    D = 3;
+                otherwise
+                    error('not implemented');
+            end
+
+            % At what index is the boundary?
+            switch name
+                case {'l', 'w', 's', 'd'}
+                    index = 1;
+                case {'r', 'e', 'n', 'u'}
+                    index = obj.m(D);
+                otherwise
+                    error('not implemented');
+            end
+
+
+
+            I = cell(1, obj.d);
+            for i = 1:obj.d
+                if i == D
+                    I{i} = index;
+                else
+                    I{i} = ':';
+                end
+            end
+
+            % Calculate size of result:
+            m = obj.m;
+            m(D) = [];
+            N = prod(m);
+
+            X = zeros(N, obj.d);
+
+            coordMat = obj.matrices();
+            for i = 1:length(coordMat)
+                Xtemp = coordMat{i}(I{:});
+                X(:,i) = reshapeRowMaj(Xtemp, [N,1]);
+            end
         end
-
     end
 end
\ No newline at end of file
--- a/+grid/CartesianTest.m	Mon Mar 21 16:33:49 2016 +0100
+++ b/+grid/CartesianTest.m	Mon Apr 04 18:23:50 2016 +0200
@@ -197,9 +197,122 @@
 
 
 function testGetBoundaryNames(testCase)
-    testCase.verifyFail();
+    in = {
+        {[1 2 3]},
+        {[1 2 3], [4 5]},
+        {[1 2 3], [4 5], [6 7 8]},
+    };
+
+    out = {
+        {'l', 'r'},
+        {'w', 'e', 's', 'n'},
+        {'w', 'e', 's', 'n', 'd', 'u'},
+    };
+
+    for i = 1:length(in)
+        g = grid.Cartesian(in{i}{:});
+        testCase.verifyEqual(g.getBoundaryNames(), out{i});
+    end
 end
 
 function testGetBoundary(testCase)
-    testCase.verifyFail();
+    grids = {
+        {[1 2 3]},
+        {[1 2 3], [4 5]},
+        {[1 2 3], [4 5], [6 7 8]},
+    };
+
+    boundaries = {
+        {'l', 'r'},
+        {'w', 'e', 's', 'n'},
+        {'w', 'e', 's', 'n', 'd', 'u'},
+    };
+
+
+    % 1d
+    out{1,1} = 1;
+    out{1,2} = 3;
+
+    % 2d
+    out{2,1} = [
+        1,4;
+        1,5;
+    ];
+    out{2,2} = [
+        3,4;
+        3,5;
+    ];
+    out{2,3} = [
+        1,4;
+        2,4;
+        3,4;
+    ];
+    out{2,4} = [
+        1,5;
+        2,5;
+        3,5;
+    ];
+
+    % 3d
+    out{3,1} = [
+        1,4,6;
+        1,4,7;
+        1,4,8;
+        1,5,6;
+        1,5,7;
+        1,5,8;
+    ];
+    out{3,2} = [
+        3,4,6;
+        3,4,7;
+        3,4,8;
+        3,5,6;
+        3,5,7;
+        3,5,8;
+    ];
+    out{3,3} = [
+        1,4,6;
+        1,4,7;
+        1,4,8;
+        2,4,6;
+        2,4,7;
+        2,4,8;
+        3,4,6;
+        3,4,7;
+        3,4,8;
+    ];
+    out{3,4} = [
+        1,5,6;
+        1,5,7;
+        1,5,8;
+        2,5,6;
+        2,5,7;
+        2,5,8;
+        3,5,6;
+        3,5,7;
+        3,5,8;
+    ];
+    out{3,5} = [
+        1,4,6;
+        1,5,6;
+        2,4,6;
+        2,5,6;
+        3,4,6;
+        3,5,6;
+    ];
+    out{3,6} = [
+        1,4,8;
+        1,5,8;
+        2,4,8;
+        2,5,8;
+        3,4,8;
+        3,5,8;
+    ];
+
+    for ig = 1:length(grids)
+        g = grid.Cartesian(grids{ig}{:});
+        for ib = 1:length(boundaries{ig})
+            testCase.verifyEqual(g.getBoundary(boundaries{ig}{ib}), out{ig,ib});
+        end
+    end
 end
--- a/+grid/Curvilinear.m	Mon Mar 21 16:33:49 2016 +0100
+++ b/+grid/Curvilinear.m	Mon Apr 04 18:23:50 2016 +0200
@@ -83,12 +83,60 @@
 
         % Return the names of all boundaries in this grid.
         function bs = getBoundaryNames(obj)
-            bs = [];
+            bs = obj.logic.getBoundaryNames();
         end
 
         % Return coordinates for the given boundary
-        function b = getBoundary(obj, name)
-            b = [];
+        function X = getBoundary(obj, name)
+              % In what dimension is the boundary?
+            switch name
+                case {'l', 'r', 'w', 'e'}
+                    D = 1;
+                case {'s', 'n'}
+                    D = 2;
+                case {'d', 'u'}
+                    D = 3;
+                otherwise
+                    error('not implemented');
+            end
+
+            % At what index is the boundary?
+            switch name
+                case {'l', 'w', 's', 'd'}
+                    index = 1;
+                case {'r', 'e', 'n', 'u'}
+                    index = obj.logic.m(D);
+                otherwise
+                    error('not implemented');
+            end
+
+
+
+            I = cell(1, obj.D);
+            for i = 1:obj.D
+                if i == D
+                    I{i} = index;
+                else
+                    I{i} = ':';
+                end
+            end
+
+            % Calculate size of result:
+            m = obj.logic.m;
+            m(D) = [];
+            N = prod(m);
+
+            X = zeros(N, obj.D);
+
+            p = obj.points;
+            for i = 1:obj.D()
+                coordMat{i} = reshapeRowMaj(p(:,i), obj.logic.m);
+            end
+
+            for i = 1:length(coordMat)
+                Xtemp = coordMat{i}(I{:});
+                X(:,i) = reshapeRowMaj(Xtemp, [N,1]);
+            end
         end
     end
 end
--- a/+grid/CurvilinearTest.m	Mon Mar 21 16:33:49 2016 +0100
+++ b/+grid/CurvilinearTest.m	Mon Apr 04 18:23:50 2016 +0200
@@ -81,10 +81,51 @@
 end
 
 function testGetBoundaryNames(testCase)
-    testCase.verifyFail();
+    in = {
+        {{1:10}, @(x) exp(x)},
+        {{1:10,1:6}, @(x,y) [exp(x+y); exp(x-y)]},
+        {{1:10,1:5,1:7}, @(x,y,z)[exp(x+y+z); exp(x-y-z); 2+x+y-z]},
+    };
+
+    out = {
+        {'l', 'r'},
+        {'w', 'e', 's', 'n'},
+        {'w', 'e', 's', 'n', 'd', 'u'},
+    };
+
+    for i = 1:length(in)
+        g = grid.Curvilinear(in{i}{2},in{i}{1}{:});
+        testCase.verifyEqual(g.getBoundaryNames(), out{i});
+    end
 end
 
 function testGetBoundary(testCase)
-    testCase.verifyFail();
+    grids = {
+        {{1:10}, @(x) exp(x)},
+        {{1:10,1:6}, @(x,y) [exp(x+y); exp(x-y)]},
+        {{1:10,1:5,1:7}, @(x,y,z)[exp(x+y+z); exp(x-y-z); 2+x+y-z]},
+    };
+
+    boundaries = {
+        {'l', 'r'},
+        {'w', 'e', 's', 'n'},
+        {'w', 'e', 's', 'n', 'd', 'u'},
+    };
+
+
+    for ig = 1:length(grids)
+        g = grid.Curvilinear(grids{ig}{2},grids{ig}{1}{:});
+
+        logicalGrid = grid.Cartesian(grids{ig}{1}{:});
+
+        for ib = 1:length(boundaries{ig})
+
+            logicalBoundary = logicalGrid.getBoundary(boundaries{ig}{ib});
+
+            x = num2cell(logicalBoundary',2);
+            expectedBoundary = grids{ig}{2}(x{:})';
+            testCase.verifyEqual(g.getBoundary(boundaries{ig}{ib}), expectedBoundary);
+        end
+    end
 end