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