unmaintained: New home for misfit Factor vocabularies.

[factor-unmaintained.git] / adsoda / solution2 / solution2.factor

diff --git a/adsoda/solution2/solution2.factor b/adsoda/solution2/solution2.factor
new file mode 100644 (file)
index 0000000..9b75e05
--- /dev/null
+++ b/adsoda/solution2/solution2.factor
@@ -0,0 +1,126 @@
+USING: kernel
+sequences
+namespaces
+
+math
+math.vectors
+math.matrices
+;
+IN: adsoda.solution2
+
+! -------------------
+! correctif solution
+! ---------------
+SYMBOL: matrix
+: MIN-VAL-adsoda ( -- x ) 0.00000001
+! 0.000000000001 
+;
+
+: zero? ( x -- ? ) 
+    abs MIN-VAL-adsoda <
+;
+
+! [ number>string string>number ] map 
+
+: with-matrix ( matrix quot -- )
+    [ swap matrix set call matrix get ] with-scope ; inline
+
+: nth-row ( row# -- seq ) matrix get nth ;
+
+: change-row ( row# quot -- seq ) ! row# quot -- | quot: seq -- seq )
+    matrix get swap change-nth ; inline
+
+: exchange-rows ( row# row# -- ) matrix get exchange ;
+
+: rows ( -- n ) matrix get length ;
+
+: cols ( -- n ) 0 nth-row length ;
+
+: skip ( i seq quot -- n )
+    over [ find-from drop ] dip length or ; inline
+
+: first-col ( row# -- n )
+    ! First non-zero column
+    0 swap nth-row [ zero? not ] skip ;
+
+: clear-scale ( col# pivot-row i-row -- n )
+    [ over ] dip nth dup zero? [
+        3drop 0
+    ] [
+        [ nth dup zero? ] dip swap [
+            2drop 0
+        ] [
+            swap / neg
+        ] if
+    ] if ;
+
+: (clear-col) ( col# pivot-row i -- )
+    [ [ clear-scale ] 2keep [ n*v ] dip v+ ] change-row ;
+
+: rows-from ( row# -- slice )
+    rows dup <slice> ;
+
+: clear-col ( col# row# rows -- )
+    [ nth-row ] dip [ [ 2dup ] dip (clear-col) ] each 2drop ;
+
+: do-row ( exchange-with row# -- )
+    [ exchange-rows ] keep
+    [ first-col ] keep
+    dup 1 + rows-from clear-col ;
+
+: find-row ( row# quot -- i elt )
+    [ rows-from ] dip find ; inline
+
+: pivot-row ( col# row# -- n )
+    [ dupd nth-row nth zero? not ] find-row 2nip ;
+
+: (echelon) ( col# row# -- )
+    over cols < over rows < and [
+        2dup pivot-row [ over do-row 1 + ] when*
+        [ 1 + ] dip (echelon)
+    ] [
+        2drop
+    ] if ;
+
+: echelon ( matrix -- matrix' )
+    [ 0 0 (echelon) ] with-matrix ;
+
+: nonzero-rows ( matrix -- matrix' )
+    [ [ zero? ] all? ] reject ;
+
+: null/rank ( matrix -- null rank )
+    echelon dup length swap nonzero-rows length [ - ] keep ;
+
+: leading ( seq -- n elt ) [ zero? not ] find ;
+
+: reduced ( matrix' -- matrix'' )
+    [
+        rows <reversed> [
+            dup nth-row leading drop
+            dup [ swap dup clear-col ] [ 2drop ] if
+        ] each
+    ] with-matrix ;
+
+: basis-vector ( row col# -- )
+    [ clone ] dip
+    [ swap nth neg recip ] 2keep
+    [ 0 spin set-nth ] 2keep
+    [ n*v ] dip
+    matrix get set-nth ;
+
+: nullspace ( matrix -- seq )
+    echelon reduced dup empty? [
+        dup first length identity-matrix [
+            [
+                dup leading drop
+                dup [ basis-vector ] [ 2drop ] if
+            ] each
+        ] with-matrix flip nonzero-rows
+    ] unless ;
+
+: 1-pivots ( matrix -- matrix )
+    [ dup leading nip [ recip v*n ] when* ] map ;
+
+: solution ( matrix -- matrix )
+    echelon nonzero-rows reduced 1-pivots ;
+