testSymbolic.m 35.4 KB
Newer Older
1
2
3
4
5
6
function [tests ] = testSymbolic()
%TESTQUANTITY Unittests for the quantity objects
%   Add all unittests for the quantity files here
tests = functiontests(localfunctions());
end

Jakob Gabriel's avatar
Jakob Gabriel committed
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
function testSubsDomain(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
zeta = quantity.Domain("zeta", linspace(0, 1, 11));
zLr = quantity.Domain("z", linspace(0, 1, 5));
eta = quantity.Domain("eta", linspace(0, 1, 11));
f = z.Symbolic() * zeta.Symbolic();
%
fofZEta = f.subsDomain("zeta", eta);	% = f(z, eta)
tc.verifyEqual(fofZEta.on(), f.on());
tc.verifyEqual([fofZEta(1).domain.name], ["z", "eta"]);
%
fofZlrZeta = f.subsDomain("z", zLr);	% = f(z, zeta)
tc.verifyEqual(fofZlrZeta.on(), f.on([zLr, zeta]));
tc.verifyEqual([fofZlrZeta(1).domain.name], ["z", "zeta"]);

%
fofZ = f.subsDomain("zeta", zLr);		% = f(z) 
tc.verifyEqual(fofZ.on(), diag(f.on([z, quantity.Domain("zeta", z.grid)])));
tc.verifyEqual([fofZ(1).domain.name], "z");
end % testSubsDomain()

function testSubsNumeric(tc)
z = quantity.Domain("z", linspace(0, 1, 5));
zeta = quantity.Domain("zeta", linspace(0, 1, 4));
t = quantity.Domain("t", linspace(0, 1, 6));

quanScalar = z.Symbolic() + 1 + zeta.Symbolic() * t.Symbolic();
tc.verifyEqual(quanScalar.subsNumeric(["z", "t"], [0, 0]).on(), ones(zeta.n, 1));
tc.verifyEqual(quanScalar.subsNumeric(["t", "z", "zeta"], [0.1, 0.2, 0.3]), 0.2 + 1 + 0.3 * 0.1);

quanMatrix = quantity.Symbolic(...
	[sym("z") + 1 + sym("zeta") * sym("t"), 1; ...
	sym("z") + sym("t"), -sym("zeta")], [zeta, t, z]);
tc.verifyEqual(quanMatrix.subsNumeric(["z", "zeta", "t"], [0.1, 0.2, 0.3]), ...
	[0.1 + 1 + 0.2 * 0.3, 1; 0.1 + 0.3, -0.2], "AbsTol", 10*eps);
tc.verifyEqual(quanMatrix.subsNumeric(["z", "zeta", "t"], [0, 0, 0]), [1, 1; 0, 0]);
tc.verifyEqual(quanMatrix.subsNumeric(["z", "zeta", "t"], [1, 1, 1]), [3, 1; 2, -1]);
tc.verifyEqual(quanMatrix(1, 1).subsNumeric(["z", "zeta", "t"], [0.1, 0.2, 0.3]), ...
	0.1 + 1 + 0.2 * 0.3, "AbsTol", 10*eps);
end % testSubsNumeric()

function testChangeDomain(tc)
zeta = quantity.Domain("zeta", linspace(0, 1, 4));
t = quantity.Domain("t", linspace(0, 1, 4));
thisSym = sym("zeta") + sym("t");

quan = quantity.Symbolic(thisSym, [zeta, t]);
newDomain = [quantity.Domain("zeta", linspace(0, 1, 11)), quantity.Domain("t", linspace(0, 1, 21))];
newQuan = quan.changeDomain(newDomain);

tc.verifyEqual(newQuan.on(), quan.on(newDomain));
tc.verifyTrue(isa(newQuan, "quantity.Symbolic"));
end % testChangeDomain()

61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
function testOn2(tc)
zeta = quantity.Domain("zeta", linspace(0, 1, 4));
t = quantity.Domain("t", linspace(0, 1, 4));
thisSym = sym("zeta") + sym("t");
quan = quantity.Symbolic(thisSym, [zeta, t]);

solution = 0.2 + 0.3; % this obtained by manual calculation
% verify solution
tc.verifyEqual(double(subs(subs(thisSym, "zeta", 0.2), "t", 0.3)), solution, "AbsTol", 10*eps);

% test on with numeric values only
tc.verifyEqual(quan.on({0.2, 0.3}), solution);

% test on with domain
tc.verifyEqual(quan.on([quantity.Domain("zeta", 0.2), quantity.Domain("t", 0.3)]), solution);
end % testOn2()

function testOn2b(tc)
zeta = quantity.Domain("zeta", linspace(0, 1, 4));
t = quantity.Domain("t", linspace(0, 1, 4));
thisSym = sym("zeta") + sym("t");
quan = quantity.Symbolic(thisSym, [t, zeta]);

solution = 0.2 + 0.3; % this obtained by manual calculation
% verify solution
tc.verifyEqual(double(subs(subs(thisSym, "zeta", 0.2), "t", 0.3)), solution, "AbsTol", 10*eps);

% test on with numeric values only
tc.verifyEqual(quan.on({0.3, 0.2}), solution);

% test on with domain
tc.verifyEqual(quan.on([quantity.Domain("zeta", 0.2), quantity.Domain("t", 0.3)]), solution);
end % testOn2()

function testOn3(tc)
z = quantity.Domain("z", linspace(0, 1, 4));
t = quantity.Domain("t", linspace(0, 1, 4));
thisSym = sym("z") + sym("t");
quan = quantity.Symbolic(thisSym, [z, t]);

solution = 0.2 + 0.3; % this obtained by manual calculation
% verify solution
tc.verifyEqual(double(subs(subs(thisSym, "z", 0.2), "t", 0.3)), solution, "AbsTol", 10*eps);

% test on with numeric values only
tc.verifyEqual(quan.on({0.2, 0.3}), solution);

% test on with domain
tc.verifyEqual(quan.on([quantity.Domain("z", 0.2), quantity.Domain("t", 0.3)]), solution);
end % testOn2()

112
113
114
115
116
117
118
119
120
121
122
123
124
125
function testDiffOfConstants(tc)
% unittest to fix a bug...
z = quantity.Domain("z", linspace(0, 1, 4));
K0 = quantity.Symbolic.zeros(1, z);%[z, zeta]);
tc.verifyEqual(K0.on(), zeros([4, 1]));

