Commit e37e9a54 authored by Jakob Gabriel's avatar Jakob Gabriel
Browse files

quantity.Discrete: simplified implementation of power()

+ wrote unittests (also for quantity.Symbolic)
parent 8518e01c
......@@ -1609,24 +1609,10 @@ classdef (InferiorClasses = {?quantity.Symbolic}) Discrete ...
end % sqrtm()
function P = power(obj, p)
% a.^p elementwise power implemented by multiplication
% #todo unittest required
P = quantity.Discrete.zeros(size(obj), obj(1).domain);
if round(p) == p
% power of integer order:
for k = 1:numel(obj)
P(k) = obj(k)^p;
end
else
% power of not integer order
for k = 1:numel(obj)
P(k) = obj(k).copy();
P(k).valueDiscrete = obj(k).valueDiscrete.^p;
P(k).name = obj(k).name + ".^" + num2str(p);
end
end
end % power
% a.^p elementwise power
P = quantity.Discrete(obj.on().^(p), obj(1).domain, ...
"name", obj(1).name + ".^{" + num2str(p) + "}");
end % power()
function P = mpower(a, p)
% Matrix power a^p is matrix or scalar a to the power p.
......
......@@ -4,6 +4,29 @@ function [tests] = testDiscrete()
tests = functiontests(localfunctions);
end
function testPower(tc)
%% scalar case
z = quantity.Domain("z", linspace(0, 1, 11));
A = quantity.Discrete(sin(z.grid * pi), z);
aa = sin(z.grid * pi) .* sin(z.grid * pi);
AA = A.^2;
tc.verifyEqual( aa, AA.on());
tc.verifyEqual( on(A.^0), ones(z.n, 1));
tc.verifyEqual( on(A.^3), on(AA * A), 'AbsTol', 1e-15);
%% matrix case
zeta = quantity.Domain("zeta", linspace(0, 1, 21));
B = [0*z.Discrete + zeta.Discrete, A + 0*zeta.Discrete; ...
z.Discrete+zeta.Discrete, A * subs(A, "z", "zeta")];
BBB = B.^3;
BBB_alt = B .* B .* B;
tc.verifyEqual(BBB.on(), BBB_alt.on(), 'AbsTol', 1e-14);
tc.verifyEqual(BBB.at({0.5, 0.6}), (B.at({0.5, 0.6})).^3, 'AbsTol', 1e-15);
tc.verifyEqual(on(B.^0), ones([B(1).domain(1).n, B(1).domain(2).n, 2, 2]), ...
'AbsTol', 1e-15);
end % testPower()
function setupOnce(testCase)
syms z zeta t sy
testCase.TestData.sym.z = z;
......
......@@ -13,12 +13,35 @@ B = [0*z.Symbolic + zeta.Symbolic, A + 0*zeta.Symbolic; ...
B2 = [zeta.Symbolic + 0*z.Symbolic, A + 0*zeta.Symbolic; ...
z.Symbolic+zeta.Symbolic, A * subs(A, "z", "zeta")];
tc.verifyEqual(B.on([z, zeta]), B2.on([z, zeta]));
% TODO: fix this error!
% TODO: fix this error! See issue #41
% this problem is maybe caused as input arguments of B2(1).valueContinuous is somehow ordered
% differently compared to B2(2:4).valueContinuous or all B(:).valueContinuous,
% i.e. B2(1).valueContinuous is (zeta, z) instead of (z, zeta).
end % testValueContinuousBug()
function testPower(tc)
%% scalar case
z = quantity.Domain("z", linspace(0, 1, 11));
A = quantity.Symbolic(sin(sym("z") * pi), z);
aa = sin(z.grid * pi) .* sin(z.grid * pi);
AA = A.^2;
tc.verifyEqual( aa, AA.on());
tc.verifyEqual( on(A.^0), ones(z.n, 1));
tc.verifyEqual( on(A.^3), on(AA * A), 'AbsTol', 1e-15);
%% matrix case
zeta = quantity.Domain("zeta", linspace(0, 1, 21));
B = [0*z.Symbolic + zeta.Symbolic, A + 0*zeta.Symbolic; ...
z.Symbolic+zeta.Symbolic, A * subs(A, "z", "zeta")];
BBB = B.^3;
BBB_alt = B .* B .* B;
tc.verifyEqual(BBB.on(), BBB_alt.on(), 'AbsTol', 1e-14);
tc.verifyEqual(BBB.at({0.5, 0.6}), (B.at({0.5, 0.6})).^3, 'AbsTol', 1e-15);
tc.verifyEqual(on(B.^0), ones([B(1).domain(1).n, B(1).domain(2).n, 2, 2]), ...
'AbsTol', 1e-15);
end % testPower()
function testMPower(tc)
%% scalar case
z = quantity.Domain("z", linspace(0, 1, 11));
......
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