~ubuntu-branches/ubuntu/hoary/scilab/hoary

« back to all changes in this revision

Viewing changes to demos/basic/intro/dem01.dem

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2005-01-09 22:58:21 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20050109225821-473xr8vhgugxxx5j
Tags: 3.0-12
changed configure.in to build scilab's own malloc.o, closes: #255869

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
mode(7)
 
2
//               SCILAB OBJECTS
 
3
//               1. SCALARS
 
4
// Copyright INRIA
 
5
a=1               //real constant
 
6
1==1              //boolean
 
7
'string'          //character string
 
8
z=poly(0,'z')     // polynomial with variable 'z' and with one root at zero
 
9
p=1+3*z+4.5*z^2   //polynomial 
 
10
r=z/p             //rational
 
11
 
 
12
//                2. MATRICES
 
13
A=[a+1 2 3
 
14
     0 0 atan(1)
 
15
     5 9 -1]      //3 x 3 constant matrix
 
16
 
 
17
b=[%t,%f]         //1 x 2 boolean matrix
 
18
 
 
19
Mc=['this','is';
 
20
    'a' ,'matrix']   //2 x 2 matrix of strings
 
21
 
 
22
Mp=[p,1-z;
 
23
    1,z*p]        //2 x 2 polynomial matrix
 
24
 
 
25
F=Mp/poly([1+%i 1-%i 1],'z')   //rational matrix
 
26
 
 
27
Sp=sparse([1,2;4,5;3,10],[1,2,3])   //sparse matrix
 
28
 
 
29
Sp(1,10)==Sp(1,1)                   //boolean sparse matrix
 
30
 
 
31
//                 3. LISTS
 
32
L=list(a,-(1:5), Mp,['this','is';'a','list'])   //list
 
33
L(2)(3)     //sub-entry in list
 
34
Lt=tlist(['mylist','color','position','weight'],'blue',[0,1],10)  //typed-list
 
35
Lt('color')      //extracting
 
36
Lt('weight')     //extracting
 
37
A=diag([2,3,4]);B=[1 0;0 1;0 0];C=[1 -1 0];D=0*C*B;x0=[0;0;0];
 
38
Sl=syslin('c',A,B,C,D,x0)    //Standard state-space linear system 
 
39
Sl("A"), Sl("C")             //Retrieving elements of a typed list
 
40
Slt=ss2tf(Sl)                // Transfer matrix
 
41
Slt('num'), Slt('den')
 
42
//                  OPERATIONS
 
43
v=1:5;W=v'*v                 //constant matrix multiplication
 
44
W(1,:)                       //extracting first row
 
45
W(:,$)                       //extracting last column
 
46
Mp'*Mp+eye()                   //polynomial matrix
 
47
Mp1=Mp(1,1)+4.5*%i           //complex
 
48
Fi=C*(z*eye()-A)^(-1)*B;       //transfer function evaluation
 
49
F(:,1)*Fi                    //operations with rationals
 
