diff --git a/+quantity/Discrete.m b/+quantity/Discrete.m
index afc6b5ac86c5ce8317fa02ecedd2fba50986e464..b8e78f9c255f67c9e70ba5c2d0a7e712f265a7cd 100644
--- a/+quantity/Discrete.m
+++ b/+quantity/Discrete.m
@@ -329,55 +329,101 @@ classdef  (InferiorClasses = {?quantity.Symbolic, ?quantity.Operator}) Discrete
 		end
 		
 		function c = horzcat(a, varargin)
-			% Before Concatenate, check if all quantites have same grid.
-			% Furthermore, copies of the input quantites are used for
-			% concatenation.
-			if nargin == 1
-				objCell = {a};
-			else
-				objCell = {a, varargin{:}};
-			end
-			[fineGrid, fineGridName] = getFinestGrid(objCell{:});
-			for it = 1 : nargin
-				objCopyTemp = copy(objCell{it});
-				objCell{it} = objCopyTemp.changeGrid(fineGrid, fineGridName);
-			end
-			assertSameGrid(objCell{:});
-			 c = builtin('horzcat', objCell{:});
+			%HORZCAT Horizontal concatenation.
+			%   [A B] is the horizontal concatenation of objects A and B from the
+			%   class quantity.Discrete. A and B must have the same
+			%   number of rows and the same grid. [A,B] is the same thing.
+			%   Any number of matrices can be concatenated within one pair
+			%   of brackets. Horizontal and vertical concatenation can be
+			%   combined together as in [1 2;3 4].
+			%
+			%   [A B; C] is allowed if the number of rows of A equals the
+			%   number of rows of B and the number of columns of A plus the
+			%   number of columns of B equals the number of columns of C.
+			%   The matrices in a concatenation expression can themselves
+			%   by formed via a concatenation as in [A B;[C D]].  These
+			%   rules generalize in a hopefully obvious way to allow fairly
+			%   complicated constructions.
+			%
+			%   N-D arrays are concatenated along the second dimension. The
+			%   first and remaining dimensions must match.
+			%
+			%   C = HORZCAT(A,B) is called for the syntax '[A  B]' when A
+			%   or B is an object.
+			%
+			%   Y = HORZCAT(X1,X2,X3,...) is called for the syntax '[X1 X2
+			%   X3 ...]' when any of X1, X2, X3, etc. is an object.
+			%
+			%	See also HORZCAT, CAT.
+			 c = a.cat(2, varargin{:});
 		end
 		function c = vertcat(a, varargin)