K1 = quantity.Symbolic.ones(1, z);%[z, zeta]);
tc.verifyEqual(K1.diff("z", 1).on(), zeros([4, 1]));

zeta = quantity.Domain("zeta", linspace(0, 1, 5));
K2 = quantity.Symbolic.ones([2, 2], [z, zeta]);
tc.verifyEqual(K2.diff("z", 1).on(), zeros([4, 5, 2, 2]));
end % testDiffOfConstants

126
function testZeroOn(tc)
127
% unittest to fix a bug related with isNumber() usage in on() of symbolic
128
z = quantity.Domain("z", linspace(0, 1, 11));
129
quan = 0 * quantity.Symbolic(1-0.5*sym("z"), z, "name", "Gamma");
130
131
132
tc.verifyEqual(quan.on({linspace(0, 1, 5)}), zeros(5, 1));
tc.verifyEqual(quan.on(linspace(0, 1, 5)), zeros(5, 1));
end % testZeroOn()
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151

function testOnes(tc)
myDomain = [quantity.Domain("z", linspace(0, 1, 11)), quantity.Domain("zeta", linspace(0, 1, 5))];
tc.verifyEqual(MAX(abs(quantity.Symbolic.ones(1, myDomain) - 1)), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.ones(2, myDomain) - ones(2))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.ones([2, 1], myDomain) - ones([2, 1]))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.ones([2, 2], myDomain(2)) - ones([2, 2]))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.ones([4, 4], myDomain) - ones(4))), 0);
end % testOnes()

function testEye(tc)
myDomain = [quantity.Domain("z", linspace(0, 1, 11)), quantity.Domain("zeta", linspace(0, 1, 5))];
tc.verifyEqual(MAX(abs(quantity.Symbolic.eye(1, myDomain) - 1)), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.eye(2, myDomain) - eye(2))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.eye([2, 1], myDomain) - eye([2, 1]))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.eye([2, 2], myDomain(2)) - eye([2, 2]))), 0);
tc.verifyEqual(MAX(abs(quantity.Symbolic.eye([4, 4], myDomain) - eye(4))), 0);
end % testEye()

152
153
154
155
function testValueContinuousBug(tc)
z = quantity.Domain("z", linspace(0, 1, 7));
zeta = quantity.Domain("zeta", linspace(0, 1, 5));
A = quantity.Symbolic(sin(sym("z") * pi), z);
156

157
B1 = [0*z.Symbolic + zeta.Symbolic, A + 0*zeta.Symbolic; ...
158
159
160
	z.Symbolic+zeta.Symbolic, A * subs(A, "z", "zeta")];
B2 = [zeta.Symbolic + 0*z.Symbolic, A + 0*zeta.Symbolic; ...
	z.Symbolic+zeta.Symbolic, A * subs(A, "z", "zeta")];
161
162
tc.verifyEqual(B1.on([z, zeta]), B2.on([z, zeta]));

163
164
165
x = A * subs(A, "z", "zeta");
C1 = [0*z.Symbolic + zeta.Symbolic; x];
C2 = [zeta.Symbolic + 0*z.Symbolic; x];
166
167
tc.verifyEqual(C1.on([z, zeta]), C2.on([z, zeta]));

168
e = zeta.Symbolic + 0*z.Symbolic;
169
e = e.sortDomain();
170
f = e.changeDomain([z z.rename("zeta")]);
171
172


173
% TODO: fix this error! See issue #41
174
175
176
177
178
% 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()

179
180
181
182
183
184
185
186
187
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));
188
tc.verifyEqual( on(A.^3), on(AA * A), "AbsTol", 1e-15);
189
190
191
192
193
194
195

%% 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;
196
197
198
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);
199
200
end % testPower()

201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
function testMPower(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));

%% 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")];
B_discrete = quantity.Discrete(B);
BBB = B^3;
BBB_alt = B * B * B;
BBB_discrete = B_discrete^3;
220
221
222
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(BBB_discrete.at({0.5, 0.6}), (B.at({0.5, 0.6}))^3, "AbsTol", 1e-15);
223
tc.verifyEqual(on(B^0), permute(repmat(eye(2), [1, 1, B(1).domain(1).n, B(1).domain(2).n]), [3, 4, 1, 2]), ...
224
	"AbsTol", 1e-15);
225
226
end % testMPower()

