[Racket] Ferti hydroponic nutrient solver, redux.
Update max-w-in-R algorithm.
services/nnls.rkt
@@ -34,7 +34,7 @@
34
34
(lawson-hanson-1974 fertilizer-product-matrix deficits error-threshold))
35
35
36
36
;; Algorithm lifted from the Wikipedia article on NNLS
37
Removed:
(define/contract (lawson-hanson-1974 A y ε)
37
Added:
(define (lawson-hanson-1974 A y ε)
38
38
;;;;;;;;;
39
39
;; Inputs
40
40
;;;;;;;;;
@@ -63,16 +63,14 @@
63
63
(define (compute-w x)
64
64
(matrix* (matrix-transpose A) (matrix- y (matrix* A x))))
65
65
66
Removed:
;; max over j in R of w_j, returning (values max-val j*)
66
Added:
;; max over j in R of w_j.
67
67
(define (max-w-in-R w)
68
Removed:
(for/fold ([max-val -inf.0]
69
Removed:
[max-j #f])
70
Removed:
([j (in-set R)])
71
Removed:
(define v (colv-ref w j))
72
Removed:
(if (> v max-val)
73
Removed:
(values v j)
74
Removed:
(values max-val max-j))))
68
Added:
(if (set-empty? R)
69
Added:
(values -inf.0 #f)
70
Added:
(let ([max-j (argmax (λ (j) (colv-ref w j)) (set->list R))])
71
Added:
(values (colv-ref w max-j) max-j))))
75
72
73
Added:
76
74
;; Build full candidate vector s from current P:
77
75
;; s_P = (A_Pᵀ A_P)⁻¹ A_Pᵀ y, s_R = 0
78
76
(define (make-s-from-P)
@@ -124,20 +122,16 @@
124
122
;; Inner loop: adjust until s_P > 0
125
123
(let inner-loop ()
126
124
(define s (make-s-from-P))
127
Removed:
128
Removed:
;; min(s_P)
129
125
(define min-sP
130
126
(if (set-empty? P)
131
127
+inf.0
132
128
(for/fold ([mn +inf.0]) ([j (in-set P)])
133
129
(min mn (colv-ref s j)))))
134
Removed:
135
130
(cond
136
131
;; If all s_P > 0 (or P empty), accept s as new x and go back to outer loop
137
132
[(or (set-empty? P) (> min-sP 0))
138
133
(set! x s)
139
134
(outer-loop)]
140
Removed:
141
135
[else
142
136
;; Compute α = min_{i in P, s_i <= 0} x_i / (x_i - s_i)
143
137
(define α