Mercurial > repos > public > sbplib

changeset 154:c7b2f645101f feature/grids

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Added classes and functions for Cartesian and equidistant grids.
author Jonatan Werpers <jonatan@werpers.com>
date Thu, 18 Feb 2016 16:46:02 +0100
parents 7aee9eba3bb8
children cf748f33cd5b
files +grid/Cartesian.m +grid/CartesianTest.m +grid/equidistant.m +grid/equidistantTest.m
diffstat 4 files changed, 253 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+grid/Cartesian.m	Thu Feb 18 16:46:02 2016 +0100
@@ -0,0 +1,78 @@
+classdef Cartesian < grid.Grid
+    properties
+        n % Number of points in the grid
+        d % Number of dimensions
+        m % Number of points in each direction
+        x % Cell array of vectors with node placement for each dimension.
+    end
+
+    % General d dimensional grid with n points
+    methods
+        % Creates a cartesian grid given vectors conatining the coordinates
+        % in each direction
+        function obj = Cartesian(varargin)
+            obj.d = length(varargin);
+            for i = 1:obj.d
+                obj.x{i} = varargin{i};
+                obj.m(i) = length(varargin{i});
+            end
+            obj.n = prod(obj.m);
+            if obj.n == 0
+                error('grid:Cartesian:EmptyGrid','Input parameter gives an empty grid.')
+            end
+        end
+        % n returns the number of points in the grid
+        function o = N(obj)
+            o = obj.n;
+        end
+
+        % d returns the spatial dimension of the grid
+        function o = D(obj)
+            o = obj.d;
+        end
+
+        % points returns a n x d matrix containing the coordianets for all points.
+        % points are ordered according to the kronecker product with X*Y*Z
+        function X = points(obj)
+            X = zeros(obj.n, obj.d);
+
+            for i = 1:obj.d
+                if iscolumn(obj.x{i})
+                    c = obj.x{i};
+                else
+                    c = obj.x{i}';
+                end
+
+                m_before = prod(obj.m(1:i-1));
+                m_after = prod(obj.m(i+1:end));
+
+                X(:,i) = kr(ones(m_before,1),c,ones(m_after,1));
+            end
+        end
+
+        % matrices returns a cell array with coordinates in matrix form.
+        % For 2d case these will have to be transposed to work with plotting routines.
+        function X = matrices(obj)
+
+            if obj.d == 1 % There is no 1d matrix data type in matlab, handle special case
+                X{1} = reshape(obj.x{1}, [obj.m 1]);
+                return
+            end
+
+            X = cell(1,obj.d);
+            for i = 1:obj.d
+                s = ones(1,obj.d);
+                s(i) = obj.m(i);
+
+                t = reshape(obj.x{i},s);
+
+                s = obj.m;
+                s(i) = 1;
+                X{i} = repmat(t,s);
+            end
+        end
+
+        % coordVectors()
+        % coordMatrices()
+    end
+end
\ No newline at end of file
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+grid/CartesianTest.m	Thu Feb 18 16:46:02 2016 +0100
@@ -0,0 +1,98 @@
+function tests = CartesianTest()
+    tests = functiontests(localfunctions);
+end
+
+
+function testWarningEmptyGrid(testCase)
+    in  = {
+        {[]},
+        {[],[1]},
+        {[1],[2], []},
+    };
+
+    for i = 1:length(in)
+        testCase.verifyError(@()grid.Cartesian(in{i}{:}),'grid:Cartesian:EmptyGrid');
+    end
+end
+
+function testN(testCase)
+    in  = {
+        {[1 2 3]},
+        {[1 2 3],[1 2]},
+        {[1 2 3],[1 2 3]},
+        {[1 2 3],[1 2 3], [1]},
+        {[1 2 3],[1 2 3], [1 3 4]},
+    };
+
+    out = [3,6,9,9,27];
+
+    for i = 1:length(in)
+        g = grid.Cartesian(in{i}{:});
+        testCase.verifyEqual(g.N(),out(i));
+    end
+end
+
+
+function testD(testCase)
+    in  = {
+        {[1 2 3]},
+        {[1 2 3],[1 2]},
+        {[1 2 3],[1 2 3]},
+        {[1 2 3],[1 2 3], [1]},
+        {[1 2 3],[1 2 3], [1 3 4]},
+    };
+
+    out = [1,2,2,3,3];
+
+    for i = 1:length(in)
+        g = grid.Cartesian(in{i}{:});
+        testCase.verifyEqual(g.D(),out(i));
+    end
+end
+
+
+function testPoints(testCase)
+    in  = {
+        {[1 2]},
+        {[1 2],[3 4]},
+        {[1 2],[3 4], [5 6]},
+    };
+
+    out = {
+        [[1; 2]],
+        [[1; 1; 2; 2],[3; 4; 3; 4]],
+        [[1; 1; 1; 1; 2; 2; 2; 2],[3; 3; 4; 4; 3; 3; 4; 4],[ 5; 6; 5; 6; 5; 6; 5; 6]],
+    };
+
+    for i = 1:length(in)
+        g = grid.Cartesian(in{i}{:});
+        testCase.verifyEqual(g.points(),out{i});
+    end
+end
+
+function testMatrices(testCase)
+    in  = {
+        {[1 2]},
+        {[1 2],[3 4]},
+        {[1 2],[3 4], [5 6]},
+    };
+
+    out{1}{1} = [1; 2];
+
+    out{2}{1} = [1, 1; 2, 2];
+    out{2}{2} = [3, 4; 3, 4];
+
+    out{3}{1}(:,:,1) = [1, 1; 2, 2];
+    out{3}{1}(:,:,2) = [1, 1; 2, 2];
+
+    out{3}{2}(:,:,1) = [3, 4; 3, 4];
+    out{3}{2}(:,:,2) = [3, 4; 3, 4];
+
+    out{3}{3}(:,:,1) = [5, 5; 5, 5];
+    out{3}{3}(:,:,2) = [6, 6; 6, 6];
+
+    for i = 1:length(in)
+        g = grid.Cartesian(in{i}{:});
+        testCase.verifyEqual(g.matrices(),out{i});
+    end
+end
\ No newline at end of file
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+grid/equidistant.m	Thu Feb 18 16:46:02 2016 +0100
@@ -0,0 +1,27 @@
+% Creates a cartesian grid of dimension length(m).
+% over the doman xlim, ylim, ...
+% Examples:
+%   g = grid.equidistant([mx, my], xlim, ylim)
+%   g = grid.equidistant([10, 15], {0,1}, {0,2})
+function g = equidistant(m, varargin)
+    if length(m) ~= length(varargin)
+        error('grid:equidistant:NonMatchingParameters','The number of provided dimensions do not match.')
+    end
+
+    for i = 1:length(m)
+        if ~iscell(varargin{i}) || numel(varargin{i}) ~= 2
+           error('grid:equidistant:InvalidLimits','The limits should be cell arrays with 2 elements.');
+        end
+
+        if varargin{i}{1} > varargin{i}{2}
+            error('grid:equidistant:InvalidLimits','The elements of the limit must be increasing.');
+        end
+    end
+
+    X = {};
+    for i = 1:length(m)
+        X{i} = util.get_grid(varargin{i}{:},m(i));
+    end
+
+    g = grid.Cartesian(X{:});
+end
\ No newline at end of file
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+grid/equidistantTest.m	Thu Feb 18 16:46:02 2016 +0100
@@ -0,0 +1,50 @@
+function tests = equidistantTest()
+    tests = functiontests(localfunctions);
+end
+
+
+function testErrorInvalidLimits(testCase)
+     in  = {
+        {10,{1}},
+        {10,[0,1]},
+        {[10,10],{0,1},{1}},
+        {[10,10],{1},{1,0}},
+        {10,{1,0}},
+        {[10, 5],{1,0}, {0,-1}},
+    };
+
+    for i = 1:length(in)
+        testCase.verifyError(@()grid.equidistant(in{i}{:}),'grid:equidistant:InvalidLimits',sprintf('in(%d) = %s',i,toString(in{i})));
+    end
+end
+
+function testErrorNonMatchingParam(testCase)
+    in  = {
+        {[],{1}},
+        {[],{1},{0,1}},
+        {[5,5],{0,1},{0,1},{0,1}},
+        {[5,5,4],{0,1},{0,1}},
+    };
+
+    for i = 1:length(in)
+        testCase.verifyError(@()grid.equidistant(in{i}{:}),'grid:equidistant:NonMatchingParameters',sprintf('in(%d) = %s',i,toString(in{i})));
+    end
+end
+
+
+function testCompiles(testCase)
+    in  = {
+        {5, {0,1}},
+        {[3 3],{0,1},{0,2}},
+    };
+
+    out = {
+        [[0; 0.25; 0.5; 0.75; 1]],
+        [[0; 0; 0; 0.5; 0.5; 0.5; 1; 1; 1;],[0; 1; 2; 0; 1; 2; 0; 1; 2;]],
+    };
+
+    for i = 1:length(in)
+        g = grid.equidistant(in{i}{:});
+        testCase.verifyEqual(g.points(),out{i});
+    end
+end
\ No newline at end of file