1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
|
function [x,w] = nnls(E,f)
// Non Negative Least Squares (nnls) for Ex=f with the constraint x>=0
// Reference
// Lawson, C.L. and R.J. Hanson, Solving Least Squares Problems,
// Prentice-Hall, 1974, Chapter 23, p. 161.
m2 = size(E,1);
n = size(E,2);
iterMax = 3*n;
x = zeros(n,1);
z = zeros(n,1);
w = zeros(n,1);
wActive = w;
// P indicates if a variable is inactive
P = x~=0; // setting all to False
// z indiactes if a variable is active
Z = x==0; // setting all to True
Ep = zeros(size(E));
// Step 2
w = E'*(f-E*x);
iter = 0;
while or(Z) & or(w(Z)>0)
// Step 4
wActive(P) = -%inf;
wActive(Z) = w(Z);
[maxval,maxindex] = max(wActive);
// Step 5
P(maxindex) = %T;
Z(maxindex) = %F;
// Step 6
z(P)=E(:,P)\f;
z(Z)=0;
while or(z(P)<=0)
iter = iter+1;
if iter>iterM`ax then
x = z;
return;
end
// Step 8
Q = (z <= 0) & P;
// Step 9
alpha = min(x(P)./(x(P)-z(P)));
// Step 10
x = x+alpha*(z-x);
indices = find(x>=0);
P(indices) = %T;
W(indices) = %F;
// Step 6
z(P)=E(:,P)\f;
z(Z)=0;
end
x = z;
// Step 2
w = E'*(f-E*x);
wActive = w(Z);
end
endfunction
|