diff --git a/+misc/Parser.m b/+misc/Parser.m
index a5545a45e98306278e4c43f6e1f538e5c4bcc15d..4df1f6a5ee236e7703eec28cfcb93e6fb38fd337 100644
--- a/+misc/Parser.m
+++ b/+misc/Parser.m
@@ -6,13 +6,15 @@ classdef Parser < inputParser
methods
- function obj = Parser()
+ function obj = Parser(varargin)
obj = obj@inputParser();
% set default properties
obj.KeepUnmatched = true;
obj.CaseSensitive = true;
obj.PartialMatching = false;
+
+
end % Parser()
function i = isDefault(obj, name)
diff --git a/+misc/ensureIsCell.m b/+misc/ensureIsCell.m
new file mode 100644
index 0000000000000000000000000000000000000000..b93e156a35a6bc0ec047105af10ddd3ab10ea612
--- /dev/null
+++ b/+misc/ensureIsCell.m
@@ -0,0 +1,11 @@
+function value = ensureIsCell(value)
+%ENSUREISCELL ensures that the value is a cell.
+% c = ensureIsCell(value) checks if value is a cell. If it is not a cell,
+% it is converted as a cell.
+
+if ~iscell(value)
+ value = {value};
+end
+
+end
+
diff --git a/+quantity/Discrete.m b/+quantity/Discrete.m
index d30e3e6489dcbbabf66d4ef52101d6392657c292..123bbc1b80fa5ed2770b1d9c521de24e669c2b4e 100644
--- a/+quantity/Discrete.m
+++ b/+quantity/Discrete.m
@@ -2,9 +2,6 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
properties (SetAccess = protected)
% Discrete evaluation of the continuous quantity
valueDiscrete double;
-
- % Name of the domains that generate the grid.
- gridName {mustBe.unique};
end
properties (Hidden, Access = protected, Dependent)
@@ -12,6 +9,21 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
properties
+ % ID of the figure handle in which the handle is plotted
+ figureID double = 1;
+
+ % Name of this object
+ name char;
+
+ % domain
+ domain;
+ end
+
+ properties ( Dependent )
+
+ % Name of the domains that generate the grid.
+ gridName {mustBe.unique};
+
% Grid for the evaluation of the continuous quantity. For the
% example with the function f(x,t), the grid would be
% {[<spatial domain>], [<temporal domain>]}
@@ -19,12 +31,6 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% spatial domain and <temporal domain> the discrete description of
% the temporal domain.
grid; % in set.grid it is ensured that, grid is a (1,:)-cell-array
-
- % ID of the figure handle in which the handle is plotted
- figureID double = 1;
-
- % Name of this object
- name char;
end
methods
@@ -32,13 +38,21 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% --- Constructor ---
%--------------------
function obj = Discrete(valueOriginal, varargin)
- % The constructor requires valueOriginal to be
+ % DISCRETE a quantity, represented by discrete values.
+ % obj = Discrete(valueOriginal, varargin) initializes a
+ % quantity. The parameters to be set are:
+ % 'valueOrigin' must be
% 1) a cell-array of double arrays with
% size(valueOriginal) == size(obj) and
% size(valueOriginal{it}) == gridSize
+ % Example: valueOrigin = { f(Z, T), g(Z, T) } is a cell array
+ % wich contains the functions f(z,t) and g(z,t) evaluated on
+ % the discrete domain (Z x T). Then, the name-value-pair
+ % parameter 'domain' must be set with quantity.Domain
+ % objects, according to the domains Z and T.
% OR
% 2) adouble-array with
- % size(valueOriginal) == [gridSize, size(quantity)]
+ % size(valueOriginal) == [gridSize, size(quantity)]
% Furthermore, 'gridName' must be part of the name-value-pairs
% in varargin. Additional parameters can be specified using
% name-value-pair-syntax in varargin.
@@ -56,7 +70,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
obj = valueOriginal;
else
% empty object. this is needed for instance, to create
- % quantity.Discrete([]), which is useful for creating default
+ % quantity.Discrete([]), which is useful for creating default
% values.
obj = quantity.Discrete.empty(size(valueOriginal));
end
@@ -64,106 +78,122 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
%% input parser
myParser = misc.Parser();
- myParser.addParameter('gridName', [], @(g) ischar(g) || iscell(g));
- myParser.addParameter('grid', [], @(g) isnumeric(g) || iscell(g));
myParser.addParameter('name', string(), @isstr);
myParser.addParameter('figureID', 1, @isnumeric);
myParser.parse(varargin{:});
- assert(all(~contains(myParser.UsingDefaults, 'gridName')), ...
- 'gridName is a mandatory input for quantity');
- if iscell(myParser.Results.gridName)
- myGridName = myParser.Results.gridName;
+ %% domain parser
+ domainParser = misc.Parser();
+ domainParser.addParameter('domain', {}, @(g) isa(g, 'quantity.Domain'));
+ domainParser.addParameter('gridName', '', @(g) ischar(g) || iscell(g));
+ domainParser.addParameter('grid', [], @(g) isnumeric(g) || iscell(g));
+ domainParser.parse(varargin{:});
+
+ if domainParser.isDefault('domain') && ...
+ domainParser.isDefault('grid') && ...
+ domainParser.isDefault('gridName')
+ % case 1: nothing about the grid is defined
+ % -> use default grid
+ assert(iscell(valueOriginal), 'If no grid is specified, valueOriginal must be a cell-array')
+ valueOriginalGridSize = size(valueOriginal{1});
+ myDomain = quantity.Domain.defaultGrid(valueOriginalGridSize);
+
+ elseif domainParser.isDefault('domain') && ...
+ domainParser.isDefault('grid')
+ % case 2: the gridNames are specified
+ % -> use default grid with specified names
+ assert(iscell(valueOriginal), 'If no grid is specified, valueOriginal must be a cell-array')
+ valueOriginalGridSize = size(valueOriginal{1});
+ myDomain = quantity.Domain.defaultGrid(valueOriginalGridSize, ...
+ misc.ensureIsCell(domainParser.Results.gridName));
+
+ elseif domainParser.isDefault('domain')
+ % case 3: the gridNames and the gridValues are defined:
+ % -> initialize quantity.Domain objects with the
+ % specified values
+
+ myGridName = misc.ensureIsCell(domainParser.Results.gridName);
+ myGrid = misc.ensureIsCell(domainParser.Results.grid);
+
+ assert(isequal(size(myGrid), size(myGridName)), ...
+ 'number of grids and gridNames must be equal');
+
+ % initialize the domain objects
+ myDomain = quantity.Domain.empty();
+ for k = 1:numel(myGrid)
+ myDomain(k) = quantity.Domain('grid', myGrid{k}, ...
+ 'name', myGridName{k});
+ end
else
- myGridName = {myParser.Results.gridName};
+ % else case: the domains are specified as domain
+ % objects.
+ myDomain = domainParser.Results.domain;
end
- %% allow initialization of empty objects:
+ % #TODO check uniqueness of gridNames
+
+ %% get the sizes of obj and grid
+ gridLength = myDomain.gridLength;
+
+ % convert double valued valueOriginal to cell-valued value
+ % original
+ if ~iscell(valueOriginal)
+ valueOriginal = quantity.Discrete.value2cell(valueOriginal, gridLength);
+ end
+
+ % Check if the grid fits to the values. In addition, catch
+ % the case if all values are empty. This is required for
+ % the initialization of quantity.Function and
+ % quantity.Symbolic objects
+ assert( numGridElements(myDomain) == numel(valueOriginal{1}) || ...
+ misc.alln( cellfun(@isempty, valueOriginal ) ), ...
+ 'grids do not fit to valueOriginal');
+
+ % allow initialization of empty objects
valueOriginalSize = size(valueOriginal);
if any(valueOriginalSize == 0)
% If the size is specified in the arguements, it should
% be chosen instead of the default size from the
% valueOriginal.
myParser = misc.Parser();
- myParser.addParameter('size', valueOriginalSize((1+numel(myGridName)):end));
+ myParser.addParameter('size', valueOriginalSize((1+ndims(myDomain)):end));
myParser.parse(varargin{:});
obj = quantity.Discrete.empty(myParser.Results.size);
return;
end
- %% get the sizes of obj and grid
- if iscell(valueOriginal)
- if isempty(valueOriginal{1})
- % if valueOriginal is a cell-array with empty
- % cells, then grid must be specified as an input
- % parameter. This case is important for
- % constructing Symbolic or Function quantities
- % without discrete values.
- assert(all(~contains(myParser.UsingDefaults, 'grid')), ...
- ['grid is a mandatory input for quantity, ', ...
- 'if no discrete values are specified']);
- if ~iscell(myParser.Results.grid)
- gridSize = numel(myParser.Results.grid);
- else
- gridSize = cellfun(@(v) numel(v), myParser.Results.grid);
- end
+ % set valueDiscrete
+ for k = 1:numel(valueOriginal)
+ if numel(myDomain) == 1
+ % TODO: Which case is this? Why does it need extra
+ % treatment?
+ obj(k).valueDiscrete = valueOriginal{k}(:); %#ok<AGROW>
else
- gridSize = size(valueOriginal{1});
- end
- objSize = size(valueOriginal);
- elseif isnumeric(valueOriginal)
- gridSize = valueOriginalSize(1 : numel(myGridName));
- objSize = [valueOriginalSize(numel(myGridName)+1 : end), 1, 1];
- end
-
- %% get grid and check size
- if any(contains(myParser.UsingDefaults, 'grid'))
- myGrid = quantity.Discrete.defaultGrid(gridSize);
- else
- myGrid = myParser.Results.grid;
- end
- if ~iscell(myGrid)
- myGrid = {myGrid};
- end
- if isempty(myGridName) || isempty(myGrid)
- if ~(isempty(myGridName) && isempty(myGrid))
- error(['If one of grid and gridName is empty, ', ...
- 'then both must be empty.']);
- end
- else
- assert(isequal(size(myGrid), size(myGridName)), ...
- 'number of grids and gridNames must be equal');
- myGridSize = cellfun(@(v) numel(v), myGrid);
- assert(isequal(gridSize(gridSize>1), myGridSize(myGridSize>1)), ...
- 'grids do not fit to valueOriginal');
- end
-
- %% set valueDiscrete
- if ~iscell(valueOriginal)
- valueOriginal = quantity.Discrete.value2cell(valueOriginal, gridSize, objSize);
- end
- for k = 1:prod(objSize)
- if numel(myGrid) == 1
- obj(k).valueDiscrete = valueOriginal{k}(:);
- else
- obj(k).valueDiscrete = valueOriginal{k};
+ obj(k).valueDiscrete = valueOriginal{k}; %#ok<AGROW>
end
end
%% set further properties
- [obj.grid] = deal(myGrid);
- [obj.gridName] = deal(myGridName);
+ [obj.domain] = deal(myDomain);
[obj.name] = deal(myParser.Results.name);
[obj.figureID] = deal(myParser.Results.figureID);
%% reshape object from vector to matrix
- obj = reshape(obj, objSize);
+ obj = reshape(obj, size(valueOriginal));
end
end% Discrete() constructor
-
+
%---------------------------
% --- getter and setters ---
- %---------------------------
+ %---------------------------
+ function gridName = get.gridName(obj)
+ gridName = {obj.domain.name};
+ end
+
+ function grid = get.grid(obj)
+ grid = {obj.domain.grid};
+ end
+
function itIs = isConstant(obj)
% the quantity is interpreted as constant if it has no grid or
% it has a grid that is only one point.
@@ -178,7 +208,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
obj.gridName = name;
end
-
+
function set.grid(obj, grid)
if ~iscell(grid)
grid = {grid};
@@ -225,7 +255,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
function o = quantity.Function(obj)
props = nameValuePair( obj(1) );
-
+
for k = 1:numel(obj)
F = griddedInterpolant(obj(k).grid{:}', obj(k).on());
o(k) = quantity.Function(@(varargin) F(varargin{:}), ...
@@ -249,7 +279,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
error('Not yet implemented')
end
end
-
+
function obj = setName(obj, newName)
% Function to set all names of all elements of the quantity obj to newName.
[obj.name] = deal(newName);
@@ -307,13 +337,13 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
value = zeros(size(obj));
indexGrid = arrayfun(@(s)linspace(1,s,s), size(obj), 'UniformOutput', false);
interpolant = numeric.interpolant(...
- [indexGrid{:}], value);
+ [indexGrid{:}], value);
else
- myGrid = obj(1).grid;
+ myGrid = obj(1).grid;
value = obj.obj2value();
indexGrid = arrayfun(@(s)linspace(1,s,s), size(obj), 'UniformOutput', false);
interpolant = numeric.interpolant(...
- [myGrid, indexGrid{:}], value);
+ [myGrid, indexGrid{:}], value);
end
end
@@ -363,15 +393,15 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
[referenceGrid, referenceGridName] = varargin{1}.getFinestGrid(varargin{2:end});
else
referenceGrid = {};
- referenceGridName = '';
+ referenceGridName = '';
end
return;
else
referenceGridName = a(1).gridName;
referenceGrid = a(1).grid;
- referenceGridSize = a(1).gridSize(referenceGridName);
+ referenceGridSize = a(1).gridSize(referenceGridName);
end
-
+
for it = 1 : numel(varargin)
if isempty(varargin{it})
continue;
@@ -411,9 +441,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
tempGrid1 = obj1(1).grid{index1};
tempGrid2 = obj2(1).grid{index2};
- if ~obj1.isConstant && ~obj2.isConstant
- assert(tempGrid1(1) == tempGrid2(1), 'Grids must have same domain for gridJoin')
- assert(tempGrid1(end) == tempGrid2(end), 'Grids must have same domain for gridJoin')
+ if ~obj1.isConstant && ~obj2.isConstant
+ assert(tempGrid1(1) == tempGrid2(1), 'Grids must have same domain for gridJoin')
+ assert(tempGrid1(end) == tempGrid2(end), 'Grids must have same domain for gridJoin')
end
if numel(tempGrid1) > numel(tempGrid2)
gridJoined{it} = tempGrid1;
@@ -425,10 +455,10 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end % gridJoin()
function obj = sort(obj, varargin)
%SORT sorts the grid of the object in a desired order
- % obj = sortGrid(obj) sorts the grid in alphabetical order.
+ % obj = sortGrid(obj) sorts the grid in alphabetical order.
% obj = sort(obj, 'descend') sorts the grid in descending
% alphabetical order.
-
+
if nargin == 2 && strcmp(varargin{1}, 'descend')
descend = 1;
else
@@ -442,23 +472,23 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% this is the default case for ascending alphabetical
% order
[sortedNames, I] = sort(gridNames);
-
+
% if descending: flip the order of the entries
- if descend
+ if descend
sortedNames = flip(sortedNames);
I = flip(I);
end
-
+
% sort the grid entries
[obj.grid] = deal(obj(1).grid(I));
-
+
% assign the new grid names
[obj.gridName] = deal(sortedNames);
-
+
% permute the value discrete
for k = 1:numel(obj)
obj(k).valueDiscrete = permute(obj(k).valueDiscrete, I);
- end
+ end
end
end% sort()
function c = horzcat(a, varargin)
@@ -488,7 +518,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% X3 ...]' when any of X1, X2, X3, etc. is an object.
%
% See also HORZCAT, CAT.
- c = cat(2, a, varargin{:});
+ c = cat(2, a, varargin{:});
end
function c = vertcat(a, varargin)
%VERTCAT Vertical concatenation.
@@ -535,18 +565,18 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% single matrix.
%
% Examples:
- % a = magic(3); b = pascal(3);
+ % 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),
+ % for i=1:length(s),
% siz{i} = size(s{i});
% end
% sizes = cat(1,siz{:})
% produces a 2-by-2 array of size vectors.
-
+
if nargin == 1
- objCell = {a};
+ objCell = {a};
else
objCell = [{a}, varargin(:)'];
@@ -565,7 +595,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% concatenated, so a new empty object is initialized.
s = cellfun(@(o) size(o), objCell, 'UniformOutput', false);
if dim == 1
- S = sum(cat(3, s{:}), 3);
+ S = sum(cat(3, s{:}), 3);
elseif dim == 2
S = s{1};
else
@@ -589,7 +619,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if isempty(value)
M = zeros([prod(obj(1).gridSize), size(value(l))]);
end
- M = reshape(M, [obj(1).gridSize, size(value)]);
+ M = reshape(M, [obj(1).gridSize, size(value)]);
o = quantity.Discrete( M, ...
'size', size(value), ...
'gridName', obj(1).gridName, ...
@@ -683,13 +713,13 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
solution = objInverseTemp.on(rhs);
-
-% solution = zeros(numel(obj), 1);
-% for it = 1 : numel(obj)
-% objInverseTemp = obj(it).invert(gridName);
-% solution(it) = objInverseTemp.on(rhs(it));
-% end
-% solution = reshape(solution, size(obj));
+
+ % solution = zeros(numel(obj), 1);
+ % for it = 1 : numel(obj)
+ % objInverseTemp = obj(it).invert(gridName);
+ % solution(it) = objInverseTemp.on(rhs(it));
+ % end
+ % solution = reshape(solution, size(obj));
end
function inverse = invert(obj, gridName)
@@ -709,7 +739,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
assert(isequal(size(obj), [1, 1]));
inverse = quantity.Discrete(repmat(obj(1).grid{obj.gridIndex(gridName)}(:), [1, size(obj)]), ...
'size', size(obj), 'grid', obj.on(), 'gridName', {[obj(1).name]}, ...
- 'name', gridName);
+ 'name', gridName);
end
@@ -732,7 +762,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
myParser.parse(varargin{:});
variableGrid = myParser.Results.variableGrid;
- myGridSize = [numel(variableGrid), ...
+ myGridSize = [numel(variableGrid), ...
numel(myParser.Results.initialValueGrid)];
% the time (s) vector has to start at 0, to ensure the IC. If
@@ -753,14 +783,14 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if numel(positiveVariableGrid) > 1
[resultGridPositive, odeSolutionPositive] = ...
ode45(@(y, z) obj(it).on(z), ...
- positiveVariableGrid, ...
- myParser.Results.initialValueGrid(icIdx), options);
+ positiveVariableGrid, ...
+ myParser.Results.initialValueGrid(icIdx), options);
end
if numel(negativeVariableGrid) >1
[resultGridNegative, odeSolutionNegative] = ...
ode45(@(y, z) obj(it).on(z), ...
- negativeVariableGrid, ...
- myParser.Results.initialValueGrid(icIdx), options);
+ negativeVariableGrid, ...
+ myParser.Results.initialValueGrid(icIdx), options);
end
if any(variableGrid == 0)
resultGrid = [flip(resultGridNegative(2:end)); 0 ; resultGridPositive(2:end)];
@@ -782,7 +812,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
{variableGrid, myParser.Results.initialValueGrid}, ...
'size', size(obj), 'name', ['solve(', obj(1).name, ')']);
end
-
+
function solution = subs(obj, gridName2Replace, values)
if nargin == 1 || isempty(gridName2Replace)
% if gridName2Replace is empty, then nothing must be done.
@@ -803,8 +833,8 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
assert(numel(values) == numel(gridName2Replace), ...
'gridName2Replace and values must be of same size');
- % here substitution starts:
- % The first (gridName2Replace{1}, values{1})-pair is
+ % here substitution starts:
+ % The first (gridName2Replace{1}, values{1})-pair is
% replaced. If there are more cell-elements in those inputs
% then subs() is called again for the remaining pairs
% (gridName2Replace{2:end}, values{2:end}).
@@ -846,10 +876,10 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
newValue = misc.diagNd(obj.on(newGridForOn), gridIndices);
newGrid = {newGridForOn{gridIndices(1)}, ...
newGridForOn{1:1:numel(newGridForOn) ~= gridIndices(1) ...
- & 1:1:numel(newGridForOn) ~= gridIndices(2)}};
+ & 1:1:numel(newGridForOn) ~= gridIndices(2)}};
newGridName = {values{1}, ...
obj(1).gridName{1:1:numel(obj(1).gridName) ~= gridIndices(1) ...
- & 1:1:numel(obj(1).gridName) ~= gridIndices(2)}};
+ & 1:1:numel(obj(1).gridName) ~= gridIndices(2)}};
else
% this is the default case. just grid name is
@@ -904,7 +934,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
function value = at(obj, point)
- % at() evaluates the object at one point and returns it as array
+ % at() evaluates the object at one point and returns it as array
% with the same size as size(obj).
value = reshape(obj.on(point), size(obj));
end
@@ -913,11 +943,11 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% ATINDEX returns the valueDiscrete at the requested index.
% value = atIndex(obj, varargin) returns the
% quantity.Discrete.valueDiscrete at the index defined by
- % varargin.
+ % varargin.
% value = atIndex(obj, 1) returns the first element of
% "valueDiscrete"
% value = atIndex(obj, ':') returns all elements of
- % obj.valueDiscrete in vectorized form.
+ % obj.valueDiscrete in vectorized form.
% value = atIndex(obj, 1, end) returns the obj.valueDiscrete
% at the index (1, end).
% If a range of index is requested, the result is returned
@@ -941,15 +971,15 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
'UniformOutput', false);
valueSize = size(value{1});
-
+
if all(cellfun(@numel, varargin) == 1) && all(cellfun(@isnumeric, varargin))
outputSize = [];
else
outputSize = valueSize(1:obj(1).nargin);
- end
-
+ end
+
value = reshape([value{:}], [outputSize, size(obj)]);
- end
+ end
end
function n = nargin(obj)
@@ -1079,7 +1109,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
xlabel(labelHelper(1), 'Interpreter','latex');
if p.Results.showTitle
- title(titleHelper(), 'Interpreter','latex');
+ title(titleHelper(), 'Interpreter','latex');
end
a = gca();
a.TickLabelInterpreter = 'latex';
@@ -1092,10 +1122,10 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if p.Results.export
matlab2tikz(p.Results.exportOptions{:});
end
-
+
function myLabel = labelHelper(gridNumber)
if ~isempty(obj(rowIdx, columnIdx, figureIdx).gridName)
- myLabel = ['$$', greek2tex(obj(rowIdx, columnIdx, figureIdx).gridName{gridNumber}), '$$'];
+ myLabel = ['$$', greek2tex(obj(rowIdx, columnIdx, figureIdx).gridName{gridNumber}), '$$'];
else
myLabel = '';
end
@@ -1149,9 +1179,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
si = num2cell( size(obj) );
s(si{:}) = struct();
for l = 1:numel(obj)
-
+
doNotCopyProperties = obj(l).doNotCopy;
-
+
for k = 1:length(properties)
if ~any(strcmp(doNotCopyProperties, properties{k}))
s(l).(properties{k}) = obj(1).(properties{k});
@@ -1159,7 +1189,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
end
else
- % Without this else-case
+ % Without this else-case
% this method is in conflict with the default matlab built-in
% calls struct(). When calling struct('myFieldName', myQuantity) this
% quantity-method is called (and errors) and not the
@@ -1247,7 +1277,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% state-transition matrix for a ordinary differential equation
% d / dt x(t) = A(t) x(t)
% with A as this object. The system matrix can be constant or
- % variable. It is the solution of
+ % variable. It is the solution of
% d / dt F(t,t0) = A(t) F(t,t0)
% F(t,t) = I.
% Optional parameters for the computation can be defined as
@@ -1267,7 +1297,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
ip = misc.Parser();
ip.addParameter('grid', []);
ip.addParameter('gridName1', '');
- ip.addParameter('gridName2', '');
+ ip.addParameter('gridName2', '');
ip.parse(varargin{:});
myGrid = ip.Results.grid;
gridName1 = ip.Results.gridName1;
@@ -1282,14 +1312,14 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
assert(numel(obj(1).gridName) == 1, 'No gridName is defined. Use name-value-pairs property "gridName1" to define a grid name');
gridName1 = obj(1).gridName{1};
end
-
+
if ip.isDefault('gridName2')
gridName2 = [gridName1 '0'];
end
assert(numel(myGrid) > 1, 'If the state transition matrix is computed for constant values, a spatial domain has to be defined!')
-
+
if obj.isConstant
% for a constant system matrix, the matrix exponential
% function can be used.
@@ -1297,7 +1327,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
zeta = sym(gridName2, 'real');
f0 = expm(obj.atIndex(1)*(z - zeta));
F = quantity.Symbolic(f0, 'grid', {myGrid, myGrid});
-
+
elseif isa(obj, 'quantity.Symbolic')
f0 = misc.fundamentalMatrix.odeSolver_par(obj.function_handle, myGrid);
F = quantity.Discrete(f0, ...
@@ -1313,7 +1343,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
'gridName', {gridName1, gridName2}, ...
'grid', {myGrid, myGrid});
end
- end
+ end
function s = sum(obj, dim)
s = quantity.Discrete(sum(obj.on(), obj.nargin + dim), ...
@@ -1367,7 +1397,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function P = mtimes(a, b)
% TODO rewrite the selection of the special cases! the
% if-then-cosntruct is pretty ugly!
-
+
% numeric computation of the matrix product
if isempty(b) || isempty(a)
% TODO: actually this only holds for multiplication of
@@ -1388,7 +1418,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
return
else
b = quantity.Discrete(b, 'size', size(b), 'grid', {}, ...
- 'gridName', {}, 'name', '');
+ 'gridName', {}, 'name', '');
end
end
@@ -1412,7 +1442,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if isscalar(b)
% do the scalar multiplication in a loop
for k = 1:numel(a)
- P(k) = innerMTimes(a(k), b);
+ P(k) = innerMTimes(a(k), b);
end
P = reshape(P, size(a));
return
@@ -1421,13 +1451,13 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if isscalar(a)
% do the scalar multiplication in a loop
for k = 1:numel(b)
- P(k) = innerMTimes(a, b(k));
+ P(k) = innerMTimes(a, b(k));
end
P = reshape(P, size(b));
return
end
-
- P = innerMTimes(a, b);
+
+ P = innerMTimes(a, b);
end
function P = innerMTimes(a, b)
@@ -1438,7 +1468,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% misc.multArray is very efficient, but requires that the
% multiple dimensions of the input are correctly arranged.
[idx, permuteGrid] = computePermutationVectors(a, b);
-
+
parameters = struct(a(1));
parameters.name = [a(1).name, ' ', b(1).name];
parameters.grid = [...
@@ -1465,7 +1495,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if b.gridIndex(parameters.gridName{it})
bGrid{b.gridIndex(parameters.gridName{it})} = ...
parameters.grid{it};
- end
+ end
end
parameters = misc.struct2namevaluepair(parameters);
@@ -1494,11 +1524,11 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
objDiscreteReshaped = permute(reshape(objDiscreteOriginal, ...
[prod(obj(1).gridSize), size(obj)]), [2, 3, 1]);
invDiscrete = zeros([prod(obj(1).gridSize), size(obj)]);
-
+
parfor it = 1 : size(invDiscrete, 1)
invDiscrete(it, :, :) = inv(objDiscreteReshaped(:, :, it));
end
-
+
y = quantity.Discrete(reshape(invDiscrete, size(objDiscreteOriginal)),...
'size', size(obj), ...
'name', ['{(', obj(1).name, ')}^{-1}'], ...
@@ -1515,9 +1545,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
objT = obj.';
objCtMat = conj(objT.on());
objCt = quantity.Discrete(objCtMat, 'grid', obj(1).grid, ...
- 'gridName', obj(1).gridName, 'name', ['{', obj(1).name, '}^{H}']);
+ 'gridName', obj(1).gridName, 'name', ['{', obj(1).name, '}^{H}']);
end % ctranspose(obj)
-
+
function y = exp(obj)
% exp() is the exponential function using obj as the exponent.
y = quantity.Discrete(exp(obj.on()), ...
@@ -1534,7 +1564,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% The integral domain is specified by integralGrid.
if obj.nargin > 1
- assert(ischar(integralGridName), 'integralGrid must specify a gridName as a char');
+ assert(ischar(integralGridName), 'integralGrid must specify a gridName as a char');
else
integralGridName = obj(1).gridName;
end
@@ -1555,7 +1585,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% xNorm = sqrt(x.' * x).
% Optionally, a weight can be defined, such that instead
% xNorm = sqrt(x.' * weight * x).
-
+
myParser = misc.Parser();
myParser.addParameter('weight', eye(size(obj, 1)));
myParser.parse(varargin{:});
@@ -1592,7 +1622,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% number of rows."
x = inv(A) * B;
end
-
+
function x = mrdivide(B, A)
% mRdivide see doc mrdivide of matlab: "x = B/A solves the
% system of linear equations x*A = B for x. The matrices A and
@@ -1752,11 +1782,11 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
return
end
- if nargin < 3 % if no diffGridName is specified, then the derivative
+ if nargin < 3 % if no diffGridName is specified, then the derivative
% w.r.t. to all gridNames is applied
diffGridName = obj(1).gridName;
end
-
+
% diff for each element of diffGridName (this is rather
% inefficient, but an easy implementation of the specification)
@@ -1770,35 +1800,35 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
end
end
else
- gridSelector = strcmp(obj(1).gridName, diffGridName);
- gridSelectionIndex = find(gridSelector);
-
- spacing = gradient(obj(1).grid{gridSelectionIndex}, 1);
- assert(numeric.near(spacing, spacing(1)), ...
- 'diff is currently only implemented for equidistant grid');
- permutationVector = 1 : (numel(obj(1).grid)+ndims(obj));
-
- objDiscrete = permute(obj.on(), ...
- [permutationVector(gridSelectionIndex), ...
- permutationVector(permutationVector ~= gridSelectionIndex)]);
-
- if size(objDiscrete, 2) == 1
- derivativeDiscrete = gradient(objDiscrete, spacing(1));
- else
- [~, derivativeDiscrete] = gradient(objDiscrete, spacing(1));
- end
-
- rePermutationVector = [2:(gridSelectionIndex), ...
- 1, (gridSelectionIndex+1):ndims(derivativeDiscrete)];
- result = quantity.Discrete(...
- permute(derivativeDiscrete, rePermutationVector), ...
- 'size', size(obj), 'grid', obj(1).grid, ...
- 'gridName', obj(1).gridName, ...
+ gridSelector = strcmp(obj(1).gridName, diffGridName);
+ gridSelectionIndex = find(gridSelector);
+
+ spacing = gradient(obj(1).grid{gridSelectionIndex}, 1);
+ assert(numeric.near(spacing, spacing(1)), ...
+ 'diff is currently only implemented for equidistant grid');
+ permutationVector = 1 : (numel(obj(1).grid)+ndims(obj));
+
+ objDiscrete = permute(obj.on(), ...
+ [permutationVector(gridSelectionIndex), ...
+ permutationVector(permutationVector ~= gridSelectionIndex)]);
+
+ if size(objDiscrete, 2) == 1
+ derivativeDiscrete = gradient(objDiscrete, spacing(1));
+ else
+ [~, derivativeDiscrete] = gradient(objDiscrete, spacing(1));
+ end
+
+ rePermutationVector = [2:(gridSelectionIndex), ...
+ 1, (gridSelectionIndex+1):ndims(derivativeDiscrete)];
+ result = quantity.Discrete(...
+ permute(derivativeDiscrete, rePermutationVector), ...
+ 'size', size(obj), 'grid', obj(1).grid, ...
+ 'gridName', obj(1).gridName, ...
'name', ['(d_{', diffGridName, '}', obj(1).name, ')']);
if k > 1
-% % if a higher order derivative is requested, call the function
-% % recursivly until the first-order derivative is reached
+ % % if a higher order derivative is requested, call the function
+ % % recursivly until the first-order derivative is reached
result = result.diff(k-1, diffGridName);
end
end
@@ -1836,15 +1866,15 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if iscell(intGridName) && numel(intGridName) > 1
obj = obj.int(intGridName{1}, a(1), b(1));
I = obj.int(intGridName{2:end}, a(2:end), b(2:end));
- return
+ return
end
[idxGrid, isGrid] = obj.gridIndex(intGridName);
if nargin == 4
I = cumInt( obj, varargin{:});
return;
-%
-% assert(all((varargin{2} == a(idxGrid)) & (varargin{3} == b(idxGrid))), ...
-% 'only integration from beginning to end of domain is implemented');
+ %
+ % assert(all((varargin{2} == a(idxGrid)) & (varargin{3} == b(idxGrid))), ...
+ % 'only integration from beginning to end of domain is implemented');
end
%% do the integration over one dimension
@@ -1855,12 +1885,14 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
[domainIdx{1:obj(1).nargin}] = deal(':');
currentGrid = obj(1).grid{ idxGrid };
domainIdx{idxGrid} = find(currentGrid >= a(idxGrid) & currentGrid <= b(idxGrid));
-
+
myGrid = currentGrid(domainIdx{idxGrid});
+
+
for kObj = 1:numel(obj)
J = numeric.trapz_fast_nDim(myGrid, obj(kObj).atIndex(domainIdx{:}), idxGrid);
-
+
% If the quantity is integrated, so that the first
% dimension collapses, this dimension has to be shifted
% away, so that the valueDiscrete can be set correctly.
@@ -1882,7 +1914,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
I(kObj).grid = newGrid;
I(kObj).name = ['int(' obj(kObj).name ')'];
end
-
+
end
function result = cumInt(obj, domain, lowerBound, upperBound)
@@ -1915,7 +1947,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
result = quantity.Discrete(F, 'grid', myGrid, ...
'gridName', obj(1).gridName);
- % int_lowerBound^upperBound f(.) =
+ % int_lowerBound^upperBound f(.) =
% F(upperBound) - F(lowerBound)
result = result.subs({domain}, upperBound) ...
- result.subs({domain}, lowerBound);
@@ -1923,7 +1955,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
result.setName(deal(['int(', obj(1).name, ')']));
end
end
-
+
function C = plus(A, B)
%% PLUS is the sum of two quantities.
assert(isequal(size(A), size(B)), 'plus() not supports mismatching sizes')
@@ -1954,7 +1986,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
[aDiscrete] = A.expandValueDiscrete(gridJoined, gridNameJoined, gridJoinedLength);
[bDiscrete] = B.expandValueDiscrete(gridJoined, gridNameJoined, gridJoinedLength);
-
+
% create result object
C = quantity.Discrete(aDiscrete + bDiscrete, ...
'grid', gridJoined, 'gridName', gridNameJoined, ...
@@ -1965,7 +1997,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% EXPANDVALUEDISCRETE
% [valDiscrete, gridNameSorted] = ...
% expandValueDiscrete(obj, gridIndex, valDiscrete)
- % Expanses the value of obj, so that
+ % Expanses the value of obj, so that
% todo
% get the index of obj.grid in the joined grid
@@ -1975,15 +2007,15 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
oldDim = ndims(valDiscrete);
valDiscrete = permute(valDiscrete, [(1:sum(~gridLogical)) + oldDim, 1:oldDim] );
valDiscrete = repmat(valDiscrete, [gridJoinedLength(~gridLogical), ones(1, ndims(valDiscrete))]);
-%
- valDiscrete = reshape(valDiscrete, ...
- [gridJoinedLength(~gridLogical), gridJoinedLength(gridLogical), size(obj)]);
-
-% % permute valDiscrete such that grids are in the order specified
-% % by gridNameJoined.
- gridIndex = 1:numel(gridLogical);
+ %
+ valDiscrete = reshape(valDiscrete, ...
+ [gridJoinedLength(~gridLogical), gridJoinedLength(gridLogical), size(obj)]);
+
+ % % permute valDiscrete such that grids are in the order specified
+ % % by gridNameJoined.
+ gridIndex = 1:numel(gridLogical);
gridOrder = [gridIndex(~gridLogical), gridIndex(gridLogical)];
- valDiscrete = permute(valDiscrete, [gridOrder, numel(gridLogical)+(1:ndims(obj))]);
+ valDiscrete = permute(valDiscrete, [gridOrder, numel(gridLogical)+(1:ndims(obj))]);
end
@@ -1997,7 +2029,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
elseif isnumeric(B)
[C.name] = deal([A(1).name, '-c']);
else
- [C.name] = deal([A(1).name, '-', B(1).name]);
+ [C.name] = deal([A(1).name, '-', B(1).name]);
end
end
@@ -2034,7 +2066,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function maxValue = max(obj)
% max returns the maximal value of all elements of A over all
% variables as a double array.
- maxValue = reshape(max(obj.on(), [], 1:obj(1).nargin), [size(obj)]);
+ maxValue = reshape(max(obj.on(), [], 1:obj(1).nargin), [size(obj)]);
end
function [i, m, s] = near(obj, B, varargin)
@@ -2078,7 +2110,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function minValue = min(obj)
% max returns the minimum value of all elements of A over all
% variables as a double array.
- minValue = reshape(min(obj.on(), [], 1:obj(1).nargin), [size(obj)]);
+ minValue = reshape(min(obj.on(), [], 1:obj(1).nargin), [size(obj)]);
end % min()
function absQuantity = abs(obj)
@@ -2086,9 +2118,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
if isempty(obj)
absQuantity = quantity.Discrete.empty(size(obj));
else
- absQuantity = quantity.Discrete(abs(obj.on()), ...
- 'grid', obj(1).grid, 'gridName', obj(1).gridName, ...
- 'size', size(obj), 'name', ['|', obj(1).name, '|']);
+ absQuantity = quantity.Discrete(abs(obj.on()), ...
+ 'grid', obj(1).grid, 'gridName', obj(1).gridName, ...
+ 'size', size(obj), 'name', ['|', obj(1).name, '|']);
end
end % abs()
@@ -2111,7 +2143,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function meanValue = mean(obj)
% mean returns the mean value of all elements of A over all
% variables as a double array.
- meanValue = reshape(mean(obj.on(), 1:obj(1).nargin), size(obj));
+ meanValue = reshape(mean(obj.on(), 1:obj(1).nargin), size(obj));
end
function medianValue = median(obj)
@@ -2120,9 +2152,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% ATTENTION:
% median over multiple dimensions is only possible from MATLAB 2019a and above!
% TODO Update Matlab!
- medianValue = reshape(median(obj.on(), 1:obj(1).nargin), size(obj));
+ medianValue = reshape(median(obj.on(), 1:obj(1).nargin), size(obj));
end
-
+
function q = obj2value(obj)
% for a case with a 3-dimensional grid, the common cell2mat
% with reshape did not work. So it has to be done manually:
@@ -2187,27 +2219,20 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function q = value2cell(value, gridSize, valueSize)
% VALUE2CELL
- gs = num2cell(gridSize(:));
- vs = arrayfun(@(n) ones(n, 1), valueSize, 'UniformOutput', false);
- s = [gs(:); vs(:)]; %{gs{:}, vs{:}};%
- q = reshape(mat2cell(value, s{:}), valueSize);
- end
-
- function g = defaultGrid(gridSize)
- g = cell(1, length(gridSize));
- for k = 1:length(gridSize)
-
- o = ones(1, length(gridSize) + 1); % + 1 is required to deal with one dimensional grids
- o(k) = gridSize(k);
- O = ones(o);
- O(:) = linspace(0, 1, gridSize(k));
-
- g{k} = O;
+ fullSize = size(value);
+
+ myGridSize = num2cell( fullSize(1:length(gridSize)) );
+
+ if nargin == 2
+ valueSize = [fullSize(length(myGridSize)+1:length(fullSize)), 1];
end
+
+ myValueSize = arrayfun(@(n) ones(n, 1), valueSize, 'UniformOutput', false);
+ s = [myGridSize(:); myValueSize(:)]; %{gs{:}, vs{:}};%
+ q = reshape(mat2cell(value, s{:}), valueSize);
end
-
function [p, q, parameters, gridSize] = inputNormalizer(a, b)
parameters = [];
@@ -2242,7 +2267,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function [idx, permuteGrid] = computePermutationVectors(a, b)
% Computes the required permutation vectors to use
% misc.multArray for multiplication of matrices
-
+
uA = unique(a(1).gridName, 'stable');
assert(numel(uA) == numel(a(1).gridName), 'Gridnames have to be unique!');
uB = unique(b(1).gridName, 'stable');
@@ -2313,7 +2338,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% quantity.Discrete class.
permuteGrid = [idxGrid(lGridA), idxGrid(lGridB), idxGrid(lGridVA), idxGrid(lGridVB)];
end
-
+
function f = setValueContinuous(obj, f)
end
@@ -2339,7 +2364,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
s = getPropertyGroups@matlab.mixin.CustomDisplay(obj);
return;
else
- s = getPropertyGroups@matlab.mixin.CustomDisplay(obj(1));
+ s = getPropertyGroups@matlab.mixin.CustomDisplay(obj(1));
end
if numel(obj) ~= 1
diff --git a/+quantity/Domain.m b/+quantity/Domain.m
new file mode 100644
index 0000000000000000000000000000000000000000..1484f5c5e9645fe3aaf11cfbb27766423c64598f
--- /dev/null
+++ b/+quantity/Domain.m
@@ -0,0 +1,103 @@
+classdef Domain
+ %DOMAIN class to describes a range of values on which a function can be defined.
+
+ % todo:
+ % * EquidistantDomain
+ % * multi dimensional
+
+ properties
+ % The discrete points of the grid for the evaluation of a
+ % continuous quantity. For an example, the function f(x) should be
+ % considered on the domain x \in X = [0, 1]. Then, a grid can be
+ % generated by X_grid = linspace(0, 1).
+ grid double {mustBeReal};
+
+ % a speaking name for this domain; Should be unique, so that the
+ % domain can be identified by the name.
+ name char;
+ end
+
+ properties (Dependent)
+ n; % number of discretization points
+ lower; % lower bound of the domain
+ upper; % upper bound of the domain
+ end
+
+ methods
+ function obj = Domain(varargin)
+ %DOMAIN initialize the domain
+ %
+ if nargin >= 1
+ parser = misc.Parser();
+ parser.addParameter('grid', [], @isvector);
+ parser.addParameter('name', '', @ischar);
+ parser.parse(varargin{:});
+
+ % todo: assertions
+ % * ascending ?
+
+ obj.grid = parser.Results.grid(:);
+ obj.name = parser.Results.name;
+ else
+ obj = quantity.Domain.empty();
+ end
+ end
+
+ function n = get.n(obj)
+ n = length(obj.grid);
+ end
+
+ function lower = get.lower(obj)
+ lower = min( obj.grid );
+ end
+
+ function upper = get.upper(obj)
+ upper = max( obj.grid );
+ end
+
+ function nd = ndims(obj)
+ %NDIMS number of dimensions of the domain specified by the
+ %object-array.
+ nd = size(obj(:), 1);
+ end
+
+ function n = numGridElements(obj)
+ % NUMGRIDLEMENTS returns the number of the elements of the grid
+ n = prod([obj.n]);
+ end
+
+ function s = gridLength(obj)
+ %GRIDLENGTH number of discretization points for each grid in the
+ %object-array.
+ s = [obj.n];
+ end
+
+ end
+
+ methods (Static)
+ function g = defaultGrid(gridSize, name)
+
+ if nargin == 1
+ % If no names are specified, chose x_1, x_2, ... as default
+ % names.
+ name = cell(1, length(gridSize));
+ for k = 1:length(gridSize)
+ name{k} = ['x_' num2str(k)];
+ end
+ end
+
+ % generate a default gird with given sizes
+ g = quantity.Domain.empty();
+ for k = 1:length(gridSize)
+
+ o = ones(1, length(gridSize) + 1); % + 1 is required to deal with one dimensional grids
+ o(k) = gridSize(k);
+ O = ones(o);
+ O(:) = linspace(0, 1, gridSize(k));
+
+ g(k) = quantity.Domain('grid', O, 'name', name{k});
+ end
+ end
+ end
+end
+
diff --git a/+quantity/EquidistantDomain.m b/+quantity/EquidistantDomain.m
new file mode 100644
index 0000000000000000000000000000000000000000..ff763973dfff6ae55c52e884a40144faf8264875
--- /dev/null
+++ b/+quantity/EquidistantDomain.m
@@ -0,0 +1,24 @@
+classdef EquidistantDomain < quantity.Domain
+ %EQUIDISTANTDOMAIN class to handle the discretization of the range of
+ %definition of a function. The discretization points are equally
+ %distributed over the domain.
+
+ properties
+ Property1
+ end
+
+ methods
+ function obj = untitled(inputArg1,inputArg2)
+ %UNTITLED Construct an instance of this class
+ % Detailed explanation goes here
+ obj.Property1 = inputArg1 + inputArg2;
+ end
+
+ function outputArg = method1(obj,inputArg)
+ %METHOD1 Summary of this method goes here
+ % Detailed explanation goes here
+ outputArg = obj.Property1 + inputArg;
+ end
+ end
+end
+
diff --git a/+unittests/+quantity/testDomain.m b/+unittests/+quantity/testDomain.m
new file mode 100644
index 0000000000000000000000000000000000000000..3f14eb8a4ca10a2f3192b285e2841e6a45a6d9f1
--- /dev/null
+++ b/+unittests/+quantity/testDomain.m
@@ -0,0 +1,45 @@
+function [tests] = testDomain()
+ tests = functiontests(localfunctions);
+end
+
+function setupOnce(testCase)
+end
+
+function testDomainInit(testCase)
+
+Z = linspace(0,pi, 3);
+d = quantity.Domain('name', 'z', 'grid', Z);
+D = [d d];
+
+testCase.verifyEqual( ndims(D), 2);
+testCase.verifyEqual( ndims(d), 1);
+end
+
+
+function testDomainGridLength(testCase)
+
+Z = linspace(0,pi, 3);
+d = quantity.Domain('name', 'z', 'grid', Z);
+D = [d d];
+
+testCase.verifyEqual( cellfun(@(v) numel(v), {Z, Z}), D.gridLength)
+
+end
+
+function testDomainNumGridElements(testCase)
+
+Z = linspace(0,pi, 3);
+d = quantity.Domain('name', 'z', 'grid', Z);
+D = [d d];
+
+testCase.verifyEqual( D.numGridElements, prod([length(Z), length(Z)]));
+
+end
+
+function testDomainEmpty(testCase)
+
+d = quantity.Domain();
+
+testCase.verifyTrue( isempty(d) )
+
+end
\ No newline at end of file