1
2
3
4
5
6
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
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
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
|
/*
Simple testsuite for diag.mac
Run this file with
batch ("diag_test", 'test)
*/
kill (all);
done$
(load ("diag"), 0);
0$
diag ([1,2,3]);
matrix ([1,0,0], [0,2,0], [0,0,3])$
diag ([matrix ([1]), matrix ([2]), matrix ([3])]);
matrix ([1,0,0], [0,2,0], [0,0,3])$
(a1: matrix ([1,2,3], [0,4,5], [0,0,6]),
a2: matrix([1,1], [1,0]),
0);
0$
diag ([a1, x, a2]);
matrix ([1,2,3,0,0,0],
[0,4,5,0,0,0],
[0,0,6,0,0,0],
[0,0,0,x,0,0],
[0,0,0,0,1,1],
[0,0,0,0,1,0])$
kill (a1, a2);
done$
JF (0, 1);
matrix ([0])$
JF (x, 3);
matrix ([x, 1, 0],
[0, x, 1],
[0, 0, x])$
dispJordan ([[0, 1], [2, 2], [3, 2, 1]]);
''(diag ([JF (0, 1), JF (2, 2), JF (3, 2), JF (3, 1)]))$
/*
For testing jordan(), use JCF examples built by dispJordan and
conjugate by invertible matrices.
*/
(make_conjugate (j_info) := block (
[n: lsum (x, x, map (lambda ([lst], lsum (x, x, rest (lst))), j_info)),
T],
T: genmatrix (lambda([i,j], if is (i = n+1-j) then 0 else 1), n, n),
T^^(-1) . dispJordan (j_info) . T),
test_jordan (j_info) := is (jordan (make_conjugate (j_info)) = j_info),
test_modematrix (j_info) := block (
[conj: make_conjugate (j_info), T],
T: ModeMatrix (conj, j_info),
is (T^^(-1) . conj . T = dispJordan (j_info))),
0);
0$
jordan (matrix ([1]));
[[1, 1]]$
test_jordan ([[%i, 2]]);
true$
test_jordan ([[1, 1], [2, 2, 1]]);
true$
minimalPoly ([[%i, 1]]);
x - %i$
minimalPoly ([[1, 2, 1]]);
(x - 1)^2$
minimalPoly ([[1, 1], [2, 1]]);
(x - 1)*(x - 2)$
diag_sorted_list_histogram ([]);
[]$
diag_sorted_list_histogram ([1]);
[[1, 1]]$
diag_sorted_list_histogram ([3, 3, 2, 2, 1]);
[[3, 2], [2, 2], [1, 1]]$
/*
Testing diag_kernel_element, diag_general_jordan_chain and
diag_find_li_chain is probably a bit difficult. Let's punt on that
and hope they're well exercised by ModeMatrix.
*/
test_modematrix ([[%i, 2]]);
true$
test_modematrix ([[1, 1], [2, 2, 1]]);
true$
diag_taylor_coefficients (exp(x), x, 3);
[exp(x), exp(x), exp(x)/2, exp(x)/6]$
diag_taylor_expand_block (diag_taylor_coefficients (exp(x), x, 3),
x, 1, 3);
matrix ([exp(1), exp(1), exp(1)/2],
[0, exp(1), exp(1)],
[0, 0, exp(1)])$
mat_function (sin, ident(3)*t);
''(ident(3)*sin(t))$
/*
Bug 2542
The "correct" answer was calculated numerically by Octave. There's
probably a better way of checking this!
*/
block ([a: matrix([6,6,-3,2],[-4,-4,2,0],[8,7,-4,4],[1,0,-1,-2]),
num_ans:
matrix ([ 2.9848e+00, 2.3771e+00, -1.2382e+00, 6.1912e-01],
[-1.2382e+00, -7.2991e-01, 6.1912e-01, 1.3645e-16],
[ 3.7147e+00, 3.2064e+00, -1.3490e+00, 1.2382e+00],
[-1.1079e-01, -2.6557e-01, -9.9383e-02, -1.1079e-01]),
err],
err: abs (float (rectform (mat_function (exp, a))) - num_ans),
is (lmax (map (lmax, args (err))) < 1e-4));
true$
|