50
M=[Mp -Mp; Mp' Mp+eye()]       //concatenation of polynomial matrices
 
51
[Fi, Fi(:,1)]                // ... or rationals
 
52
F=syslin('c',F);             
 
53
Num=F('num');Den=F('den');           //operation on transfer matrix
 
54
 
 
55
//                  SOME NUMERICAL PRIMITIVES
 
56
inv(A)                       //Inverse
 
57
inv(Mp)                      //Inverse
 
58
inv(Sl*Sl')                  //Product of two linear systems and inverse
 
59
w=ss2tf(ans)                 //Transfer function representation
 
60
w1=inv(ss2tf(Sl)*ss2tf(Sl'))    //Product of two transfer functions and inverse
 
61
clean(w-w1)                 
 
62
A=rand(3,3);;B=rand(3,1);n=contr(A,B)                 //Controllability
 
63
K=ppol(A,B,[-1-%i -1+%i -1])        //Pole placement
 
64
poly(A-B*K,'z')-poly([-1-%i -1+%i -1],'z')    //Check...
 
65
 
 
66
s=sin(0:0.1:5*%pi);
 
67
ss=fft(s(1:128),-1);        //FFT
 
68
xbasc();
 
69
plot2d3("enn",1,abs(ss)');      //simple plot
 
70
 
 
71
 
 
72
//              ON LINE DEFINITION OF FUNCTION
 
73
deff('[x]=fact(n)','if n==0 then x=1,else x=n*fact(n-1),end')
 
74
10+fact(5)
 
75
//                    OPTIMIZATION
 
76
deff('[f,g,ind]=rosenbro(x,ind)', 'a=x(2)-x(1)^2 , b=1-x(2) ,...
 
77
f=100.*a^2 + b^2 , g(1)=-400.*x(1)*a , g(2)=200.*a -2.*b ');
 
78
[f,x,g]=optim(rosenbro,[2;2],'qn')
 
79
 
 
80
//                   SIMULATION
 
81
a=rand(3,3)
 
82
e=expm(a)
 
83
deff('[ydot]=f(t,y)','ydot=a*y');
 
84
e(:,1)-ode([1;0;0],0,1,f)
 
85
 
 
86
//                  SYSTEM DEFINITION
 
87
s=poly(0,'s');
 
88
h=[1/s,1/(s+1);1/s/(s+1),1/(s+2)/(s+2)]
 
89
w=tf2ss(h);
 
90
ss2tf(w)
 
91
h1=clean(ans)
 
92
 
 
93
//             EXAMPLE: SECOND ORDER SYSTEM ANALYSIS
 
94
sl=syslin('c',1/(s*s+0.2*s+1))
 
95
instants=0:0.05:20;
 
96
//             step response:
 
97
y=csim('step',instants,sl);
 
98
xbasc();plot2d(instants',y')
 
99
//             Delayed step response
 
100
deff('[in]=u(t)','if t<3 then in=0;else in=1;end');
 
101
y1=csim(u,instants,sl);plot2d(instants',y1');
 
102
 
 
103
//             Impulse response;
 
104
yi=csim('imp',instants,sl);xbasc();plot2d(instants',yi');
 
105
yi1=csim('step',instants,s*sl);plot2d(instants',yi1');
 
106
 
 
107
//              Discretization
 
108
dt=0.05;
 
109
sld=dscr(tf2ss(sl),0.05);
 
110
 
 
111
//               Step response
 
112
u=ones(instants);
 
113
yyy=flts(u,sld);
 
114
xbasc();plot(instants,yyy)
 
115
 
 
116
//              Impulse response
 
117
u=0*ones(instants);u(1)=1/dt;
 
118
yy=flts(u,sld);
 
119
xbasc();plot(instants,yy)
 
120
 
 
121
//            system interconnexion
 
122
w1=[w,w];
 
123
clean(ss2tf(w1))
 
124
w2=[w;w];
 
125
clean(ss2tf(w2))
 
126
 
 
127
//               change of variable
 
128
z=poly(0,'z');
 
129
horner(h,(1-z)/(1+z))  //bilinear transform
 
130
 
 
131
//                 PRIMITIVES
 
132
H=[1.    1.    1.    0.;
 
133
   2.   -1.    0.    1;
 
134
   1.    0.    1.    1.;
 
135
   0.    1.    2.   -1];
 
136
 
 
137
ww=spec(H)
 
138
 
 
139
//             STABLE SUBSPACES
 
140
[X,d]=schur(H,'cont');
 
141
X'*H*X
 
142
 
 
143
[X,d]=schur(H,'disc');
 
144
X'*H*X
 
145
 
 
146
//Selection of user-defined eigenvalues (# 3 and 4 here);
 
147
function [flg]=sel(x)
 
148
  flg=%f
 
149
  if abs(x-ww(3))<0.0001|abs(x-ww(4))<0.00001 then flg=%t,end 
 
150
endfunction
 
151
 
 
152
[X,d]=schur(H,sel)
 
153
 
 
154
X'*H*X
 
155
 
 
156
//               With matrix pencil
 
157
function [flg]=sel(x,y)
 
158
  flg=%f
 
159
  if abs(x/y-ww(3))<0.0001|abs(x/y-ww(4))<0.00001 then flg=%t,end 
 
160
endfunction
 
161
 
 
162
[X,d]=schur(H,eye(H),sel)
 
163
X'*H*X
 
164
 
 
165
//            block diagonalization
 
166
[ab,x,bs]=bdiag(H);
 
167
 
 
168
inv(x)*H*x
 
169
 
 
170
//                     Matrix pencils
 
171
E=rand(3,2)*rand(2,3);
 
172
A=rand(3,2)*rand(2,3);
 
173
s=poly(0,'s');
 
174
 
 
175
w=det(s*E-A)   //determinant
 
176
[al,be]=spec(A,E);
 
177
al./(be+%eps*ones(be))
 
178
roots(w)
 
179
[Ns,d]=coffg(s*E-A);    //inverse of polynomial matrix;
 
180
clean(Ns/d*(s*E-A))
 
181
[Q,M,i1]=pencan(E,A);   // Canonical form;
 
182
clean(M*E*Q)
 
183
clean(M*A*Q)
 
184
 
 
185
//           PAUSE-RESUME
 
186
write(%io(2),'pause command...');
 
187
write(%io(2),'TO CONTINUE...');
 
188
write(%io(2),'ENTER ''resume (or return) or click on resume!!''');
 
189
//pause;
 
190
 
 
191
//           CALLING EXTERNAL ROUTINE
 
192
 
 
193
foo=['void foo(double *a,double *b,double *c)';
 
194
     '{ *c = *a + *b; }'  ];
 
195
 
 
196
// we use TMPDIR for compilation 
 
197
        
 
198
if ~c_link('foo') then
 
199
  path = getcwd(); 
 
200
  chdir(TMPDIR); 
 
201
  mputl(foo,'foo.c');
 
202
  ilib_for_link(['foo'],'foo.o',[],"c");
 
203
  exec loader.sce 
 
204
  chdir(path) 
 
205
end     
 
206
 
 
207
//5+7 by C function
 
208
call('foo',5,1,'d',7,2,'d','out',[1,1],3,'d')