← Back to branch summary

~ubuntu-branches/ubuntu/karmic/scilab/karmic

« back to all changes in this revision

Viewing changes to tests/calpol.tst

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2002-03-21 16:57:43 UTC
  • Revision ID: james.westby@ubuntu.com-20020321165743-e9mv12c1tb1plztg
Tags: upstream-2.6
Import upstream version 2.6

Show diffs side-by-side

added added

removed removed

Lines of Context:
tests/calpol.tst
 
1
// Copyright INRIA
 
2
//
 
3
deff('[ok]=cmpr(h1,h2,eps)',['h1=h1-h2;';
 
4
         'if norm(coeff(h1(2)))>eps then ok=0,else ok=1,end'])
 
5
 
 
6
s=poly(0,'s');
 
7
//
 
8
//rationals
 
9
//
 
10
num=1;den=1+s;
 
11
eps=5000*%eps;
 
12
rtyp=['r','num','den','dt'];
 
13
if cmpr(num/den,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
14
if cmpr(den\num,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
15
if cmpr(num./den,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
16
if cmpr(den.\num,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
17
num=1.5+s**3;
 
18
if cmpr(num/den,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
19
if cmpr(den\num,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
20
if cmpr(num./den,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
21
if cmpr(den.\num,tlist(rtyp,num,den,[]),eps)<>1 then pause,end
 
22
//
 
23
h1=num/den;x=1.5;
 
24
if cmpr(h1+x,tlist(rtyp,num+x*den,den,[]),eps)<>1 then pause,end
 
25
if cmpr(x+h1,tlist(rtyp,num+x*den,den,[]),eps)<>1 then pause,end
 
26
if cmpr(h1-x,tlist(rtyp,num-x*den,den,[]),eps)<>1 then pause,end
 
27
if cmpr(x-h1,tlist(rtyp,-num+x*den,den,[]),eps)<>1 then pause,end
 
28
 
 
29
x=1.5+3*s;
 
30
if cmpr(h1+x,tlist(rtyp,num+x*den,den,[]),eps)<>1 then pause,end
 
31
if cmpr(x+h1,tlist(rtyp,num+x*den,den,[]),eps)<>1 then pause,end
 
32
if cmpr(h1-x,tlist(rtyp,num-x*den,den,[]),eps)<>1 then pause,end
 
33
if cmpr(x-h1,tlist(rtyp,-num+x*den,den,[]),eps)<>1 then pause,end
 
34
y=s**3;h2=x/y;
 
35
if cmpr(h1+h2,tlist(rtyp,num*y+x*den,den*y,[]),eps)<>1 then pause,end
 
36
if cmpr(h1-h2,tlist(rtyp,num*y-x*den,den*y,[]),eps)<>1 then pause,end
 
37
//
 
38
// concatenations
 
39
//
 
40
h1=num/den;x=1.5;
 
41
if cmpr([h1,x],tlist(rtyp,[num,x],[den,1],[]),eps)<>1 then pause,end
 
42
if cmpr([x,h1],tlist(rtyp,[x,num],[1,den],[]),eps)<>1 then pause,end
 
43
if cmpr([h1;x],tlist(rtyp,[num;x],[den;1],[]),eps)<>1 then pause,end
 
44
if cmpr([x;h1],tlist(rtyp,[x;num],[1;den],[]),eps)<>1 then pause,end
 
45
x=1.5+3*s;
 
46
if cmpr([h1,x],tlist(rtyp,[num,x],[den,1],[]),eps)<>1 then pause,end
 
47
if cmpr([x,h1],tlist(rtyp,[x,num],[1,den],[]),eps)<>1 then pause,end
 
48
if cmpr([h1;x],tlist(rtyp,[num;x],[den;1],[]),eps)<>1 then pause,end
 
49
if cmpr([x;h1],tlist(rtyp,[x;num],[1;den],[]),eps)<>1 then pause,end
 
50
y=-0.5+s**3;h2=x/y;
 
51
if cmpr([h1,h2],tlist(rtyp,[num,x],[den,y],[]),eps)<>1 then pause,end
 
52
if cmpr([h1;h2],tlist(rtyp,[num;x],[den;y],[]),eps)<>1 then pause,end
 
53
h1=[num/den,den/num];x=[0.3 1.5];
 
54
if cmpr([h1,x],tlist(rtyp,[h1(2),x],[h1(3),ones(x)],[]),eps)<>1 then pause,end
 
55
if cmpr([x,h1],tlist(rtyp,[x,h1(2)],[ones(x),h1(3)],[]),eps)<>1 then pause,end
 
56
if cmpr([h1;x],tlist(rtyp,[h1(2);x],[h1(3);ones(x)],[]),eps)<>1 then pause,end
 
57
if cmpr([x;h1],tlist(rtyp,[x;h1(2)],[ones(x);h1(3)],[]),eps)<>1 then pause,end
 
58
h1=[num/den;den/num];x=[0.3;-1.5];
 
59
if cmpr([h1,x],tlist(rtyp,[h1(2),x],[h1(3),ones(x)],[]),eps)<>1 then pause,end
 
60
if cmpr([x,h1],tlist(rtyp,[x,h1(2)],[ones(x),h1(3)],[]),eps)<>1 then pause,end
 
61
if cmpr([h1;x],tlist(rtyp,[h1(2);x],[h1(3);ones(x)],[]),eps)<>1 then pause,end
 
62
if cmpr([x;h1],tlist(rtyp,[x;h1(2)],[ones(x);h1(3)],[]),eps)<>1 then pause,end
 
63
x=[1.5+3*s;-1+s**3];
 
64
if cmpr([h1,x],tlist(rtyp,[h1(2),x],[h1(3),ones(x)],[]),eps)<>1 then pause,end
 
65
if cmpr([x,h1],tlist(rtyp,[x,h1(2)],[ones(x),h1(3)],[]),eps)<>1 then pause,end
 
66
if cmpr([h1;x],tlist(rtyp,[h1(2);x],[h1(3);ones(x)],[]),eps)<>1 then pause,end
 
67
if cmpr([x;h1],tlist(rtyp,[x;h1(2)],[ones(x);h1(3)],[]),eps)<>1 then pause,end
 
68
h1=[num/den,den/num];x=[1.5+3*s,-1+s**2];
 
69
if cmpr([h1,x],tlist(rtyp,[h1(2),x],[h1(3),ones(x)],[]),eps)<>1 then pause,end
 
70
if cmpr([x,h1],tlist(rtyp,[x,h1(2)],[ones(x),h1(3)],[]),eps)<>1 then pause,end
 
71
if cmpr([h1;x],tlist(rtyp,[h1(2);x],[h1(3);ones(x)],[]),eps)<>1 then pause,end
 
72
if cmpr([x;h1],tlist(rtyp,[x;h1(2)],[ones(x);h1(3)],[]),eps)<>1 then pause,end
 
73
h1=[num/den;den/num];y=-0.5+s**3;h2=[num/y;y*y/(y+1)];
 
74
if cmpr([h1,h2],tlist(rtyp,[h1(2),h2(2)],[h1(3),h2(3)],[]),eps)<>1 then pause,end
 
75
if cmpr([h1;h2],tlist(rtyp,[h1(2);h2(2)],[h1(3);h2(3)],[]),eps)<>1 then pause,end
 
76
h1=[num/den,den/num];y=-0.5+s**3;h2=[num/y,y*y/(y+1)];
 
77
if cmpr([h1,h2],tlist(rtyp,[h1(2),h2(2)],[h1(3),h2(3)],[]),eps)<>1 then pause,end
 
78
if cmpr([h1;h2],tlist(rtyp,[h1(2);h2(2)],[h1(3);h2(3)],[]),eps)<>1 then pause,end
 
79
//
 
80
// extraction
 
81
//
 
82
h1=[num/den,den/num];
 
83
if cmpr(h1(1,1),num/den,eps)<>1 then pause,end
 
84
if cmpr(h1(1,1:2),h1,eps)<>1 then pause,end
 
85
if cmpr(h1(1,[2 1]),[den/num,num/den],eps)<>1 then pause,end
 
86
h1=[num/den;den/num];
 
87
if cmpr(h1(2,1),den/num,eps)<>1 then pause,end
 
88
if cmpr(h1(1:2,1),h1,eps)<>1 then pause,end
 
89
if cmpr(h1([2 1],1),[den/num;num/den],eps)<>1 then pause,end
 
90
y=-0.5+s**3;h1=[num/den,den/num];h2=[num/den,den/num;num/y,y*y/(y+1)];
 
91
if cmpr(h2(2,1),num/y,eps)<>1 then pause,end
 
92
if cmpr(h2(1:2,1:2),h2,eps)<>1 then pause,end
 
93
if cmpr(h2([2 1],1),[num/y;num/den],eps)<>1 then pause,end
 
94
if cmpr(h2(:,1),[num/den;num/y],eps)<>1 then pause,end
 
95
if cmpr(h2(2,:),[num/y,y*y/(y+1)],eps)<>1 then pause,end
 
96
if cmpr(h2(:,:),h2,eps)<>1 then pause,end
 
97
//
 
98
// insertions
 
99
//
 
100
h1=[num/den,den/num];x=0.33;
 
101
hh=h1;hh(1,1)=x;if cmpr(hh,[x,h1(1,2)],eps)<>1 then pause,end;
 
102
x=[-2.67 0.8];
 
103
hh=h1;hh(1,1:2)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
104
hh=h1;hh(1,[2 1])=x;if cmpr(hh,x([2,1]),eps)<>1 then pause,end;
 
105
hh=h1;hh(1,[2;1])=x;if cmpr(hh,x([2,1]),eps)<>1 then pause,end;
 
106
hh=h1;hh(1,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
107
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
108
h1=[num/den;den/num];x=0.33;
 
109
hh=h1;hh(1,1)=x;if cmpr(hh,[x;h1(2,1)],eps)<>1 then pause,end;
 
110
x=[-2.67;0.8];
 
111
hh=h1;hh(1:2,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
112
hh=h1;hh(:,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
113
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
114
hh=h1;hh([2 1],1)=x;if cmpr(hh,x([2 1]),eps)<>1 then pause,end;
 
115
hh=h1;hh([2;1],1)=x;if cmpr(hh,x([2 1]),eps)<>1 then pause,end;
 
116
 
 
117
h1=[num/den,den/num];x=0.33*s+1;
 
118
hh=h1;hh(1,1)=x;if cmpr(hh,[x,h1(1,2)],eps)<>1 then pause,end;
 
119
x=[-2.67 0.8+s**3];
 
120
hh=h1;hh(1,1:2)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
121
hh=h1;hh(1,[2 1])=x;if cmpr(hh,x([2,1]),eps)<>1 then pause,end;
 
122
hh=h1;hh(1,[2;1])=x;if cmpr(hh,x([2,1]),eps)<>1 then pause,end;
 
123
hh=h1;hh(1,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
124
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
125
h1=[num/den;den/num];x=-0.33+38*s;
 
126
hh=h1;hh(1,1)=x;if cmpr(hh,[x;h1(2,1)],eps)<>1 then pause,end;
 
127
x=[-2.67-s*8;0.8];
 
128
hh=h1;hh(1:2,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
129
hh=h1;hh(:,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
130
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
131
hh=h1;hh([2 1],1)=x;if cmpr(hh,x([2 1]),eps)<>1 then pause,end;
 
132
hh=h1;hh([2;1],1)=x;if cmpr(hh,x([2 1]),eps)<>1 then pause,end;
 
133
h1=[num/den,den/num];y=0.33*s+1;x=y*y/(y+1);
 
134
hh=h1;hh(1,1)=x;if cmpr(hh,[x,h1(1,2)],eps)<>1 then pause,end;
 
135
x=[num/y,y*y/(y+1)];
 
136
hh=h1;hh(1,1:2)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
137
hh=h1;hh(1,[2 1])=x;if cmpr(hh,x(1,[2,1]),eps)<>1 then pause,end;
 
138
hh=h1;hh(1,[2;1])=x;if cmpr(hh,x(1,[2,1]),eps)<>1 then pause,end;
 
139
hh=h1;hh(1,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
140
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
141
h1=[num/den;den/num];x=[num/y;y*y/(y+1)];
 
142
hh=h1;hh(1:2,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
143
hh=h1;hh(:,1)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
144
hh=h1;hh(:,:)=x;if cmpr(hh,x,eps)<>1 then pause,end;
 
145
hh=h1;hh([2 1],1)=x;if cmpr(hh,x([2 1],1),eps)<>1 then pause,end;
 
146
hh=h1;hh([2;1],1)=x;if cmpr(hh,x([2 1],1),eps)<>1 then pause,end;
 
147
//
 
148
// matrix operations 
 
149
//
 
150
h1=[num/den,den/num];x=[0.3 1.5];
 
151
if cmpr(-h1,(-1)*h1,eps)<>1 then pause,end
 
152
if cmpr(h1+x,[h1(1,1)+x(1,1) h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
153
if cmpr(x+h1,[h1(1,1)+x(1,1) h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
154
if cmpr(h1-x,[h1(1,1)-x(1,1) h1(1,2)-x(1,2)],eps)<>1 then pause,end
 
155
if cmpr(x-h1,[-h1(1,1)+x(1,1) -h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
156
h1=[num/den;den/num];x=[0.3;1.5];
 
157
if cmpr(-h1,(-1)*h1,eps)<>1 then pause,end
 
158
if cmpr(h1+x,[h1(1,1)+x(1,1); h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
159
if cmpr(x+h1,[h1(1,1)+x(1,1); h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
160
if cmpr(h1-x,[h1(1,1)-x(1,1); h1(2,1)-x(2,1)],eps)<>1 then pause,end
 
161
if cmpr(x-h1,[-h1(1,1)+x(1,1); -h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
162
 
 
163
 
 
164
h1=[num/den,den/num];x=[0.3+s, 1.5];
 
165
if cmpr(-h1,(-1)*h1,eps)<>1 then pause,end
 
166
if cmpr(h1+x,[h1(1,1)+x(1,1) h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
167
if cmpr(x+h1,[h1(1,1)+x(1,1) h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
168
if cmpr(h1-x,[h1(1,1)-x(1,1) h1(1,2)-x(1,2)],eps)<>1 then pause,end
 
169
if cmpr(x-h1,[-h1(1,1)+x(1,1) -h1(1,2)+x(1,2)],eps)<>1 then pause,end
 
170
h1=[num/den;den/num];x=[0.3;1.5-3*s];
 
171
if cmpr(-h1,(-1)*h1,eps)<>1 then pause,end
 
172
if cmpr(h1+x,[h1(1,1)+x(1,1); h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
173
if cmpr(x+h1,[h1(1,1)+x(1,1); h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
174
if cmpr(h1-x,[h1(1,1)-x(1,1); h1(2,1)-x(2,1)],eps)<>1 then pause,end
 
175
if cmpr(x-h1,[-h1(1,1)+x(1,1); -h1(2,1)+x(2,1)],eps)<>1 then pause,end
 
176
//
 
177
//
 
178
h1=[num/den,den/num];
 
179
if cmpr([num,den]/den,[num/den,1],eps)<>1 then pause,end
 
180
if cmpr(den\[num,den],[num/den,1],eps)<>1 then pause,end
 
181
if cmpr([num,den]./[den,num],h1,eps)<>1 then pause,end
 
182
if cmpr([den,num].\[num,den],h1,eps)<>1 then pause,end
 
183
h1=[num/den,den/num];x=[0.3 1.5];
 
184
if cmpr(h1/x(1),[h1(1,1)/x(1),h1(1,2)/x(1)],eps)<>1 then pause,end
 
185
if cmpr(x(1)\h1,[h1(1,1)/x(1),h1(1,2)/x(1)],eps)<>1 then pause,end
 
186
if cmpr(h1./x,[h1(1,1)/x(1),h1(1,2)/x(2)],eps)<>1 then pause,end
 
187
if cmpr(x.\h1,[h1(1,1)/x(1),h1(1,2)/x(2)],eps)<>1 then pause,end
 
188
if cmpr(h1*x(1),[h1(1,1)*x(1),h1(1,2)*x(1)],eps)<>1 then pause,end
 
189
if cmpr(h1.*x,[h1(1,1)*x(1),h1(1,2)*x(2)],eps)<>1 then pause,end
 
190
h1=[num/den;den/num];
 
191
if cmpr([num;den]/den,[num/den;1],eps)<>1 then pause,end
 
192
if cmpr(den\[num;den],[num/den;1],eps)<>1 then pause,end
 
193
if cmpr([num;den]./[den;num],h1,eps)<>1 then pause,end
 
194
if cmpr([den;num].\[num;den],h1,eps)<>1 then pause,end
 
195
 
 
196
x=[0.3;1.5];
 
197
if cmpr(h1/x(1),[h1(1,1)/x(1);h1(2,1)/x(1)],eps)<>1 then pause,end
 
198
if cmpr(x(1)\h1,[h1(1,1)/x(1);h1(2,1)/x(1)],eps)<>1 then pause,end
 
199
if cmpr(h1./x,[h1(1,1)/x(1);h1(2,1)/x(2)],eps)<>1 then pause,end
 
200
if cmpr(x.\h1,[h1(1,1)/x(1);h1(2,1)/x(2)],eps)<>1 then pause,end
 
201
if cmpr(h1*x(1),[h1(1,1)*x(1);h1(2,1)*x(1)],eps)<>1 then pause,end
 
202
if cmpr(h1.*x,[h1(1,1)*x(1);h1(2,1)*x(2)],eps)<>1 then pause,end
 
203
if cmpr(h1'*x,h1(1,1)*x(1)+h1(2,1)*x(2),eps)<>1 then pause,end
 
204
if cmpr(h1*x',[h1(1,1)*x(1),h1(1,1)*x(2);h1(2,1)*x(1),h1(2,1)*x(2)],eps)<>1,
 
205
        then pause,end
 
206
h1=[num/den,den/num];x=[0.3+3*s, 1.5];
 
207
if cmpr(h1/x(1),[h1(1,1)/x(1),h1(1,2)/x(1)],eps)<>1 then pause,end
 
208
if cmpr(x(1)\h1,[h1(1,1)/x(1),h1(1,2)/x(1)],eps)<>1 then pause,end
 
209
if cmpr(h1./x,[h1(1,1)/x(1),h1(1,2)/x(2)],eps)<>1 then pause,end
 
210
if cmpr(x.\h1,[h1(1,1)/x(1),h1(1,2)/x(2)],eps)<>1 then pause,end
 
211
if cmpr(h1*x(1),[h1(1,1)*x(1),h1(1,2)*x(1)],eps)<>1 then pause,end
 
212
if cmpr(h1.*x,[h1(1,1)*x(1),h1(1,2)*x(2)],eps)<>1 then pause,end
 
213
if cmpr(h1*x',h1(1,1)*x(1)+h1(1,2)*x(2),eps)<>1 then pause,end
 
214
if cmpr(h1'*x,[h1(1,1)*x(1),h1(1,1)*x(2);h1(1,2)*x(1),h1(1,2)*x(2)],eps)<>1,
 
215
        then pause,end
 
216
h1=[num/den;den/num];x=[0.3+3*s; 1.5];
 
217
if cmpr(h1/x(1),[h1(1,1)/x(1);h1(2,1)/x(1)],eps)<>1 then pause,end
 
218
if cmpr(x(1)\h1,[h1(1,1)/x(1);h1(2,1)/x(1)],eps)<>1 then pause,end
 
219
if cmpr(h1./x,[h1(1,1)/x(1);h1(2,1)/x(2)],eps)<>1 then pause,end
 
220
if cmpr(x.\h1,[h1(1,1)/x(1);h1(2,1)/x(2)],eps)<>1 then pause,end
 
221
if cmpr(h1*x(1),[h1(1,1)*x(1);h1(2,1)*x(1)],eps)<>1 then pause,end
 
222
if cmpr(h1.*x,[h1(1,1)*x(1);h1(2,1)*x(2)],eps)<>1 then pause,end
 
223
if cmpr(h1'*x,h1(1,1)*x(1)+h1(2,1)*x(2),eps)<>1 then pause,end
 
224
if cmpr(h1*x',[h1(1,1)*x(1),h1(1,1)*x(2);h1(2,1)*x(1),h1(2,1)*x(2)],eps)<>1,
 
225
        then pause,end
 
226
h1=[num/den,den/num];x=[0.3/s,1.5-s**2/(1+s**2)];
 
227
if cmpr(h1/x(1,1),[h1(1,1)/x(1,1),h1(1,2)/x(1,1)],eps)<>1 then pause,end
 
228
if cmpr(x(1,1)\h1,[h1(1,1)/x(1,1),h1(1,2)/x(1,1)],eps)<>1 then pause,end
 
229
if cmpr(h1./x,[h1(1,1)/x(1,1),h1(1,2)/x(1,2)],eps)<>1 then pause,end
 
230
if cmpr(x.\h1,[h1(1,1)/x(1,1),h1(1,2)/x(1,2)],eps)<>1 then pause,end
 
231
if cmpr(h1*x(1,1),[h1(1,1)*x(1,1),h1(1,2)*x(1,1)],eps)<>1 then pause,end
 
232
if cmpr(h1.*x,[h1(1,1)*x(1,1),h1(1,2)*x(1,2)],eps)<>1 then pause,end
 
233
 
 
234
h1=[num/den;den/num];x=[0.3/s;1.5-s**2/(1+s**2)];
 
235
if cmpr(h1/x(1,1),[h1(1,1)/x(1,1);h1(2,1)/x(1,1)],eps)<>1 then pause,end
 
236
if cmpr(x(1,1)\h1,[h1(1,1)/x(1,1);h1(2,1)/x(1,1)],eps)<>1 then pause,end
 
237
if cmpr(h1./x,[h1(1,1)/x(1,1);h1(2,1)/x(2,1)],eps)<>1 then pause,end
 
238
if cmpr(x.\h1,[h1(1,1)/x(1,1);h1(2,1)/x(2,1)],eps)<>1 then pause,end
 
239
if cmpr(h1*x(1,1),[h1(1,1)*x(1,1);h1(2,1)*x(1,1)],eps)<>1 then pause,end
 
240
if cmpr(h1.*x,[h1(1,1)*x(1,1);h1(2,1)*x(2,1)],eps)<>1 then pause,end
 
241
if cmpr(h1'*x,h1(1,1)*x(1,1)+h1(2,1)*x(2,1),eps)<>1 then pause,end
 
242
if cmpr(h1*x',[h1(1,1)*x(1,1),h1(1,1)*x(2,1)
 
243
               h1(2,1)*x(1,1),h1(2,1)*x(2,1)],eps)<>1,
 
244
        then pause,end 
 
245
//
 
246
 
 
247
h1=h1*x';x=[1 2;3 4];
 
248
if cmpr(h1/x,h1*inv(x),eps)<>1 then pause,end
 
249
if cmpr(x\h1,inv(x)*h1,eps)<>1 then pause,end
 
250
x=[s s*s+s-1;1 s+1];xi=[1+s,1-s-s*s;-1 s];
 
251
if norm(coeff(invr(x)-xi))>eps then pause,end
 
252
if cmpr(h1/x,h1*xi,eps)<>1 then pause,end
 
253
if cmpr(x\h1,xi*h1,eps)<>1 then pause,end
 
254
if cmpr(x**(-1),xi,eps)<>1 then pause,end
 
255
x=[1/(1+s) 1/(s*(s+1))-1;1 1/s];xi=[1/s,-1/(s*(s+1))+1;-1,1/(1+s)];
 
256
if cmpr(invr(x),xi,eps)<>1 then pause,end
 
257
if cmpr(h1/x,h1*xi,eps)<>1 then pause,end
 
258
if cmpr(x\h1,xi*h1,eps)<>1 then pause,end
 
259
if cmpr((1/(1+s))**3,1/((1+s)**3),eps)<>1 then pause,end
 
260
x=[1/(1+s),1;0 1/s];x3=[x(1,1)**3,x(1,1)**2+x(2,2)*(x(2,2)+x(1,1));0,x(2,2)**3];
 
261
if cmpr(x**3,x3,eps)<>1 then pause,end
 
262
 
 
263
 
 
264
//Bezout
 
265
 
 
266
mode(5)
 
267
//test
 
268
un=poly(1,'s','c');
 
269
zer=0*un;
 
270
s=poly(0,'s');
 
271
[p,q]=bezout(un,s);
 
272
if norm(coeff([un s]*q-[p 0]))>10*%eps then pause,end
 
273
[p,q]=bezout(s,un);
 
274
if norm(coeff([s un]*q-[p 0]))>10*%eps then pause,end
 
275
[p,q]=bezout(un,un);
 
276
if norm(coeff([un un]*q-[p 0]))>10*%eps then pause,end
 
277
//
 
278
[p,q]=bezout(zer,s);
 
279
if norm(coeff([zer s]*q-[p 0]))>10*%eps then pause,end
 
280
[p,q]=bezout(s,zer);
 
281
if norm(coeff([s zer]*q-[p 0]))>10*%eps then pause,end
 
282
[p,q]=bezout(zer,zer);
 
283
if norm(coeff([zer zer]*q-[p 0]))>10*%eps then pause,end
 
284
//
 
285
[p,q]=bezout(zer,un);
 
286
if norm(coeff([zer un]*q-[p 0]))>10*%eps then pause,end
 
287
[p,q]=bezout(un,zer);
 
288
if norm(coeff([un zer]*q-[p 0]))>10*%eps then pause,end
 
289
//
 
290
 
 
291
//simple
 
292
a=poly([1 2 3],'z');
 
293
b=poly([4 1],'z');
 
294
 
 
295
[p q]=bezout(a,b);
 
296
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
297
if norm(coeff(dt/dt0-1))>10*%eps then pause,end
 
298
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
299
if norm(coeff([p 0]*qi-[a b]))>100*%eps then pause,end
 
300
 
 
301
//more difficult
 
302
 
 
303
b1=poly([4 1+1000*%eps],'z');del=10*norm(coeff(b1-b));
 
304
[p,q]=bezout(a,b1);
 
305
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
306
if norm(coeff(dt/dt0-1))>del then pause,end
 
307
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
308
if norm(coeff([p 0]*qi-[a b1]))>del then pause,end
 
309
 
 
310
b1=poly([4 1+.001],'z');del=10*norm(coeff(b1-b));
 
311
[p,q]=bezout(a,b1);
 
312
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
313
if norm(coeff(dt/dt0-1))>100000*%eps then pause,end
 
314
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
315
if norm(coeff([p 0]*qi-[a b1]))>100000*%eps then pause,end
 
316
 
 
317
 
 
318
b1=poly([4 1+10000*%eps],'z');del=10*norm(coeff(b1-b));
 
319
[p,q]=bezout(a,b1);
 
320
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
321
if norm(coeff(dt/dt0-1))>del then pause,end
 
322
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
323
if norm(coeff([p 0]*qi-[a b1]))>del then pause,end
 
324
//
 
325
z=poly(0,'z');
 
326
num = 0.99999999999999922+z*(-4.24619123578530730+..
 
327
      z*(10.0503215227460350+z*(-14.6836461849931740+..
 
328
      z*(13.924822877119892+z*(-5.63165998008533460+..
 
329
      z*(-5.63165998008530710+z*(13.9248228771198730+..
 
330
      z*(-14.683646184993167+z*(10.0503215227460330+..
 
331
      z*(-4.24619123578530910+z*(0.99999999999999989)))))))))));
 
332
den = -0.17323463717888873+z*(1.91435457459735380+..
 
333
      z*(-9.90126732768255560+z*(31.6286096652930410+..
 
334
      z*(-69.3385546880304280+z*(109.586435800377690+..
 
335
      z*(-127.516160100808290+z*(109.388684898145950+..
 
336
      z*(-67.92645394857864+z*(29.1602681026148110+..
 
337
      z*(-7.8212498781094952+z*(0.99999999999999989)))))))))));
 
338
//
 
339
[p,q]=bezout(num,den);del=1.d-4;
 
340
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
341
if norm(coeff(dt/dt0-1))>del then pause,end
 
342
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
343
// JPC 
 
344
del=3*del
 
345
if norm(coeff([p 0]*qi-[num den]))>del then pause,end
 
346
if degree(p)>0 then pause,end
 
347
// une autre "solution" telle que l'erreur directe est petite mais l'erreur
 
348
// inverse grande
 
349
q1=[]
 
350
q1(1,1)=poly([0.011533588674 , 3.570224012502 , -28.78817740957 ,...
 
351
            102.9479419903, -210.8258579715 , 266.2028963639 ,...
 
352
           -207.427018299 , 92.74478640319, -18.5259652457],'z','c');
 
353
q1(1,2)=poly([0.000270220664 , -0.002465565223 , 0.010039706635 ,...
 
354
           -0.023913827007, 0.036387144032 , -0.036175176058 ,...
 
355
            0.022954475171 , -0.008514798968,  0.001417382492],'z','c');
 
356
q1(2,1)=poly([0.000639302018 , 20.502472606 , -26.66621972106 ,...
 
357
            39.74001534015,  -5.945830824342 , -7.973226036298 ,...
 
358
            39.84118622788 , -26.51337424414, 18.5259652457],'z','c');
 
359
q1(2,2) =poly( [0.001562930385 , -0.003589174974 , 0.005076237479 ,...
 
360
            -0.003483568682,  -0.000135940266 , 0.003651509121 ,...
 
361
            -0.005059502869 , 0.003447573440, -0.001417382492],'z','c');
 
362
p1 =poly( [0.011422839421 , -0.029264363070 , 0.030070175223 ,...
 
363
         -0.012596066108],'z','c');
 
364
 
 
365
//
 
366
//simplification
 
367
num =[0.03398330733500143,-0.20716057008572933,0.64660689206696986,...
 
368
     -1.97665462134021790,3.38751027286891300,-3.58940006392108120,...
 
369
      5.09956159043231680,5.2514918861834694,1.00000000000000020];
 
370
den = [0.03398330733500360,-0.20716057008816283,0.64660689204312,...
 
371
      -1.97665462079896660,3.38751027286891300,-3.58940006392108350,...
 
372
       5.09956159043232040,5.2514918861834712,1];
 
373
num=poly(num,'z','c');
 
374
den=poly(den,'z','c');
 
375
[p,q]=bezout(num,den);del=1.d-8;
 
376
dt=q(1,1)*q(2,2)-q(1,2)*q(2,1);dt0=coeff(dt,0);
 
377
if norm(coeff(dt/dt0-1))>del then pause,end
 
378
qi=[q(2,2) -q(1,2);-q(2,1) q(1,1)]/dt0;
 
379
if norm(coeff([p 0]*qi-[num den]))>del then pause,end
 
380
if degree(p)<8 then pause,end