227
228
229
function testNorm(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
quan = Symbolic(z) * [-1; 2; 3];
230
231
tc.verifyEqual(quan.norm.on(), on(z.Discrete*sqrt(1+4+9)), "AbsTol", 1e-15);
tc.verifyEqual(quan.norm.on(), on(quan.norm(2)), "AbsTol", 1e-15);
232
end % testNorm()
233

234
235
236
237
function testLog(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
quan = Symbolic(z) * ones(2, 2) + [1, 1; 0, 0]+0.1*ones(2,2);

238
239
tc.verifyEqual(log(quan.on()), on(log(quan)), "AbsTol", 1e-15);
tc.verifyEqual(on(quan), on(log(exp(quan))), "AbsTol", 1e-15);
240
241
242
243
244
245
end % testLog()

function testLog10(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
quan = Symbolic(z) * ones(2, 2) + [1, 1; 0, 0]+0.1*ones(2,2);

246
tc.verifyEqual(log10(quan.on()), on(log10(quan)), "AbsTol", 1e-15);
247
248
end % testLog10()

249
250
function testDet(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
251
zSymbolic = quantity.Symbolic(sym("z"), z);
252
253
254
255
256
257
258
myMatrix2x2 = [1, 2; 3, 4] * zSymbolic;
myMatrix3x3 = [1, 2, 3; 2, -3, 4; 3, 4, 5] * zSymbolic;
myMatrix5x5 = [1, 2, 3, 4, 5; ...
				2, 3, 4, 5, 6; ...
				3, 4, 5, 6, 7; ...
				4, -6, 7, 8, 9; ...
				1, 1, 1, 1, 1] * zSymbolic;
259
260
261
262
tc.verifyEqual(zSymbolic.det.on(), zSymbolic.on(), "AbsTol", 10*eps);
tc.verifyEqual(myMatrix2x2.det.at(0.5), det(myMatrix2x2.at(0.5)), "AbsTol", 10*eps);
tc.verifyEqual(myMatrix3x3.det.at(0.5), det(myMatrix3x3.at(0.5)), "AbsTol", 10*eps);
tc.verifyEqual(myMatrix5x5.det.at(0.5), det(myMatrix5x5.at(0.5)), "AbsTol", 10*eps);
263
264
265
266
end % testDet()

function testSum(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
267
myOnes = ones(3, 4) * (1 + 0*quantity.Symbolic(sym("z"), z));
268
269
270
271
tc.verifyEqual(on(sum(myOnes)), 3*4*ones(11, 1), "AbsTol", 10*eps);
tc.verifyEqual(on(sum(myOnes, 2)), 4*ones(11, 3), "AbsTol", 10*eps);
tc.verifyEqual(on(sum(myOnes, 1)), 3*ones(11, 1, 4), "AbsTol", 10*eps);
tc.verifyEqual(on(sum(myOnes, [1, 2])), 3*4*ones(11, 1), "AbsTol", 10*eps);
272

273
myTwos3D = ones(2, 3, 4) * 2 * (1 + 0*quantity.Symbolic(sym("z"), z));
274
275
276
tc.verifyEqual(on(sum(myTwos3D)), 2*3*4*ones(11, 1)*2, "AbsTol", 10*eps);
tc.verifyEqual(on(sum(myTwos3D, 2)), 3*ones(11, 2, 1, 4)*2, "AbsTol", 10*eps);
tc.verifyEqual(on(sum(myTwos3D, [1, 3])), 2*4*ones(11, 1, 3)*2, "AbsTol", 10*eps);
277
278
279
280
end % testSum()

function testProd(tc)
z = quantity.Domain("z", linspace(0, 1, 11));
281
myTwos = 2 * ones(3, 4) * (1 + 0*quantity.Symbolic(sym("z"), z));
282
283
284
285
tc.verifyEqual(on(prod(myTwos)), 2^(3*4)*ones(11, 1), "AbsTol", 10*eps);
tc.verifyEqual(on(prod(myTwos, 2)), 2^4*ones(11, 3), "AbsTol", 10*eps);
tc.verifyEqual(on(prod(myTwos, 1)), 2^3*ones(11, 1, 4), "AbsTol", 10*eps);
tc.verifyEqual(on(prod(myTwos, [1, 2])), 2^(3*4)*ones(11, 1), "AbsTol", 10*eps);
286

287
myTwos3D = quantity.Symbolic(2*(1+0*sym("z")) * ones(2, 3, 4), z);
288
289
290
tc.verifyEqual(on(prod(myTwos3D)), 2^(2*3*4)*ones(11, 1), "AbsTol", 10*eps);
tc.verifyEqual(on(prod(myTwos3D, 2)), 2^3*ones(11, 2, 1, 4), "AbsTol", 10*eps);
tc.verifyEqual(on(prod(myTwos3D, [1, 3])), 2^(2*4)*ones(11, 1, 3), "AbsTol", 10*eps);
291
292
end % testProd()

293
294
295
function testVec2Diag(testCase)
% quantity.Symbolic
syms z
296
Z = quantity.Domain("z", linspace(0, 1, 3));
297
298
myMatrixSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1, 0; 0, 0.9-z/2], Z);
myVectorSymbolic = quantity.Symbolic([sin(0.5*z*pi)+1; 0.9-z/2], Z);
299
300
301
302

testCase.verifyTrue( myVectorSymbolic.vec2diag.near(myMatrixSymbolic) );
end

303
304
305
306
function testCastSymbolic2Function(testCase)

z = linspace(0, 2*pi)';
Z = sym("z");
307
s = quantity.Symbolic([sin(Z); cos(Z)], quantity.Domain("z", z));
308
309
310
311
312
313
314

f = quantity.Function(s);

testCase.verifyTrue(all(size(f) == size(s)));

end

315
316
317
function testSymbolicEvaluation(tc)
syms z
myGrid = linspace(0, 1, 7);
318
319
myDomain = quantity.Domain("z", myGrid);
fS = quantity.Symbolic(z * sinh(z * 1e5) / cosh(z * 1e5), myDomain, ...
320
	"symbolicEvaluation", true);
321
fF = quantity.Symbolic(z * sinh(z * 1e5) / cosh(z * 1e5), myDomain, ...
322
	"symbolicEvaluation", false);
323
324
325
326
327

tc.verifyTrue(any(isnan(fF.on)))
tc.verifyFalse(any(isnan(fS.on)))
end

328
329
function testCtranspose(tc)
syms z zeta
330
331
332
333

dom.z = quantity.Domain("z", linspace(0, 1, 21));
dom.zeta = quantity.Domain("zeta", linspace(0, 1, 41));

334
qSymbolic = quantity.Symbolic( [1+zeta, -zeta; -z, z^2], [dom.z, dom.zeta], "name", "q");
335
336
337
338
339
340
qSymbolicCtransp = qSymbolic';
tc.verifyEqual(qSymbolic(1,1).on(), qSymbolicCtransp(1,1).on());
tc.verifyEqual(qSymbolic(2,2).on(), qSymbolicCtransp(2,2).on());
tc.verifyEqual(qSymbolic(1,2).on(), qSymbolicCtransp(2,1).on());
tc.verifyEqual(qSymbolic(2,1).on(), qSymbolicCtransp(1,2).on());

341
qSymbolic2 = quantity.Symbolic(sym("z") * 2i + sym("zeta"), [dom.z, dom.zeta], "name", "im");
342
343
344
345
346
qSymolic2Ctransp = qSymbolic2';
tc.verifyEqual(qSymbolic2.real.on(), qSymolic2Ctransp.real.on());
tc.verifyEqual(qSymbolic2.imag.on(), -qSymolic2Ctransp.imag.on());
end % testCtranspose

347
348
function testTranspose(tc)
syms z zeta
349
350
dom(1) = quantity.Domain("z", linspace(0, 1, 21));
dom(2) = quantity.Domain("zeta", linspace(0, 1, 41));
351
qSymbolic = quantity.Symbolic( [1+z*zeta, -zeta; -z, z^2], dom, "name", "q");
352
353
354
355
356
357
358
qSymbolicTransp = qSymbolic.';
tc.verifyEqual(qSymbolic(1,1).on(), qSymbolicTransp(1,1).on());
tc.verifyEqual(qSymbolic(2,2).on(), qSymbolicTransp(2,2).on());
tc.verifyEqual(qSymbolic(1,2).on(), qSymbolicTransp(2,1).on());
tc.verifyEqual(qSymbolic(2,1).on(), qSymbolicTransp(1,2).on());
end % testTranspose

359
function testFlipDomain(tc)
360
361
syms z zeta
myGrid = linspace(0, 1, 11);
362
363
364
365
366

myDomain(1) = quantity.Domain("z", myGrid);
myDomain(2) = quantity.Domain("zeta", myGrid);

f = quantity.Symbolic([1+z+zeta; 2*zeta+sin(z)] + zeros(2, 1)*z*zeta, myDomain);
367
% flip zeta
368
fReference = quantity.Symbolic([1+z+(1-zeta); 2*(1-zeta)+sin(z)] + zeros(2, 1)*z*zeta, myDomain);
369
fFlipped = f.flipDomain("zeta");
370
tc.verifyEqual(fReference.on, fFlipped.on, "AbsTol", 10*eps);
371
372
373

% flip z and zeta
fReference2 = quantity.Symbolic([1+(1-z)+(1-zeta); 2*(1-zeta)+sin(1-z)] + zeros(2, 1)*z*zeta, ...
374
	myDomain);
375
fFlipped2 = f.flipDomain(["z", "zeta"]);
376
tc.verifyEqual(fReference2.on, fFlipped2.on, "AbsTol", 10*eps);
377
378


379
380
end % testFlipGrid();

381
function testInt(testCase)
382
383
t = sym("t");
T = quantity.Domain.defaultDomain(10, "t");
384

385
f = quantity.Symbolic(sin(t), T);
386
387
F = int(f);

388
testCase.verifyEqual(F, double( int( sin(t), t, 0, 1 ) ), "AbsTol", 6e-17);
389
390
391
392


end

Jakob Gabriel's avatar
Jakob Gabriel committed
393
394
function testCumInt(testCase)
tGrid = linspace(pi, 1.1*pi, 51)';
395
396
s = sym("s");
t = sym("t");
Jakob Gabriel's avatar
Jakob Gabriel committed
397
398
399
400
401
402

a = [ 1, s; t, 1];
b = [ s; 2*s];

%% int_0_t a(t,s) * b(s) ds
% compute symbolic version of the volterra integral
403
myDomain(1) = quantity.Domain("t", tGrid);
Ferdinand Fischer's avatar
Ferdinand Fischer committed
404
myDomain(2) = quantity.Domain("s", tGrid);
405
integrandSymbolic = quantity.Symbolic(a*b, myDomain);
Jakob Gabriel's avatar
Jakob Gabriel committed
406
integrandDiscrete = quantity.Discrete(integrandSymbolic);
407
V = cumInt(integrandSymbolic, "s", tGrid(1), "t");
Ferdinand Fischer's avatar
Ferdinand Fischer committed
408
f = cumInt(integrandDiscrete, "s", tGrid(1), "t");
Jakob Gabriel's avatar
Jakob Gabriel committed
409

410
testCase.verifyEqual(V.on(), f.on(), "AbsTol", 1e-5);
Jakob Gabriel's avatar
Jakob Gabriel committed
411
412

%% int_s_t a(t,s) * b(s) ds
413
414
V = cumInt(integrandSymbolic, "s", "s", "t");
f = cumInt(integrandDiscrete, "s", "s", "t");
Jakob Gabriel's avatar
Jakob Gabriel committed
415

416
testCase.verifyEqual( f.on(f(1).domain), V.on(f(1).domain), "AbsTol", 1e-5 );
Jakob Gabriel's avatar
Jakob Gabriel committed
417
418

%% testCumIntWithLowerAndUpperBoundSpecified
419
420
421
fCumInt2Bcs = cumInt(integrandSymbolic, "s", "zeta", "t");
fCumInt2Cum = cumInt(integrandSymbolic, "s", tGrid(1), "t") ...
	- cumInt(integrandSymbolic, "s", tGrid(1), "zeta");
422

423
424
testCase.verifyEqual(fCumInt2Bcs.on(fCumInt2Bcs(1).domain), ...
	fCumInt2Cum.on(fCumInt2Bcs(1).domain), "AbsTol", 100*eps);
425
426

%% test not an integrable function
427
f = quantity.Symbolic([ sin( sin( s ) ); cos( cos( t ) )], myDomain);
428

429
430
F_sym = cumInt(f, "s", 0, 1);
F_dic = cumInt( quantity.Discrete(f), "s", 0, 1);
431
432
433

testCase.verifyEqual( F_sym.on(), F_dic.on())

Jakob Gabriel's avatar
Jakob Gabriel committed
434
435
end

436
437
438
function testScalarPlusMinusQuantity(testCase)
syms z
myGrid = linspace(0, 1, 7);
439
f = quantity.Symbolic([1; 2] + zeros(2, 1)*z, quantity.Domain("z", myGrid));
440
441
testCase.verifyError(@() 1-f-1, "");
testCase.verifyError(@() 1+f+1, "");
442
443
end

444
445
446
function testNumericVectorPlusMinusQuantity(testCase)
syms z
myGrid = linspace(0, 1, 7);
447
448
f = quantity.Symbolic([1+z; 2+sin(z)] + zeros(2, 1)*z, quantity.Domain("z", myGrid));

449
a = ones(size(f));
450
451
452
453
testCase.verifyEqual(on(a-f), 1-f.on(), "RelTol", 10*eps);
testCase.verifyEqual(on(f-a), f.on()-1, "RelTol", 10*eps);
testCase.verifyEqual(on(a+f), f.on()+1, "RelTol", 10*eps);
testCase.verifyEqual(on(f+a), f.on()+1, "RelTol", 10*eps);
454
455
end

Jakob Gabriel's avatar
Jakob Gabriel committed
456
457
458
function testDiffWith2Variables(testCase)
syms z zeta
myGrid = linspace(0, 1, 7);
459
460
461
462
463
464
myDomain(1) = quantity.Domain("z", myGrid);
myDomain(2) = quantity.Domain("zeta", myGrid);

f = quantity.Symbolic([z*zeta, z; zeta, 1], myDomain);
fdz = quantity.Symbolic([zeta, 1; 0, 0], myDomain);
fdzeta = quantity.Symbolic([z, 0; 1, 0], myDomain);
Jakob Gabriel's avatar
Jakob Gabriel committed
465

466
467
testCase.verifyEqual(on(diff(f, "z")), fdz.on());
testCase.verifyEqual(on(diff(f, "zeta")), fdzeta.on());
Jakob Gabriel's avatar
Jakob Gabriel committed
468
469
470
471
472
473

end

function testGridName2variable(testCase)
syms z zeta eta e a
myGrid = linspace(0, 1, 7);
474
475

myDomain(1) = quantity.Domain("z", myGrid);
476
477
478
479
myDomain(2) = myDomain(1).rename("zeta");
myDomain(3) = myDomain(1).rename("eta");
myDomain(4) = myDomain(1).rename("e");
myDomain(5) = myDomain(1).rename("a");
480
481

obj = quantity.Symbolic([z*zeta, eta*z; e*a, a], myDomain);
Jakob Gabriel's avatar
Jakob Gabriel committed
482
483
484

thisGridName = {'eta', 'e', 'z'};
thisVariable = gridName2variable(obj, thisGridName);
485
testCase.verifyEqual(string(thisGridName).', string(thisVariable));
Jakob Gabriel's avatar
Jakob Gabriel committed
486
487
488
489
490

thisGridName = 'z';
thisVariable = gridName2variable(obj, thisGridName);
testCase.verifyEqual(thisGridName, char(thisVariable));

491
492
testCase.verifyError(@() gridName2variable(obj, 't'), "");
testCase.verifyError(@() gridName2variable(obj, 123), "");
Jakob Gabriel's avatar
Jakob Gabriel committed
493
494
495

end

496
497
function testCat(testCase)
syms z zeta
498
myDomain = quantity.Domain.defaultDomain(21, "z");
499
500
f1 = quantity.Symbolic(1+z*z, myDomain, "name", "f1");
f2 = quantity.Symbolic(sin(z), myDomain, "name", "f2");
501

502
% vertical concatenation
503
504
505
506
F = [f1; f2];
testCase.verifyTrue(F(1) == f1)
testCase.verifyTrue(F(2) == f2)

507
% horizontal concatenation
508
509
510
511
F = [f1, f2];
testCase.verifyTrue(F(1) == f1)
testCase.verifyTrue(F(2) == f2)

512
% combined concatenation
513
514
515
516
517
518
F = [F; f1, f2];
testCase.verifyTrue(F(1,1) == f1);
testCase.verifyTrue(F(1,2) == f2);
testCase.verifyTrue(F(2,1) == f1);
testCase.verifyTrue(F(2,2) == f2);

519
% concatenation on different grids
520
myDomain = quantity.Domain.defaultDomain(13, "z");
521
f3 = quantity.Symbolic(cos(z), myDomain, "name", "f1");
522
523

F = [f1, f3];
524
testCase.verifyEqual(F(2).domain.gridLength, f1(1).domain.gridLength)
525

526
527
528
529
530
531
% concatenate quantitySymbolics and other values
F = [f1, 0];
testCase.verifyEqual(F(2).valueSymbolic, sym(0))
F = [f1, sym(zeros(1, 10))];
testCase.verifyEqual([F(2:end).valueSymbolic], sym(zeros(1,10)))

532
533
534
% 
F = [0, f1];
testCase.verifyEqual(F(1).valueSymbolic, sym(0))
535
end
536
537
538
539

function testExp(testCase)
% 1 spatial variable
syms z zeta
540
myDomain = quantity.Domain.defaultDomain(21, "z");
541
s1d = quantity.Symbolic(1+z*z, myDomain, "name", "s1d");
542
543
544
testCase.verifyEqual(s1d.exp.on(), exp(s1d.on()));

% diagonal matrix
545
546
547
548
myDomain(2) = quantity.Domain.defaultDomain(41, "zeta");



549
s2dDiag = quantity.Symbolic([1+z*zeta, 0; 0, z^2], myDomain, "name", "s2dDiag");
550
551
552
553
554
testCase.verifyEqual(s2dDiag.exp.on(), exp(s2dDiag.on()));
end

function testExpm(testCase)
syms z zeta
555
556
557
558

myDomain = quantity.Domain.defaultDomain(21, "z");
myDomain(2) = quantity.Domain.defaultDomain(41, "zeta");

559
mat2d = quantity.Symbolic(ones(2, 2) + [1+z*zeta, 3*zeta; 2+5*z+zeta, z^2], myDomain, "name", "s2d");
560
561
562
563
564
565
566
mat2dMat = mat2d.on();
mat2dExpm = 0 * mat2d.on();
for zIdx = 1 : 21
	for zetaIdx = 1 : 41
		mat2dExpm(zIdx, zetaIdx, :, :) = expm(reshape(mat2dMat(zIdx, zetaIdx, :, :), [2, 2]));
	end
end
567
testCase.verifyEqual(mat2d.expm.on(), mat2dExpm, "RelTol", 100*eps);
568
569
570
571
572
573
574
575

end


function testOn(testCase)
syms z zeta
gridZ = linspace(0, 1, 40);
gridZeta = linspace(0, 1, 21);
576
577
578
579
580

myDomain = quantity.Domain("z", gridZ);
myDomain(2) = quantity.Domain("zeta", gridZeta);

A = quantity.Symbolic(2+[z, z^2, z+zeta; -sin(zeta), z*zeta, 1], myDomain, 'name', 'A');
581
582
583
584
585
586
587
588
589
590

%%
AOnPure = A.on();
testCase.verifyEqual(2+gridZ(:)*ones(1, numel(gridZeta)), AOnPure(:,:,1,1));
testCase.verifyEqual(2+gridZ(:)*gridZeta, AOnPure(:,:,2,2));
testCase.verifyEqual(3+0*gridZ(:)*gridZeta, AOnPure(:,:,2,3));

%%
gridZ2 = linspace(0, 1, 15);
gridZeta2 = linspace(0, 1, 31);
591
AOnChangedGrid = A.on({gridZ2, gridZeta2}, ["z", "zeta"]);
592
593
594
595
596
597
598
599
600
601
602
603
testCase.verifyEqual(2+gridZ2(:)*ones(1, numel(gridZeta2)), AOnChangedGrid(:,:,1,1));
testCase.verifyEqual(2+gridZ2(:)*gridZeta2, AOnChangedGrid(:,:,2,2));
testCase.verifyEqual(3+0*gridZ2(:)*gridZeta2, AOnChangedGrid(:,:,2,3));

%%
AOnChangedGrid2 = A.on({gridZeta2, gridZ2}, {'zeta', 'z'});
testCase.verifyEqual((2+gridZ2(:)*ones(1, numel(gridZeta2))).', AOnChangedGrid2(:,:,1,1));
testCase.verifyEqual((2+gridZ2(:)*gridZeta2).', AOnChangedGrid2(:,:,2,2));
testCase.verifyEqual((3+0*gridZ2(:)*gridZeta2).', AOnChangedGrid2(:,:,2,3));
% testCase.verifyEqual(ones([21, 31, 2, 3]), ...
% 	A.on({linspace(0, 1, 21), linspace(0, 1, 31)}));
% testCase.verifyEqual(ones([21, 31, 2, 3]), ...
604
% 	A.on({linspace(0, 1, 21), linspace(0, 1, 31)}, ["z", "zeta"]));
605
606
607
608
609
610
611
% testCase.verifyEqual(ones([21, 31, 2, 3]), ...
% 	A.on({linspace(0, 1, 21), linspace(0, 1, 31)}, {'zeta', 'z'}));
end

function testSqrt(testCase)
% 1 spatial variable
syms z zeta
612
myDomain = quantity.Domain("z", linspace(0,1,21));
613
s1d = quantity.Symbolic(1+z*z, myDomain, "name", "s1d");
614
615
616
testCase.verifyEqual(s1d.sqrt.on(), sqrt(s1d.on()));

% 2 spatial variables
617
myDomain(2) = quantity.Domain("zeta", linspace(0,1,41));
618
s2d = quantity.Symbolic(1+z*zeta, myDomain, "name", "s2d");
619
620
621
testCase.verifyEqual(s2d.sqrt.on(), sqrt(s2d.on()));

% diagonal matrix
622
s2dDiag = quantity.Symbolic([1+z*zeta, 0; 0, z^2], myDomain, "name", "s2dDiag");
623
624
625
626
627
628
testCase.verifyEqual(s2dDiag.sqrt.on(), sqrt(s2dDiag.on()));

end

function testSqrtm(testCase)
syms z zeta
629
630
631
632

myDomain = quantity.Domain("z", linspace(0,1,21));
myDomain(2) = quantity.Domain("zeta", linspace(0,1,41));

633
mat2d = quantity.Symbolic(ones(2, 2) + [1+z*zeta, 3*zeta; 2+5*z+zeta, z^2], myDomain, "name", "s2d");
634
635
636
637
638
639
640
mat2dMat = mat2d.on();
mat2dSqrtm = 0 * mat2d.on();
for zIdx = 1 : 21
	for zetaIdx = 1 : 41
		mat2dSqrtm(zIdx, zetaIdx, :, :) = sqrtm(reshape(mat2dMat(zIdx, zetaIdx, :, :), [2, 2]));
	end
end
641
testCase.verifyEqual(mat2d.sqrtm.on(), mat2dSqrtm, "AbsTol", 10*eps);
642
mat2drootquad = sqrtm(mat2d)*sqrtm(mat2d);
643
testCase.verifyEqual(mat2drootquad.on(), mat2d.on(), "AbsTol", 100*eps);
644
645
646
647
648
649
650

end

function testInv(testCase)
%% scalar example
syms z
zGrid = linspace(0, 1, 21);
651
652
653

myDomain = quantity.Domain("z", zGrid);

654
v = quantity.Symbolic(1+z*z, myDomain, "name", "s1d");
655
656
657
658
659
660
661
662
vInvReference = 1 ./ (1 + zGrid.^2);
vInv = v.inv();

%% matrix example
syms zeta
zetaGrid = linspace(0, pi);
[zNdGrid, zetaNdGrid] = ndgrid(zetaGrid, zGrid);

663
664
665
myDomain(2) = quantity.Domain("zeta", zetaGrid);

K = quantity.Symbolic([1+z*z, 0; 0, 1+z*zeta], myDomain);
666
667
668
669
670
671
672
673
674
675
676
kInv = inv(K);

%%
testCase.verifyEqual(vInv.on(), vInvReference(:))
testCase.verifyEqual(kInv(1, 1).on(), repmat(vInvReference(:), [1, 100]))
testCase.verifyEqual(kInv(2, 2).on(), 1./(1+zNdGrid .* zetaNdGrid).')

end

function testZero(testCase)

677
x = quantity.Symbolic(0, quantity.Domain.empty());
678
679

verifyEqual(testCase, x.on(), 0);
680
verifyTrue(testCase, all(on(x) * magic(5) == zeros(5), "all"));
681
682
683
684
685
686

end

function testCast2Quantity(testCase)

syms z t
687
688
d = [quantity.Domain("z", linspace(0, 1, 5)), ...
	 quantity.Domain("t", linspace(0, 2, 2))];
689

690
x = quantity.Symbolic(sin(z * t * pi), d);
691
692
693
694
695
696
697
698

X = quantity.Discrete(x);

verifyTrue(testCase, numeric.near(X.on(), x.on()));
end

function testMrdivide(testCase)
syms z t
699
700
d = [quantity.Domain("z", linspace(0, 1, 5)), ...
	 quantity.Domain("t", linspace(0, 2, 2))];
701
x = quantity.Symbolic(sin(z * t * pi), d);
702
703
o = x / x;

704
verifyTrue(testCase, all(o.on() == 1, "all"));
705
706
707

% test constant values
invsin = 1 / ( x + 8);
708
709
ndGridNT = x(1).domain.ndgrid;
verifyTrue(testCase, numeric.near(1 ./ (sin(ndGridNT{1} .* ndGridNT{2} * pi) + 8), invsin.on()));
710
711
712
713
714

end

function testMTimesConstants(testCase)
syms z t
715
716
d = [quantity.Domain("z", linspace(0, 1, 5)), ...
	 quantity.Domain("t", linspace(0, 2, 2))];
717
x = quantity.Symbolic(sin(z * t * pi), d);
718
719
720

% test multiplicaiton iwth a constant scalar
x5 = x * 5;
721
722
ndGridNT = x(1).domain.ndgrid;
X = sin(ndGridNT{1} .* ndGridNT{2} * pi);
723
724
725
726
727
728
729
730
731
732
733
734
735

verifyTrue(testCase, numeric.near(X*5, x5.on()));

% test multiplication with a constant vector
v = [2, 3];
x11 = x * v;
X11On = x11.on();

verifyTrue(testCase, numeric.near(X*v(1), X11On(:,:,:,1)));
verifyTrue(testCase, numeric.near(X*v(2), X11On(:,:,:,2)));

% test multiplication with zero:
x0 = x*0;
736
testCase.verifyEqual(x0.on(), 0*ndGridNT{1});
737
738
739
740
741
742
end

function testMLRdivide2ConstantSymbolic(testCase)

syms z
assume(z>0 & z<1);
743
myDomain = quantity.Domain.defaultDomain(100, "z");
744
745
Lambda = quantity.Symbolic(eye(2), myDomain, "name", "Lambda");
A = quantity.Symbolic(eye(2), myDomain, "name", "A");
746
747
748
749
750
751
752
753
754
C1 = A / Lambda;
C2 = A \ Lambda;
testCase.verifyEqual(C1.on(), A.on());
testCase.verifyEqual(C2.on(), A.on());

end
function testGridNamesNotAllowed(testCase)

try
755
756
757
	myDomain(1) = quantity.Domain.defaultDomain(7, "y");
	myDomain(2) = quantity.Domain.defaultDomain(5, "t");
	myDomain(3) = quantity.Domain.defaultDomain(3, "asdf");
758
	quantity.Symbolic(sym("z"), myDomain);
759
760
761
762
763
	
	% the initialization above has to throw an exception. If the code is
	% here, the test is not valid.
	verifyTrue(testCase, false);
catch ex
764
	verifyTrue(testCase, isa(ex, "MException"));
765
766
767
768
end

end

769
function testCastSymbolic2Discrete(testCase)
770
771

syms z t
772
x = quantity.Symbolic(sin(z * t * pi), ...
773
	[quantity.Domain.defaultDomain(100, "z"), quantity.Domain("t", linspace(0, 2, 200))]);
774
775
776
777
778
779
780
781
782
783

X1 = quantity.Discrete(x);

%%
verifyEqual(testCase, X1.on, x.on);

end

function testDiff(testCase)
%%
784
z = sym("z", "real");
785
786
787
f1 = sinh(z * pi);
f2 = cosh(z * pi);

788
f = quantity.Symbolic([f1 ; f2 ],  quantity.Domain("z", linspace(0,1,11)), "name", "f");
789

790
791
d2f = f.diff("z", 2);
d1f = f.diff("z", 1);
792
793
794
795
796
797

%%
F1 = symfun([f1; f2], z);
D1F = symfun([diff(f1,1); diff(f2,1)], z);
D2F = symfun([diff(f1,2); diff(f2,2)], z);

798
z = f(1).domain.grid;
799
800
801
802
803
804
805

R1 = F1(z);
R2 = D1F(z);
R3 = D2F(z);

%%
verifyTrue(testCase, numeric.near(double([R1{:}]), f.on(), 1e-12));
806
807
verifyTrue(testCase, numeric.near(double([R2{:}]), f.diff("z", 1).on(), 1e-12));
verifyTrue(testCase, numeric.near(double([R3{:}]), f.diff("z", 2).on(), 1e-12));
808
809
810
811
812

end

function testInit(testCase)

813
z = linspace(0,1).';
814
t = linspace(0,1,3);
815
816
817
818

syms x y

v = [sin(x * y * pi); cos(x * y * pi)];
819

820
try
821
822
823
824
825
	
	myDomain(1) = quantity.Domain("z", z);
	myDomain(2) = quantity.Domain("t", t);
	
	b = quantity.Symbolic(v, myDomain);
826
827
828
	testCase.verifTrue(false);
end
	
829
830
831
832
833
834

myDom2(1) = quantity.Domain("x", z);
myDom2(2) = quantity.Domain("y", t);

b = quantity.Symbolic(v, myDom2);
c = quantity.Symbolic(v, myDom2);
835

836
%%
837
verifyTrue(testCase, all(b.on() == c.on(), "all"));
838
839
840
end

function testPlus(testCase)
Jakob Gabriel's avatar
Jakob Gabriel committed
841
842
syms z zeta

843
844
Z = quantity.Domain("z", linspace(0, 1, 3));
Zeta = quantity.Domain("zeta", linspace(0, 1, 4));
845

Jakob Gabriel's avatar
Jakob Gabriel committed
846
%% 1 + 2 variables
847
848
fSymVec3 = quantity.Symbolic([sinh(pi * z)], Z, "name", "sinh");
fSymVec4 = quantity.Symbolic([cosh(zeta * pi * z)], [Z Zeta], "name", "cosh");
Jakob Gabriel's avatar
Jakob Gabriel committed
849
850
851
852

result = fSymVec3 - fSymVec4;
resultDiscrete = quantity.Discrete(fSymVec3) - quantity.Discrete(fSymVec4);

853
%%
854
testCase.verifyEqual(result.on(), resultDiscrete.on(), "AbsTol", 10*eps);
855

Jakob Gabriel's avatar
Jakob Gabriel committed
856
%% quantity + constant
857
858

myQuan = quantity.Symbolic([sinh(pi * z), cosh(zeta * pi * z); 1, z^2], ...
859
	[Z Zeta], "name", "myQuan");
Jakob Gabriel's avatar
Jakob Gabriel committed
860
861
myQuanPlus1 = myQuan + ones(size(myQuan));
my1PlusQuan =  ones(size(myQuan)) + myQuan;
862

Jakob Gabriel's avatar
Jakob Gabriel committed
863
%%
864
865
testCase.verifyEqual(myQuan.on()+1, myQuanPlus1.on(), "AbsTol", 10*eps);
testCase.verifyEqual(myQuan.on()+1, my1PlusQuan.on(), "AbsTol", 10*eps);
Jakob Gabriel's avatar
Jakob Gabriel committed
866
867
868


%% 1 variable
869
fSymVec = quantity.Symbolic([sinh(z * pi) ; cosh(z * pi)], Z, "name", "sinhcosh");
Jakob Gabriel's avatar
Jakob Gabriel committed
870
871

F = fSymVec + fSymVec;
872
873
val = F.on();

Jakob Gabriel's avatar
Jakob Gabriel committed
874
875
f1fun = @(z) sinh(z * pi);
f2fun = @(z) cosh(z * pi) + 0*z;
876

877
zGrid = fSymVec(1).domain.grid;
878
879

%%
Jakob Gabriel's avatar
Jakob Gabriel committed
880
881
882
883
884
885
886
887
verifyTrue(testCase, numeric.near(2 * f1fun(zGrid) , val(:, 1), 1e-12));
verifyTrue(testCase, numeric.near(2 * f2fun(zGrid) , val(:, 2), 1e-12));

%% 2 variable
syms zeta
f1fun = @(z, zeta) sinh(z * pi);
f2fun = @(z, zeta) cosh(zeta * pi) + 0*z;

888
fSymVec = quantity.Symbolic([sinh(z * pi); cosh(zeta * pi)], [Z Zeta], "name", "sinhcosh");
Jakob Gabriel's avatar
Jakob Gabriel committed
889
890
891
892
fSymVec2 = [fSymVec(2); fSymVec(1)];

F = fSymVec + fSymVec;
F1Minus2 = fSymVec - fSymVec2;
893

894
895
zZetaNdgrid = fSymVec(1).domain.ndgrid();
zGrid = zZetaNdgrid{1}; zetaGrid = zZetaNdgrid{2};
896
%%
Jakob Gabriel's avatar
Jakob Gabriel committed
897
898
899
900
901
902
verifyTrue(testCase, numeric.near(2 * f1fun(zGrid, zetaGrid) , F(1).on(), 1e-12));
verifyTrue(testCase, numeric.near(2 * f2fun(zGrid, zetaGrid) , F(2).on(), 1e-12));
verifyTrue(testCase, numeric.near(f1fun(zGrid, zetaGrid) - f2fun(zGrid, zetaGrid), ...
	F1Minus2(1).on(), 1e-12));
verifyTrue(testCase, numeric.near(-F1Minus2(2).on(), ...
	F1Minus2(1).on(), 1e-12));
903
904
905
906
907
908
909
910
end

function testMTimes(testCase)
%%
syms z
f1 = sinh(z * pi);
f2 = cosh(z * pi);

911
f = quantity.Symbolic([f1 ; f2 ], quantity.Domain.defaultDomain(3, "z"), "name", "sinhcosh");
912
913
914
915
916
917
918
919

F = f * f.';
val = F.on();

f1 = matlabFunction(f1);
f2 = matlabFunction(f2);


920
z = f(1).domain.grid;
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936

%%
verifyTrue(testCase, numeric.near(f1(z).^2, val(:, 1, 1), 1e-12));
verifyTrue(testCase, numeric.near(f2(z).^2, val(:, 2, 2), 1e-12));
verifyTrue(testCase, numeric.near(f1(z) .* f2(z) , val(:, 1, 2), 1e-12));


%% test multiplication with a scalar
ff = f * 2;
verifyTrue(testCase, numeric.near(ff.on(), [f1(z), f2(z)] * 2, 1e-12));


end

function testVariable(testCase)
syms z
937
q = quantity.Symbolic(z, quantity.Domain('z', 1));
938
939
940
941
942

verifyEqual(testCase, q.variable, z);

end

Ferdinand Fischer's avatar
tmp    
Ferdinand Fischer committed
943
function testUMinusAndUPlus(testCase)
944
945
946
947
948
%%
syms z;
f1 = (sin(z .* pi) );
f2 =  z;
f3 = (cos(z .* pi) );
949
d = quantity.Domain("z", linspace(0, 1, 3).');
Ferdinand Fischer's avatar
tmp    
Ferdinand Fischer committed
950

951
f = quantity.Symbolic([f1, f2; 0, f3], d);
Ferdinand Fischer's avatar
tmp    
Ferdinand Fischer committed
952
953
fDisc = quantity.Discrete(f);
values = f.on();
954
955

%%
Ferdinand Fischer's avatar
tmp    
Ferdinand Fischer committed
956
testCase.verifyEqual(on( -f ), -values);
957
testCase.verifyEqual(on( + ( - f ) ), -values );
Ferdinand Fischer's avatar
tmp    
Ferdinand Fischer committed
958
959
960
961
testCase.verifyEqual(on( + fDisc), on(+f));
testCase.verifyEqual(on( - fDisc), on(-f));


962
963
964
end

function testConstantValues(testCase)
965
%% evaluation of a constant on the same grid
966
syms z
967
968
myDomain = quantity.Domain.defaultDomain(4, "z");
q = quantity.Symbolic(magic(3), myDomain);
969

970
971
magic7 = magic(3);
magic700 = zeros(4, 3, 3);
972
for k = 1:q(1).domain.gridLength()
973
974
975
	magic700(k, :,:) = magic7;
end

Ferdinand Fischer's avatar
Ferdinand Fischer committed
976
977
testCase.verifyEqual(magic700, q.on());

978
979
980
981
982
983
984
985
986
987
988
989
%% evaluation of a constant on a different grid
myDomain2 = quantity.Domain.defaultDomain(7, "z");
magic7 = magic(3);
magic700 = zeros(myDomain2.n, 3, 3);
for k = 1:myDomain2.n
	magic700(k, :,:) = magic7;
end
q.on(myDomain2);

testCase.verifyEqual(magic700, q.on(myDomain2));

%%
990
p = quantity.Symbolic([0 0 0 0 z], myDomain);
991
testCase.verifyFalse(p.isNumber());
992

993
994
995
996
997
998
999
1000

%%
myDomains = [myDomain, myDomain.rename("t")];
Q = quantity.Symbolic(sym(magic(2)), myDomains);
Qon = Q.on();
magic2 = magic(2);

for k = 1:numel(magic2)