Commit e42f833b authored by Ferdinand Fischer's avatar Ferdinand Fischer
Browse files

Fixed changeGrid to work with quantity.EquidistantDomains

parent d359c4f8
......@@ -835,7 +835,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
% solves the first order ODE
% dvar / ds = obj(var(s))
% var(0) = ic
% to obtain var(s, ic) depending on both the argument s and the initial
% to obtain var(s, ic) depending on both the argument s and the initial
% condition ic. Herein, obj may only depend on one variable / gridName / ...
% domain.
assert(numel(obj(1).domain) == 1, ...
......@@ -1354,7 +1354,6 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
end
if obj(1).isConstant
newDomain(1:length( gridNew )) = quantity.Domain();
for it = 1 : length(gridNew)
newDomain(it) = ...
......@@ -1434,7 +1433,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
assert(numel(myGrid) > 1, 'If the state transition matrix is computed for constant values, a spatial domain has to be defined!')
myDomain = [quantity.Domain(gridName1, myGrid), ...
quantity.Domain(gridName2, myGrid)];
quantity.Domain(gridName2, myGrid)];
if obj.isConstant
% for a constant system matrix, the matrix exponential
......@@ -1624,7 +1623,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
% the function "gridJoin". It returns the finer grid.
% If the domain names do not coincide, they are just
% appended. At first, the domains of qunatity A, then the
% domains of quantity B.
% domains of quantity B.
joinedDomain = [ gridJoin(domainA(idx.A.common), domainB(idx.B.common)), ...
domainA(~idx.A.common), domainB(~idx.B.common) ];
......@@ -1933,9 +1932,9 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
end
if isa(diffGridName, 'quantity.Domain')
% a quantity.Domain is used instead of a grid name
% a quantity.Domain is used instead of a grid name
% -> get the grid name from the domain object
%
%
% #todo@domain: make the default case to call with a
% quantity.Domain instead of a grid name. Then, the
% section about the grid selection and so can be simplified
......@@ -1965,7 +1964,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
% obj.int() -> integrate over all dimensions of this
% quantity.
intDomain = obj(1).domain;
elseif nargin == 2
elseif nargin == 2
% obj.int(z) OR obj.int(z, a, b)
% integrate over the domain
if isa(varargin{1}, 'quantity.Domain')
......@@ -2060,7 +2059,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
end
% combine both domains with finest grid
joinedDomain = join(A(1).domain, B(1).domain);
joinedDomain = join(A(1).domain, B(1).domain);
[aDiscrete] = A.expandValueDiscrete(joinedDomain);
[bDiscrete] = B.expandValueDiscrete(joinedDomain);
......
......@@ -68,23 +68,42 @@ classdef (InferiorClasses = {?quantity.EquidistantDomain}) Domain < ...
methods (Access = public)
function d = split(obj, splittingPoints)
function d = split(obj, splittingPoints, optArgs)
%SPLIT splits the domain into subdomains
% d = split(OBJ, SPLITTINGPOINTS) splits the domain into the
% subdomains specified by SPLITTINGPOINTS. To be particular, a
% domain (0,1) called with the SPLITTINGPOINTS 0.5 will be
% split into two domains (0,0.5) and (0.5,1). To consider the
% case if the SPLITTINGPOINTS are not exactly contained in the
% grid, the closest grid point of OBJ.grid will be chosen as
% splitting point.
arguments
obj
splittingPoints double;
optArgs.AbsTol (1,1) double = 1e-3;
optArgs.warning (1,1) logical = true;
end
d(1:numel(splittingPoints)+1) = quantity.Domain();
idx = zeros(numel(splittingPoints)+1,1);
idx(1) = 1;
for k = 1:numel(splittingPoints)
% Search the grid point closest to the k-th splittingPoint
[~, idx(k+1)] = min(abs((obj.grid - splittingPoints(k))));
[val, idx(k+1)] = min( abs( (obj.grid - splittingPoints(k)) ));
delta = abs( obj.grid(idx(k+1)) - splittingPoints(k));
if delta > 1e-3
warning('The deviation of the grid from desired splitting point %d is %d\n', splittingPoints(k), delta);
if optArgs.warning
% verify that the distance of the chosen grid point to the
% k-th splitting point is within some tolerances:
delta = abs( obj.grid(idx(k+1)) - splittingPoints(k));
if delta > optArgs.AbsTol
warning('DOMAIN:split', ...
'The deviation of the grid from desired splitting point %d is %d\n', splittingPoints(k), delta);
end
end
d(k) = quantity.Domain(obj.name, obj.grid(idx(k):idx(k+1)));
d(k) = quantity.Domain(obj.name, obj.grid(idx(k):idx(k+1)));
end
d(k+1) = quantity.Domain(obj.name, obj.grid(idx(k+1):end));
......@@ -397,7 +416,6 @@ classdef (InferiorClasses = {?quantity.EquidistantDomain}) Domain < ...
end % methods (Public)
methods (Access = protected)
function s = getPropertyGroups(obj)
% Function to display the correct values
......@@ -426,7 +444,6 @@ classdef (InferiorClasses = {?quantity.EquidistantDomain}) Domain < ...
[sprintf('%i, ', obj.n) sprintf('\b\b\b')];
end
end % getPropertyGroups
end % methods (Access = protected)
methods (Static)
......@@ -449,19 +466,24 @@ classdef (InferiorClasses = {?quantity.EquidistantDomain}) Domain < ...
end
end % gridCells2domain()
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
function d = defaultDomain(gridSize, name)
%DEFAULTDOMAIN generate default domain
% d = defaultDomain() will create a default domain with a grid
% from 0 to 1 with 100 points and the default name x_1
% d = defaultDomain(GRIDSIZE) will create a default domain with
% a default grid of size specified in GRIDSIZE. The domains
% will have the default names
% x_1, x_2, ...
% d = defaultDomain(GRIDSIZE, NAME) will create a default
% domain with a default grid of size specified in GRIDSIZE. The
% domains will have the names specified in NAME.
arguments
gridSize double = 100;
name string = quantity.Domain.defaultGridNames(gridSize);
end
% generate a default gird with given sizes
g = quantity.Domain.empty();
d(1:length(gridSize)) = quantity.Domain();
for k = 1:length(gridSize)
o = ones(1, length(gridSize) + 1); % + 1 is required to deal with one dimensional grids
......@@ -469,10 +491,28 @@ classdef (InferiorClasses = {?quantity.EquidistantDomain}) Domain < ...
O = ones(o);
O(:) = linspace(0, 1, gridSize(k));
g(k) = quantity.Domain(name{k}, O);
d(k) = quantity.Domain(name{k}, O);
end
end
function g = defaultGrid(varargin)
%DEFAULTGRID generate default grid for a quantity
warning('DEPRICATED! Use quantity.Domain.defaultDomain instead.');
g = quantity.Domain.defaultDomain(varargin{:});
end
function name = defaultGridNames(gridSize)
%DEFAULTGRIDNAMES generates default grid names
% name = defaultGridNames(GRIDSIZE) generates default grid
% names for a grid of size GRIDSIZE. The default names will be
% x_1, x_2, ...
name = cell(1, length(gridSize));
for k = 1:length(gridSize)
name{k} = ['x_' num2str(k)];
end
end
function [myDomain, unmatched] = parser(varargin)
%% domain parser
domainParser = misc.Parser();
......
......@@ -1374,12 +1374,28 @@ quan = quantity.Discrete({sin(z * t * pi); cos(z * t * pi)}, 'grid', {z, t}, 'gr
gridSampled = {linspace(0, 1, 11), linspace(0, 1, 21)};
quanCopy = copy(quan);
% test the change of two grids:
quanSampled = quanCopy.changeGrid(gridSampled, {'z', 't'});
testCase.verifyEqual(quanSampled.on(), quan.on(gridSampled))
% test the change of only one grid:
newZ = linspace(0,1,5);
qNewZ = quanCopy.changeGrid( newZ, 'z')
testCase.verifyEqual(qNewZ.on(), quan.on({newZ, t}));
% test the change of grids in wrong order
quanCopy2 = copy(quan);
quanSampled2 = quanCopy2.changeGrid(gridSampled, {'t', 'z'});
testCase.verifyEqual(quanSampled2.on(), permute(quan.on(gridSampled, {'t', 'z'}), [2, 1, 3]));
% test the change of a quantity.EquidistantDomain with a quantity.Domain
e = quantity.EquidistantDomain( 'z', 0, 1);
d = quantity.Domain('z', linspace(0,1, 3));
E = quantity.Discrete( sin(e.grid), 'domain', e);
E_ = E.changeGrid(d);
testCase.verifyEqual(E_.on(), E.on(d))
end
function testSubs(testCase)
......
......@@ -29,7 +29,10 @@ testCase.verifyEqual(b(2).grid', -1:0)
testCase.verifyEqual(b(3).grid', 0:4)
testCase.verifyEqual(b(4).grid', 4:15)
c = a.split([0.1, 2.6]);
testCase.verifyWarning(@() a.split([0.1, 2.6]), 'DOMAIN:split')
testCase.verifyWarningFree(@() a.split([0.1, 2.6], 'AbsTol', 0.4), 'DOMAIN:split')
c = a.split([0.1, 2.6], 'warning', false);
testCase.verifyEqual(c(1).grid', -5:0)
testCase.verifyEqual(c(2).grid', 0:3)
......@@ -177,10 +180,15 @@ tc.verifyTrue( isempty([o d]) )
end
function testDefaultGrid(tc)
g = quantity.Domain.defaultGrid([100, 42]);
tc.verifyEqual(g(1).grid, linspace(0, 1, 100).');
tc.verifyEqual(g(2).grid, linspace(0, 1, 42).');
function testDefaultGrid(testCase)
g = quantity.Domain.defaultDomain([100, 42]);
testCase.verifyEqual(g(1).grid, linspace(0, 1, 100).');
testCase.verifyEqual(g(2).grid, linspace(0, 1, 42).');
g = quantity.Domain.defaultDomain();
testCase.verifyEqual(g.grid, linspace(0, 1).');
testCase.verifyEqual(g.name, "x_1")
end
function testGridJoin(tc)
......
......@@ -36,4 +36,8 @@ D = [d; e];
testCase.verifyEqual(E(1), quantity.Domain(e))
testCase.verifyEqual(D(2), quantity.Domain(e))
EE(1:2) = quantity.Domain();
EE(1) = e;
EE(2) = d;
end
\ No newline at end of file
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment