diff --git a/+quantity/Discrete.m b/+quantity/Discrete.m
index 537fd1370fe8bca4f59d3c72f15e4d4d2f1c9d27..027ed5598ec90b1750e31139493d92d1aeb3e72a 100644
--- a/+quantity/Discrete.m
+++ b/+quantity/Discrete.m
@@ -235,6 +235,8 @@ classdef  (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
 				o(k) = quantity.Function(@(varargin) F(varargin{:}), ...
 					props{:});
 			end
+			
+			o = reshape(o, size(obj));
 		end
 		function o = quantity.Operator(obj)
 			A = cell(size(obj, 3), 1);
diff --git a/+quantity/Symbolic.m b/+quantity/Symbolic.m
index ecf0c034cafcaabcfd4c8895c6b36b592d15c35e..179b37c4de48a4c11de426bf2f51c32d4c8fa8b9 100644
--- a/+quantity/Symbolic.m
+++ b/+quantity/Symbolic.m
@@ -179,6 +179,22 @@ classdef Symbolic < quantity.Function
 		function f = function_handle(obj)
 			f = matlabFunction(sym(obj));
 		end
+		function F = quantity.Function(obj, varargin)
+			myParser = misc.Parser();
+			myParser.addParameter('grid', obj(1).grid);
+			myParser.addParameter('gridName', obj(1).gridName);
+			myParser.addParameter('name', obj(1).name);
+			myParser.parse(varargin{:});
+			assert(isequal(myParser.Results.gridName, obj(1).gridName))
+			for k = 1:numel(obj)
+				F(k) = quantity.Function(obj(k).function_handle(), ...
+					'grid', myParser.Results.grid, ...
+					'gridName', myParser.Results.gridName, ...
+					'name', myParser.Results.name);			
+			end
+			
+			F = reshape(F, size(obj));
+		end
 	end % methods
 	
 	%% mathematical operations
@@ -603,7 +619,7 @@ classdef Symbolic < quantity.Function
 					'name', ['c ', B(1).name]);
 				return
 			end
-			if isa(B, 'quantity.Function')
+			if ~(isa(B, 'quantity.Symbolic')) && isa(B, 'quantity.Function')
 				C = (B' * A')';
 				return;
 			end
diff --git a/+unittests/+quantity/testDiscrete.m b/+unittests/+quantity/testDiscrete.m
index 2e9ec10d75af926da0177dd6c7a78776c74556e7..a53959b712d94ba2cabd071fe311bc38f86617cc 100644
--- a/+unittests/+quantity/testDiscrete.m
+++ b/+unittests/+quantity/testDiscrete.m
@@ -4,6 +4,16 @@ function [tests] = testDiscrete()
 tests = functiontests(localfunctions);
 end
 
+function testCastDiscrete2Function(testCase)
+
+z = linspace(0, 2*pi)';
+d = quantity.Discrete({sin(z); cos(z)}, 'grid', z, 'gridName', 'z');
+f = quantity.Function(d);
+
+testCase.verifyTrue(all(size(f) == size(d)));
+
+end
+
 function testQuadraticNorm(tc)
 blub = quantity.Discrete(ones(11, 4), 'grid', linspace(0, 1, 11), ...
 	'gridName', 'z', 'name', 'b');
@@ -282,16 +292,19 @@ testCase.verifyEqual(myMatrixDiscrete(2,2).on(), myVectorDiscrete(2,1).on());
 testCase.verifyEqual(numel(myVectorDiscrete), size(myMatrixDiscrete, 1));
 end
 function testVec2Diag(testCase)
-% quantity.Symbolic
-syms z
-myMatrixSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1, 0; 0, 0.9-z/2]);
-myVectorSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1; 0.9-z/2]);
-testCase.verifyEqual(myMatrixSymbolic.on(), myVectorSymbolic.vec2diag.on());
-
 % quantity.Discrete
-myMatrixDiscrete = quantity.Discrete(myMatrixSymbolic);
-myVectorDiscrete = quantity.Discrete(myVectorSymbolic);
-testCase.verifyEqual(myMatrixDiscrete.on(), myVectorDiscrete.vec2diag.on());
+n = 7;
+z = linspace(0,1,n)';
+myMatrixDiscrete = quantity.Discrete(...
+	{sin(0.5*z*pi)+1, zeros(n,1); zeros(n,1), 0.9-z/2}, ...
+	'grid', z, ...
+	'gridName', 'z');
+myVectorDiscrete = quantity.Discrete(...
+	{sin(0.5*z*pi)+1; 0.9-z/2}, ...
+	'grid', z, ...
+	'gridName', 'z');
+
+testCase.verifyTrue( myVectorDiscrete.vec2diag.near(myMatrixDiscrete) );
 end
 
 function testInvert(testCase)
diff --git a/+unittests/+quantity/testFunction.m b/+unittests/+quantity/testFunction.m
index 326592108aaf093039d1dca776810b1780479774..b5e80cb58e0ce89c6598265fd8deb81bf39317f5 100644
--- a/+unittests/+quantity/testFunction.m
+++ b/+unittests/+quantity/testFunction.m
@@ -4,6 +4,18 @@ function [tests] = testFunction()
 tests = functiontests(localfunctions);
 end
 
+function testCastSymbolic2Function(testCase)
+
+z = linspace(0, 2*pi)';
+Z = sym("z");
+s = quantity.Symbolic([sin(Z); cos(Z)], 'grid', z, 'gridName', 'z');
+
+f = quantity.Function(s);
+
+testCase.verifyTrue(all(size(f) == size(s)));
+
+end
+
 function testTimesSymbolic(testCase)
 
 z = linspace(0, 2*pi)';
diff --git a/+unittests/+quantity/testSymbolic.m b/+unittests/+quantity/testSymbolic.m
index 9010cb76b020d72a68d16b27d8db0df4350f6920..5634c445b328a6105b88f0d24e20ed1299677aa2 100644
--- a/+unittests/+quantity/testSymbolic.m
+++ b/+unittests/+quantity/testSymbolic.m
@@ -4,6 +4,15 @@ function [tests ] = testSymbolic()
 tests = functiontests(localfunctions());
 end
 
+function testVec2Diag(testCase)
+% quantity.Symbolic
+syms z
+myMatrixSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1, 0; 0, 0.9-z/2]);
+myVectorSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1; 0.9-z/2]);
+
+testCase.verifyTrue( myVectorSymbolic.vec2diag.near(myMatrixSymbolic) );
+end
+
 function testSymbolicEvaluation(tc)
 syms z
 myGrid = linspace(0, 1, 7);