Mercurial > octave-nkf

comparison scripts/linear-algebra/housh.m @ 13048:c5c94b63931f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
codesprint: linear algebra tests: cross, housh, planerot, qzhess, rref
author Roman Belov <romblv@gmail.com>
date Sat, 03 Sep 2011 18:40:46 +0400
parents fd0a3ac60b0e
children 72c96de7a403
comparison
equal deleted inserted replaced
13047:69a4609e61e2 13048:c5c94b63931f
21 ## Compute Householder reflection vector @var{housv} to reflect @var{x} 21 ## Compute Householder reflection vector @var{housv} to reflect @var{x}
22 ## to be the j-th column of identity, i.e., 22 ## to be the j-th column of identity, i.e.,
23 ## 23 ##
24 ## @example 24 ## @example
25 ## @group 25 ## @group
26 ## (I - beta*housv*housv')x = norm(x)*e(j) if x(1) < 0, 26 ## (I - beta*housv*housv')x = norm(x)*e(j) if x(j) < 0,
27 ## (I - beta*housv*housv')x = -norm(x)*e(j) if x(1) >= 0 27 ## (I - beta*housv*housv')x = -norm(x)*e(j) if x(j) >= 0
28 ## @end group 28 ## @end group
29 ## @end example 29 ## @end example
30 ## 30 ##
31 ## @noindent 31 ## @noindent
32 ## Inputs 32 ## Inputs
89 endif 89 endif
90 zer = (beta == 0); 90 zer = (beta == 0);
91 endif 91 endif
92 92
93 endfunction 93 endfunction
94
95 %!test
96 %! x = [1 2 3]';
97 %! j = 3;
98 %! [hv, b, z] = housh(x, j, 0);
99 %! r = (eye(3) - b*hv*hv') * x;
100 %! d = - norm(x) * [0 0 1]';
101 %! assert(r, d, 2e-8);
102 %! assert(z, 0, 2e-8);
103
104 %!test
105 %! x = [7 -3 1]';
106 %! j = 2;
107 %! [hv, b, z] = housh(x, j, 0);
108 %! r = (eye(3) - b*hv*hv') * x;
109 %! d = norm(x) * [0 1 0]';
110 %! assert(r, d, 2e-8);
111 %! assert(z, 0, 2e-8);
112
113 %!test
114 %! x = [1 0 0]';
115 %! j = 1;
116 %! [hv, b, z] = housh(x, j, 10);
117 %! r = (eye(3) - b*hv*hv') * x;
118 %! d = norm(x) * [1 0 0]';
119 %! assert(r, d, 2e-8);
120 %! assert(z, 1, 2e-8);
121
122 %!test
123 %! x = [5 0 4 1]';
124 %! j = 2;
125 %! [hv, b, z] = housh(x, j, 0);
126 %! r = (eye(4) - b*hv*hv') * x;
127 %! d = - norm(x) * [0 1 0 0]';
128 %! assert(r, d, 2e-8);
129 %! assert(z, 0, 2e-8);
130
131 %!error housh([0]);
132 %!error housh();
133