Commit 56f57e74 authored by Ferdinand Fischer's avatar Ferdinand Fischer
Browse files

Added symbolic evaluation of quantity.Symbolic. (This is necessary for some...

Added symbolic evaluation of quantity.Symbolic. (This is necessary for some strange functions, see e.g. the corresponding unittest)
Simplified the "isempty" check for quantities
parent fd389244
classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mixin.Copyable
classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mixin.Copyable & matlab.mixin.CustomDisplay
properties (SetAccess = protected)
% Discrete evaluation of the continuous quantity
valueDiscrete double;
......@@ -8,7 +8,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% In this cell, already computed derivatives can be stored to avoid
% multiple computations of the same derivative.
derivatives cell = {};
derivatives;
end
properties (Hidden, Access = protected, Dependent)
......@@ -196,7 +196,7 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
function valueDiscrete = get.valueDiscrete(obj)
% check if the value discrete for this object
% has already been computed.
empty = reshape(cellfun('isempty', {obj(:).valueDiscrete}), size(obj));
empty = isempty(obj.valueDiscrete);
if any(empty(:))
obj.valueDiscrete = obj.on(obj.grid);
end
......@@ -1595,12 +1595,12 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete < handle & matlab.mi
% check if there is already an object
empty = any(size(obj) == 0);% || isempty(obj(1).grid);
% if there is an object, check if the object has a nonempty
% valueContinuous
if ~empty
empty = cellfun('isempty', {obj(:).valueDiscrete});
empty = any(empty(:));
end
% if there is an object, check if the object has a nonempty
% % valueContinuous
% if ~empty
% empty = cellfun('isempty', {obj(:).valueDiscrete});
% empty = any(empty(:));
% end
end % isempty()
function P = ztzTimes(a, b)
......
......@@ -7,6 +7,9 @@ classdef Symbolic < quantity.Function
properties (Constant)
defaultSymVar = sym('z', 'real');
end
properties
symbolicEvaluation = false;
end
methods
......@@ -22,6 +25,7 @@ classdef Symbolic < quantity.Function
variableParser = misc.Parser();
variableParser.addParameter('variable', quantity.Symbolic.getVariable(valueContinuous));
variableParser.addParameter('symbolicEvaluation', false);
variableParser.parse(varargin{:});
variable = variableParser.Results.variable;
numVar = numel(variable);
......@@ -94,6 +98,8 @@ classdef Symbolic < quantity.Function
for k = 1:numel(valueContinuous)
obj(k).variable = variable; %#ok<AGROW>
obj(k).valueSymbolic = symb{k}; %#ok<AGROW>
obj(k).symbolicEvaluation = ...
variableParser.Results.symbolicEvaluation; %#ok<AGROW>
end
end
end
......@@ -147,13 +153,6 @@ classdef Symbolic < quantity.Function
'grid', myParser.Results.grid, 'gridName', myParser.Results.gridName, ...
'name', myParser.Results.name);
end
function Q = quantity.SymbolicII(obj)
for i = 1:numel(obj)
argin = misc.struct2namevaluepair(struct(obj(i)));
Q(i) = quantity.SymbolicII(obj(i).sym, argin{:});
end
Q = reshape(Q, size(obj));
end
function f = function_handle(obj)
f = matlabFunction(sym(obj));
end
......@@ -773,6 +772,17 @@ classdef Symbolic < quantity.Function
methods (Access = protected)
function v = evaluateFunction(obj, varargin)
% Evaluates the symbolic expression of this quantity. If the
% flag symbolicEvaluation is set, the expression is evaluated
% with subs(...) otherwise a function handle is used.
if obj.symbolicEvaluation
v = double(subs(sym(obj), obj.variable, varargin{:}));
else
v = evaluateFunction@quantity.Function(obj, varargin{:});
end
end
function p = getCopyParameters(obj)
p = nameValuePair(obj, 'gridName');
end
......
classdef SymbolicII < quantity.Symbolic
methods
function obj = SymbolicII(valueContinuous, varargin)
if nargin == 0
parentArgin = {};
else
parentArgin = [{valueContinuous}, varargin(:)'];
end
obj@quantity.Symbolic(parentArgin{:});
end
end % methods
%% mathematical operations
methods (Access = protected)
function v = evaluateFunction(obj, varargin)
v = double(subs(sym(obj), obj.variable, varargin{:}));
end
end
methods (Static)
function P = zeros(valueSize, varargin)
%ZEROS initializes an zero quantity.Discrete object of the size
%specified by the input valueSize.
P = quantity.SymbolicII(zeros(valueSize), varargin{:});
end
end
end
......@@ -604,13 +604,13 @@ end
function testCumInt(testCase)
tt = linspace(pi, 1.1*pi, 51)';
s = tt;
[T, S] = ndgrid(tt, tt);
tGrid = linspace(pi, 1.1*pi, 51)';
s = tGrid;
[T, S] = ndgrid(tGrid, tGrid);
syms sy tt
syms sy t
a = [ 1, sy; tt, 1];
a = [ 1, sy; t, 1];
b = [ sy; 2*sy];
A = zeros([size(T), size(a)]);
......@@ -618,7 +618,7 @@ B = zeros([length(s), size(b)]);
for i = 1:size(a,1)
for j = 1:size(a,2)
A(:,:,i,j) = subs(a(i,j), {sy, tt}, {S, T});
A(:,:,i,j) = subs(a(i,j), {sy, t}, {S, T});
end
end
......@@ -628,12 +628,11 @@ end
%% int_0_t a(t,s) * b(s) ds
% compute symbolic version of the volterra integral
v = int(a*b, sy, tt(1), tt);
V = quantity.Symbolic(v, 'grid', tt, 'gridName', 'tt');
V = V.subs('tt', 't');
v = int(a*b, sy, tGrid(1), t);
V = quantity.Symbolic(v, 'grid', tGrid, 'gridName', 't');
% compute the numeric version of the volterra integral
k = quantity.Discrete(A, 'size', size(a), 'grid', {tt, s}, 'gridName', {'t', 's'});
k = quantity.Discrete(A, 'size', size(a), 'grid', {tGrid, s}, 'gridName', {'t', 's'});
x = quantity.Discrete(B, 'size', size(b), 'grid', {s}, 'gridName', 's');
f = cumInt(k * x, 's', s(1), 't');
......@@ -641,63 +640,62 @@ f = cumInt(k * x, 's', s(1), 't');
testCase.verifyEqual(V.on(), f.on(), 'AbsTol', 1e-5);
%% int_s_t a(t,s) * b(s) ds
v = int(a*b, sy, sy, tt);
V = quantity.Symbolic(subs(v, {tt, sy}, {'t', 's'}), 'grid', {tt, s});
v = int(a*b, sy, sy, t);
V = quantity.Symbolic(subs(v, {t, sy}, {'t', 's'}), 'grid', {tGrid, s});
f = cumInt(k * x, 's', 's', 't');
testCase.verifyEqual( f.on(f(1).grid, f(1).gridName), V.on(f(1).grid, f(1).gridName), 'AbsTol', 1e-5 );
%% int_s_t a(t,s) * c(t,s) ds
c = [1, sy+1; tt+1, 1];
c = [1, sy+1; t+1, 1];
C = zeros([size(T), size(c)]);
for i = 1:numel(c)
C(:,:,i) = double(subs(c(i), {tt sy}, {T S}));
C(:,:,i) = double(subs(c(i), {t sy}, {T S}));
end
y = quantity.Discrete(C, 'size', size(c), 'grid', {tt s}, 'gridName', {'t' 's'});
y = quantity.Discrete(C, 'size', size(c), 'grid', {tGrid s}, 'gridName', {'t' 's'});
v = int(a*c, sy, sy, tt);
V = quantity.Symbolic(subs(v, {tt, sy}, {'t', 's'}), 'grid', {tt, s});
v = int(a*c, sy, sy, t);
V = quantity.Symbolic(subs(v, {t, sy}, {'t', 's'}), 'grid', {tGrid, s});
f = cumInt( k * y, 's', 's', 't');
testCase.verifyEqual( f.on(f(1).grid, f(1).gridName), V.on(f(1).grid, f(1).gridName), 'AbsTol', 1e-5 );
%% testCumIntWithLowerAndUpperBoundSpecified
tt = linspace(pi, 1.1*pi, 51)';
s = tt;
[T, S] = ndgrid(tt, tt);
tGrid = linspace(pi, 1.1*pi, 51)';
s = tGrid;
[T, S] = ndgrid(tGrid, tGrid);
syms sy tt
a = [ 1, sy; tt, 1];
a = [ 1, sy; t, 1];
A = zeros([size(T), size(a)]);
for i = 1:size(a,1)
for j = 1:size(a,2)
A(:,:,i,j) = subs(a(i,j), {sy, tt}, {S, T});
A(:,:,i,j) = subs(a(i,j), {sy, t}, {S, T});
end
end
% compute the numeric version of the volterra integral
k = quantity.Discrete(A, 'size', size(a), 'grid', {tt, s}, 'gridName', {'t', 's'});
k = quantity.Discrete(A, 'size', size(a), 'grid', {tGrid, s}, 'gridName', {'t', 's'});
fCumInt2Bcs = cumInt(k, 's', 'zeta', 't');
fCumInt2Cum = cumInt(k, 's', tt(1), 't') ...
- cumInt(k, 's', tt(1), 'zeta');
fCumInt2Cum = cumInt(k, 's', tGrid(1), 't') ...
- cumInt(k, 's', tGrid(1), 'zeta');
testCase.verifyEqual(fCumInt2Bcs.on(), fCumInt2Cum.on(), 'AbsTol', 10*eps);
%% testCumInt with one independent variable
tt = linspace(0, pi, 51)';
tGrid = linspace(0, pi, 51)';
syms t
a = [ 1, sin(t); t, 1];
% compute the numeric version of the volterra integral
f = quantity.Symbolic(a, 'grid', {tt}, 'variable', {'t'});
f = quantity.Symbolic(a, 'grid', {tGrid}, 'variable', {'t'});
f = quantity.Discrete(f);
intK = cumInt(f, 't', 0, 't');
K = quantity.Symbolic( int(a, 0, t), 'grid', {tt}, 'variable', {'t'});
K = quantity.Symbolic( int(a, 0, t), 'grid', {tGrid}, 'variable', {'t'});
testCase.verifyEqual(intK.on(), K.on(), 'AbsTol', 1e-3);
......
......@@ -4,6 +4,18 @@ function [tests ] = testSymbolic()
tests = functiontests(localfunctions());
end
function testSymbolicEvaluation(tc)
syms z
myGrid = linspace(0, 1, 7);
fS = quantity.Symbolic(z * sinh(z * 1e5) / cosh(z * 1e5), ...
'variable', {z}, 'grid', {myGrid}, 'symbolicEvaluation', true);
fF = quantity.Symbolic(z * sinh(z * 1e5) / cosh(z * 1e5), ...
'variable', {z}, 'grid', {myGrid}, 'symbolicEvaluation', false);
tc.verifyTrue(any(isnan(fF.on)))
tc.verifyFalse(any(isnan(fS.on)))
end
function testCtranspose(tc)
syms z zeta
qSymbolic = quantity.Symbolic(...
......
......@@ -30,7 +30,6 @@ verifyEqual(testCase, si.on(), SI.on());
end
function testGevreyBump(testCase)
% Test script to check if the class for the bumb function returns the right
% signals.
......
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