Fix http help lint
[factor/jcg.git] / extra / adsoda / solution2 / solution2.factor
blob3e0648128de9746937e1e4b4a87b6f33212693be
1 USING: kernel\r
2 sequences\r
3 namespaces\r
4 \r
5 math\r
6 math.vectors\r
7 math.matrices\r
8 ;\r
9 IN: adsoda.solution2\r
11 ! -------------------\r
12 ! correctif solution\r
13 ! ---------------\r
14 SYMBOL: matrix\r
15 : MIN-VAL-adsoda ( -- x ) 0.00000001\r
16 ! 0.000000000001 \r
17 ;\r
19 : zero? ( x -- ? ) \r
20     abs MIN-VAL-adsoda <\r
21 ;\r
23 ! [ number>string string>number ] map \r
25 : with-matrix ( matrix quot -- )\r
26     [ swap matrix set call matrix get ] with-scope ; inline\r
28 : nth-row ( row# -- seq ) matrix get nth ;\r
30 : change-row ( row# quot -- seq ) ! row# quot -- | quot: seq -- seq )\r
31     matrix get swap change-nth ; inline\r
33 : exchange-rows ( row# row# -- ) matrix get exchange ;\r
35 : rows ( -- n ) matrix get length ;\r
37 : cols ( -- n ) 0 nth-row length ;\r
39 : skip ( i seq quot -- n )\r
40     over [ find-from drop ] dip length or ; inline\r
42 : first-col ( row# -- n )\r
43     #! First non-zero column\r
44     0 swap nth-row [ zero? not ] skip ;\r
46 : clear-scale ( col# pivot-row i-row -- n )\r
47     [ over ] dip nth dup zero? [\r
48         3drop 0\r
49     ] [\r
50         [ nth dup zero? ] dip swap [\r
51             2drop 0\r
52         ] [\r
53             swap / neg\r
54         ] if\r
55     ] if ;\r
57 : (clear-col) ( col# pivot-row i -- )\r
58     [ [ clear-scale ] 2keep [ n*v ] dip v+ ] change-row ;\r
60 : rows-from ( row# -- slice )\r
61     rows dup <slice> ;\r
63 : clear-col ( col# row# rows -- )\r
64     [ nth-row ] dip [ [ 2dup ] dip (clear-col) ] each 2drop ;\r
66 : do-row ( exchange-with row# -- )\r
67     [ exchange-rows ] keep\r
68     [ first-col ] keep\r
69     dup 1+ rows-from clear-col ;\r
71 : find-row ( row# quot -- i elt )\r
72     [ rows-from ] dip find ; inline\r
74 : pivot-row ( col# row# -- n )\r
75     [ dupd nth-row nth zero? not ] find-row 2nip ;\r
77 : (echelon) ( col# row# -- )\r
78     over cols < over rows < and [\r
79         2dup pivot-row [ over do-row 1+ ] when*\r
80         [ 1+ ] dip (echelon)\r
81     ] [\r
82         2drop\r
83     ] if ;\r
85 : echelon ( matrix -- matrix' )\r
86     [ 0 0 (echelon) ] with-matrix ;\r
88 : nonzero-rows ( matrix -- matrix' )\r
89     [ [ zero? ] all? not ] filter ;\r
91 : null/rank ( matrix -- null rank )\r
92     echelon dup length swap nonzero-rows length [ - ] keep ;\r
94 : leading ( seq -- n elt ) [ zero? not ] find ;\r
96 : reduced ( matrix' -- matrix'' )\r
97     [\r
98         rows <reversed> [\r
99             dup nth-row leading drop\r
100             dup [ swap dup clear-col ] [ 2drop ] if\r
101         ] each\r
102     ] with-matrix ;\r
104 : basis-vector ( row col# -- )\r
105     [ clone ] dip\r
106     [ swap nth neg recip ] 2keep\r
107     [ 0 spin set-nth ] 2keep\r
108     [ n*v ] dip\r
109     matrix get set-nth ;\r
111 : nullspace ( matrix -- seq )\r
112     echelon reduced dup empty? [\r
113         dup first length identity-matrix [\r
114             [\r
115                 dup leading drop\r
116                 dup [ basis-vector ] [ 2drop ] if\r
117             ] each\r
118         ] with-matrix flip nonzero-rows\r
119     ] unless ;\r
121 : 1-pivots ( matrix -- matrix )\r
122     [ dup leading nip [ recip v*n ] when* ] map ;\r
124 : solution ( matrix -- matrix )\r
125     echelon nonzero-rows reduced 1-pivots ;\r