Annotation of /branches/vis12/src/compiler/fields/eval-poly.sml

Parent Directory | Revision Log

---

Revision 1950 - (view) (download)

1 :	jhr	1950	(* support for symbolic evaluation of polynomial to generate the coefficients for kernels *)
2 :
3 :			use "rational.sml";
4 :
5 :			datatype poly = P of Rational.rat list
6 :
7 :			local
8 :			fun add' ([], coeffs2) = coeffs2
9 :			| add' (coeffs1, []) = coeffs1
10 :			| add' (a::coeffs1, b::coeffs2) = Rational.+(a, b) :: add'(coeffs1, coeffs2)
11 :			in
12 :			fun add (P coeffs1, P coeffs2) = P(add' (coeffs1, coeffs2))
13 :
14 :			fun subtract (P coeffs1, P coeffs2) = let
15 :			fun sub ([], coeffs2) = List.map Rational.~ coeffs2
16 :			| sub (coeffs1, []) = coeffs1
17 :			| sub (a::coeffs1, b::coeffs2) = Rational.-(a, b) :: sub(coeffs1, coeffs2)
18 :			in
19 :			P(sub (coeffs1, coeffs2))
20 :			end
21 :
22 :			fun multiply (P coeffs1, P coeffs2) = let
23 :			fun mulOne (coeff, poly) = List.map (fn c => Rational.*(coeff, c)) poly
24 :			fun shift poly = Rational.zero :: poly
25 :			fun mul ([], _, acc) = acc
26 :			| mul (coeff::coeffs, poly, acc) =
27 :			mul (coeffs, shift poly, add'(acc, mulOne(coeff, poly)))
28 :			in
29 :			P(mul (coeffs1, coeffs2, []))
30 :			end
31 :
32 :			fun toString (P coeffs) = let
33 :			val d = length coeffs - 1
34 :			fun cat (x, []) = [x]
35 :			| cat (x, l) = x :: " + " :: l
36 :			fun toS (_, [], l) = String.concat(List.rev l)
37 :			| toS (d, c::r, l) = if Rational.isZero c
38 :			then toS(d+1, r, l)
39 :			else let
40 :			val pow = if (d = 0) then [""]
41 :			else if (d = 1) then ["*x"]
42 :			else ["*x^", Int.toString d]
43 :			in
44 :			toS(d+1, r, cat(concat(Rational.toString c :: pow), l))
45 :			end
46 :			in
47 :			toS (0, coeffs, [])
48 :			end
49 :			end (* local *)
50 :
51 :			local
52 :			val op + = add
53 :			val op - = subtract
54 :			val op * = multiply
55 :			fun op / (a, b) = P[Rational./(a, b)]
56 :
57 :			val x = P[Rational.zero, Rational.fromInt 1] (* x^1 *)
58 :			in
59 :
60 :			val t1 = (x + 1/1)
61 :			val c4hexic = [
62 :			69/80 + x*x*(~23/16 + x*x*(19/16 + x*(~7/12 + x*(1/16)))),
63 :			3/160 + x*(35/8 + x*(~341/32 + x*(10/1 + x*(~147/32 + x*(25/24 - x*(3/32)))))),
64 :			1539/160 + x*(~189/8 + x*(747/32 + x*(~12/1 + x*(109/32 + x*(~61/120 + x*(1/32))))))
65 :			]
66 :
67 :			end
68 :
69 :			(*
70 :			#define _C4HEXIC(x) \
71 :			(x >= 3 \
72 :			? 0 \
73 :			: (x >= 2 \
74 :			? 1539/160 + x*(~189/8 + x*(747/32 + x*(~12 + x*(109/32 + x*(~61/120 + x/32))))) \
75 :			: (x >= 1 \
76 :			? 3/160 + x*(35/8 + x*(~341/32 + x*(10 + x*(~147/32 + x*(25/24 - x*3/32))))) \
77 :			: 69/80 + x*x*(~23/16 + x*x*(19/16 + x*(~7/12 + x/16))) )))
78 :			*)

---

root@smlnj-gforge.cs.uchicago.edu	ViewVC Help
Powered by ViewVC 1.0.0