-			% Before Concatenate, check if all quantites have same grid.
-			% Furthermore, copies of the input quantites are used for
-			% concatenation.
-			if nargin == 1
-				objCell = {a};
-			else
-				objCell = {a, varargin{:}};
-			end
-			[fineGrid, fineGridName] = getFinestGrid(objCell{:});
-			for it = 1 : nargin
-				objCopyTemp = copy(objCell{it});
-				objCell{it} = objCopyTemp.changeGrid(fineGrid, fineGridName);
-			end
-			assertSameGrid(objCell{:});
-			c = builtin('vertcat', objCell{:});
-		end
-		function c = cat(a, varargin)
-			% Before Concatenate, check if all quantites have same grid.
-			% Furthermore, copies of the input quantites are used for
-			% concatenation.
+			%VERTCAT Vertical concatenation.
+			%   [A;B] is the vertical concatenation of objects A and B from
+			%   the class quantity.Discrete. A and B must have the same
+			%   number of columns and the same grid. Any number of matrices
+			%   can be concatenated within one pair of brackets.
+			%   Horizontal and vertical concatenation can be combined
+			%   together as in [1 2;3 4].
+			%
+			%   [A B; C] is allowed if the number of rows of A equals the
+			%   number of rows of B and the number of columns of A plus the
+			%   number of columns of B equals the number of columns of C.
+			%   The matrices in a concatenation expression can themselves
+			%   by formed via a concatenation as in [A B;[C D]].  These
+			%   rules generalize in a hopefully obvious way to allow fairly
+			%   complicated constructions.
+			%
+			%   N-D arrays are concatenated along the first dimension. The
+			%   remaining dimensions must match.
+			%
+			%   C = VERTCAT(A,B) is called for the syntax '[A; B]' when A
+			%   or B is an object.
+			%
+			%   Y = VERTCAT(X1,X2,X3,...) is called for the syntax '[X1;
+			%   X2; X3; ...]' when any of X1, X2, X3, etc. is an object.
+			%
+			%   See also HORZCAT, CAT.
+			c = a.cat(1, varargin{:});
+		end
+		function c = cat(a, dim, varargin)
+			%CAT Concatenate arrays.
+			%   CAT(DIM,A,B) concatenates the arrays of objects A and B
+			%   from the class quantity.Discrete along the dimension DIM.
+			%   CAT(2,A,B) is the same as [A,B].
+			%   CAT(1,A,B) is the same as [A;B].
+			%
+			%   B = CAT(DIM,A1,A2,A3,A4,...) concatenates the input
+			%   arrays A1, A2, etc. along the dimension DIM.
+			%
+			%   When used with comma separated list syntax, CAT(DIM,C{:}) or 
+			%   CAT(DIM,C.FIELD) is a convenient way to concatenate a cell or
+			%   structure array containing numeric matrices into a single matrix.
+			%
+			%   Examples:
+			%     a = magic(3); b = pascal(3); 
+			%     c = cat(4,a,b)
+			%   produces a 3-by-3-by-1-by-2 result and
+			%     s = {a b};
+			%     for i=1:length(s), 
+			%       siz{i} = size(s{i});
+			%     end
+			%     sizes = cat(1,siz{:})
+			%   produces a 2-by-2 array of size vectors.
+			%     
+			%   See also NUM2CELL.
+
+			%   Copyright 1984-2005 The MathWorks, Inc.
+			%   Built-in function.
 			if nargin == 1
 				objCell = {a};
 			else
-				objCell = {a, varargin{:}};
-			end
-			[fineGrid, fineGridName] = getFinestGrid(objCell{:});
-			for it = 1 : nargin
-				objCopyTemp = copy(objCell{it});
-				objCell{it} = objCopyTemp.changeGrid(fineGrid, fineGridName);
+				objCell = [{a}, varargin(:)'];
 			end
+			
 			assertSameGrid(objCell{:});
-			c = builtin('cat', objCell{:});
+			argin = [{dim}, objCell(:)'];
+			c = builtin('cat', argin{:});
 		end
 		
 		function solution = solveAlgebraic(obj, rhs, gridName, objLimit)
diff --git a/+unittests/+quantity/testDiscrete.m b/+unittests/+quantity/testDiscrete.m
index d4e846bad62453be2fc09da402b71e268a76ef85..38ccd91ef522dbd442ba503e0fc6d8d5b1bb336c 100644
--- a/+unittests/+quantity/testDiscrete.m
+++ b/+unittests/+quantity/testDiscrete.m
@@ -4,6 +4,34 @@ function [tests] = testDiscrete()
 tests = functiontests(localfunctions);
 end
 
+function testConcatenate(testCase)
+
+t = linspace(0, pi)';
+A = quantity.Discrete({sin(t); cos(t)}, 'grid', {t}, 'gridName', 't');
+B = quantity.Discrete({tan(t); exp(t)}, 'grid', {t}, 'gridName', 't');
+
+AB = [A, B];
+AB_ = [A', B'];
+ABA = [A, B, A];
+
+testCase.verifyTrue(all(size(AB) == [2, 2]));
+testCase.verifyTrue(all(size(AB_) == [1, 4]));
+testCase.verifyTrue(all(size(ABA) == [2, 3]));
+
+t = linspace(0, pi, 13)';
+C = quantity.Discrete({sin(t); cos(t)}, 'grid', {t}, 'gridName', 't');
+didNotWork = false;
+
+try
+	AC = [A, C];
+catch ex
+	didNotWork = true;
+end
+
+testCase.verifyTrue(didNotWork);
+
+end
+
 function testExp(testCase)
 % 1 spatial variable
 syms z zeta