Mercurial > forge

comparison main/special-matrix/pascal.m @ 0:6b33357c7561 octave-forge

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Initial revision
author pkienzle
date Wed, 10 Oct 2001 19:54:49 +0000
parents
children 2ac2777b30bc
comparison
equal deleted inserted replaced
-1:000000000000 0:6b33357c7561
1 ## Copyright (C) 1999 Peter Ekberg
2 ##
3 ## This program is free software; you can redistribute it and/or modify it
4 ## under the terms of the GNU General Public License as published by
5 ## the Free Software Foundation; either version 2, or (at your option)
6 ## any later version.
7 ##
8 ## This program is distributed in the hope that it will be useful, but
9 ## WITHOUT ANY WARRANTY; without even the implied warranty of
10 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 ## General Public License for more details.
12 ##
13 ## You should have received a copy of the GNU General Public License
14 ## along with this program; see the file COPYING. If not, write to the Free
15 ## Software Foundation, 59 Temple Place - Suite 330, Boston, MA
16 ## 02111-1307, USA.
17
18 ## usage: pascal (n, t)
19 ##
20 ## Return the Pascal matrix of order n if t=0.
21 ## t defaults to 0.
22 ## Return lower triangular Cholesky factor of the Pascal matrix if t=1.
23 ## Return a transposed and permuted version of pascal(n,1) if t=2.
24 ##
25 ## pascal(n,1)^2 == eye(n)
26 ## pascal(n,2)^3 == eye(n)
27 ##
28 ## See also: hankel, vander, sylvester_matrix, hilb, invhilb, toeplitz
29 ## hadamard, wilkinson, rosser, compan
30
31 ## Author: peda
32
33 function retval = pascal (n, t)
34
35 if (nargin > 2) || (nargin == 0)
36 usage ("pascal (n, t)");
37 endif
38
39 if (nargin == 1)
40 t = 0;
41 endif
42
43 if !is_scalar (n) || !is_scalar (t)
44 error ("pascal expecting scalar arguments, found something else");
45 endif
46
47 retval = diag((-1).^[0:n-1]);
48 retval(:,1) = ones(n, 1);
49
50 for j=2:n-1
51 for i=j+1:n
52 retval(i,j) = retval(i-1,j) - retval(i-1,j-1);
53 endfor
54 endfor
55
56 if (t==0)
57 retval = retval*retval';
58 elseif (t==2)
59 retval = retval';
60 retval = retval(n:-1:1,:);
61 retval(:,n) = -retval(:,n);
62 retval(n,:) = -retval(n,:);
63 if (rem(n,2) != 1)
64 retval = -retval;
65 endif
66 endif
67
68 endfunction