UNPKG

3.52 kBtext/coffeescriptView Raw
1#-----------------------------------------------------------------------------
2#
3# Generate all divisors of a term
4#
5# Input: Term on stack (factor * factor * ...)
6#
7# Output: Divisors on stack
8#
9#-----------------------------------------------------------------------------
10
11
12
13
14divisors = ->
15 i = 0
16 h = 0
17 n = 0
18 save()
19 h = tos - 1
20 divisors_onstack()
21 n = tos - h
22
23 #qsort(stack + h, n, sizeof (U *), __cmp)
24 subsetOfStack = stack.slice(h,h+n)
25 subsetOfStack.sort(cmp_expr)
26 stack = stack.slice(0,h).concat(subsetOfStack).concat(stack.slice(h+n))
27
28 p1 = alloc_tensor(n)
29 p1.tensor.ndim = 1
30 p1.tensor.dim[0] = n
31 for i in [0...n]
32 p1.tensor.elem[i] = stack[h + i]
33 moveTos h
34 push(p1)
35 restore()
36
37divisors_onstack = ->
38 h = 0
39 i = 0
40 k = 0
41 n = 0
42
43 save()
44
45 p1 = pop()
46
47 h = tos
48
49 # push all of the term's factors
50
51 if (isnum(p1))
52 push(p1)
53 factor_small_number()
54 else if (car(p1) == symbol(ADD))
55 push(p1)
56 __factor_add()
57 #printf(">>>\n")
58 #for (i = h; i < tos; i++)
59 #print(stdout, stack[i])
60 #printf("<<<\n")
61 else if (car(p1) == symbol(MULTIPLY))
62 p1 = cdr(p1)
63 if (isnum(car(p1)))
64 push(car(p1))
65 factor_small_number()
66 p1 = cdr(p1)
67 while (iscons(p1))
68 p2 = car(p1)
69 if (car(p2) == symbol(POWER))
70 push(cadr(p2))
71 push(caddr(p2))
72 else
73 push(p2)
74 push(one)
75 p1 = cdr(p1)
76 else if (car(p1) == symbol(POWER))
77 push(cadr(p1))
78 push(caddr(p1))
79 else
80 push(p1)
81 push(one)
82
83 k = tos
84
85 # contruct divisors by recursive descent
86
87 push(one)
88
89 gen(h, k)
90
91 # move
92
93 n = tos - k
94
95 for i in [0...n]
96 stack[h + i] = stack[k + i]
97
98 moveTos h + n
99
100 restore()
101
102#-----------------------------------------------------------------------------
103#
104# Generate divisors
105#
106# Input: Base-exponent pairs on stack
107#
108# h first pair
109#
110# k just past last pair
111#
112# Output: Divisors on stack
113#
114# For example, factor list 2 2 3 1 results in 6 divisors,
115#
116# 1
117# 3
118# 2
119# 6
120# 4
121# 12
122#
123#-----------------------------------------------------------------------------
124
125#define ACCUM p1
126#define BASE p2
127#define EXPO p3
128
129gen = (h,k) ->
130 expo = 0
131 i = 0
132
133 save()
134
135 p1 = pop()
136
137 if (h == k)
138 push(p1)
139 restore()
140 return
141
142 p2 = stack[h + 0]
143 p3 = stack[h + 1]
144
145 push(p3)
146 expo = pop_integer()
147
148 if (!isNaN(expo))
149 for i in [0..Math.abs(expo)]
150 push(p1)
151 push(p2)
152 push_integer(sign(expo) * i)
153 power()
154 multiply()
155 gen(h + 2, k)
156
157 restore()
158
159#-----------------------------------------------------------------------------
160#
161# Factor ADD expression
162#
163# Input: Expression on stack
164#
165# Output: Factors on stack
166#
167# Each factor consists of two expressions, the factor itself followed
168# by the exponent.
169#
170#-----------------------------------------------------------------------------
171
172__factor_add = ->
173 save()
174
175 p1 = pop()
176
177 # get gcd of all terms
178
179 p3 = cdr(p1)
180 push(car(p3))
181 p3 = cdr(p3)
182 while (iscons(p3))
183 push(car(p3))
184 gcd()
185 p3 = cdr(p3)
186
187 # check gcd
188
189 p2 = pop()
190 if (isplusone(p2))
191 push(p1)
192 push(one)
193 restore()
194 return
195
196 # push factored gcd
197
198 if (isnum(p2))
199 push(p2)
200 factor_small_number()
201 else if (car(p2) == symbol(MULTIPLY))
202 p3 = cdr(p2)
203 if (isnum(car(p3)))
204 push(car(p3))
205 factor_small_number()
206 else
207 push(car(p3))
208 push(one)
209 p3 = cdr(p3)
210 while (iscons(p3))
211 push(car(p3))
212 push(one)
213 p3 = cdr(p3)
214 else
215 push(p2)
216 push(one)
217
218 # divide each term by gcd
219
220 push(p2)
221 inverse()
222 p2 = pop()
223
224 push(zero)
225 p3 = cdr(p1)
226 while (iscons(p3))
227 push(p2)
228 push(car(p3))
229 multiply()
230 add()
231 p3 = cdr(p3)
232
233 push(one)
234
235 restore()
236