Permalink
Browse files

lev: fix for multiple right-hand-side dimensions

  • Loading branch information...
1 parent e13d184 commit 5af2b0bd67716ec070804b33e7f40d7989567a8a @drizzd committed Jan 16, 2010
Showing with 40 additions and 39 deletions.
  1. +5 −36 blk.m
  2. +32 −0 blk_row.m
  3. +3 −3 solve_lev_block.m
View
41 blk.m
@@ -1,39 +1,8 @@
-function y = blk(x, d, a, b, step)
-% BLK extract block
+function y = blk(varargin)
+% BLK_ROW extract block
% y = blk(x, d, a = 1, b = 1, step = 1)
%
-if nargin < 5
- step = 1;
-end
-if nargin < 4
- b = 1;
-end
-
-n = size(x, 1);
-if mod(n, d) ~= 0
- error('size of x is not a multiple of d')
-end
-
-if step ~= -1 && step ~= 1
- error('invalid step size')
-end
-
-if step == -1
- a = a - b + 1;
-end
-
-y = x((a-1)*d+1:(a+b-1)*d, 1:min(size(x, 2), d));
-
-if step == -1
- yr = zeros(size(y));
- S = size(y, 1);
- if mod(S, d) ~= 0
- error('this should not happen');
- end
- S = S/d;
- for k = 1:S
- yr((k-1)*d+1:k*d, :) = y((S-k)*d+1:(S-k+1)*d, :);
- end
- y = yr;
-end
+y = blk_row(varargin{:});
+d = varargin{2};
+y = y(:, 1:d);
View
@@ -0,0 +1,32 @@
+function y = blk_row(x, d, a, m, step)
+% BLK_ROW extract block rows
+% y = blk_row(x, d, a = 1, m = 1, step = 1)
+%
+
+if nargin < 5
+ step = 1;
+end
+if nargin < 4
+ m = 1;
+end
+
+n = size(x, 1);
+if mod(n, d) ~= 0
+ error('size of x is not a multiple of d')
+end
+
+m = m * abs(step);
+
+if step < 0
+ a = a - m + 1;
+ u = m;
+ v = 1;
+else
+ u = 1;
+ v = m;
+end
+
+idx = reshape((a-1)*d+1:(a+m-1)*d, d, m);
+idx = idx(:, u:step:v);
+
+y = x(idx, :);
View
@@ -32,16 +32,16 @@
c = b;
a = -b\blk(A.', d, 2).';
v = -c\blk(A, d, 2);
-x = b\blk(y, d, 1);
+x = b\blk_row(y, d, 1);
Y = a;
W = v;
N = n/d;
for k = 1:N-1
b = b - b * a * v;
c = c - c * v * a;
- u = blk(y, d, k+1);
+ u = blk_row(y, d, k+1);
for l = 2:k+1
- u = u - blk(A, d, l) * blk(x, d, k - l + 2);
+ u = u - blk(A, d, l) * blk_row(x, d, k - l + 2);
end
mu = c\u;
x = [x + blk(Y, d, k, k, -1) * mu; mu];

0 comments on commit 5af2b0b

Please sign in to comment.