-
Notifications
You must be signed in to change notification settings - Fork 6
/
Copy pathengfmt.lisp
279 lines (230 loc) · 8.52 KB
/
engfmt.lisp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
(in-package "ENGFMT")
(declaim (optimize (speed 3)
(safety 0)))
;; -----------------------------------------------------------------
(defun generate-output-to-stream (stream sign digits nsig expon exponent-printer)
(declare (type float sign)
(type list digits)
(type fixnum nsig expon))
(destructuring-bind (w &rest ds) digits
(declare (type fixnum w))
(format stream
(if (or (>= nsig 3)
(and (= nsig 2) ;; sufficient to show decimal point
(< w 100))
(and (= nsig 1) ;; ditto
(< w 10)))
"~D.~{~D~}"
"~D")
(if (minusp sign)
(- w)
w)
ds))
(funcall exponent-printer expon stream))
(defun generate-output (stream &rest args)
(cond ((eq stream t) (apply 'generate-output-to-stream *standard-output* args))
((null stream) (with-output-to-string (s)
(apply 'generate-output-to-stream s args)))
(t (apply 'generate-output-to-stream stream args))
))
;; -----------------------------------------------------------------
(defun finitep (v)
(or (zerop v)
(not (zerop (/ v)))))
(defvar *rnders*
#(1 ;; 1d0 - 1d15
10
100
1000
10000
100000
1000000
10000000
100000000
1000000000
10000000000
100000000000
1000000000000
10000000000000
100000000000000
1000000000000000))
(defun rnd-nsig (v nsig)
;; round a value according to nsig digits
;; return rounded scaled integer, scale factor, and number of dp
(declare (float v)
(fixnum nsig))
(um:bind*
((ns (the fixnum (- nsig 1 (floor (log v 10))))))
(case ns
(-2 (values (* 100 (round v 100)) 1 -2))
(-1 (values (* 10 (round v 10)) 1 -1))
(otherwise
(um:bind*
((rfact (aref *rnders* ns)))
(values (round (* rfact v)) rfact ns)))
)))
(defun decode-number (v nsig base lbase)
;; return: sign, list-of-digits, exponent
(declare (type float v)
(type fixnum nsig base lbase))
(if (zerop v)
;; sign, digits, exponent
(values 0.0 (loop for ix fixnum from 0 below nsig collect 0) 0)
(um:bind*
((:values (frac e2 sign) (decode-float v))
(:declare
(type float frac sign)
(type fixnum e2))
;; get the whole part of the exponent in specified base
;; and the fractional part for scaling the decoded base-2 fraction
(:values (e10w e10f) (floor (* e2 (log 2d0 base))))
(:declare
(fixnum e10w)
(float e10f))
;; get the fraction in the indicated base, the rounded value according to nsig,
;; and the new exponent in the indicated base.
(m10 (the float (* frac (expt base e10f))))
(e10 (the fixnum (* e10w lbase)))
(:declare
(float m10)
(fixnum e10))
(:values (rm10 rfact ns) (rnd-nsig m10 nsig)))
;; if rounded value is less than 1, then mult by base and decrease exponent
(loop while (< rm10 rfact)
do
(setf m10 (the float (* m10 base)))
(decf e10 lbase)
(um:bind*
((:values (rm rf n) (rnd-nsig m10 nsig)))
(setf rm10 rm
rfact rf
ns n)))
;; if rounded value is > base, then div by base and increase exponent
(loop while (>= rm10 (* rfact base))
do
(setf m10 (the float (/ m10 base)))
(incf e10 lbase)
(um:bind*
((:values (rm rf n) (rnd-nsig m10 nsig)))
(setf rm10 rm
rfact rf
ns n)))
(values sign
(if (plusp ns)
(um:bind*
((digits nil))
(loop for ix from 0 below ns
for (w f) = (multiple-value-list (truncate rm10 10))
then (multiple-value-list (truncate w 10))
do (push f digits)
finally (push w digits))
digits)
(list rm10))
e10)
)))
(defun check-nsig (nsig)
(assert (typep nsig '(integer 1 15))
(nsig)
"Number of significant digits (~S) must be an integer in [1,15]" nsig))
(defun format-with-computed-predigits (stream v nsig base lbase exponent-printer)
(declare (type float v)
(type fixnum nsig))
(um:bind*
((:values (sign digits expon) (decode-number v nsig base lbase))
(:declare
(type float sign)
(type list digits)
(type fixnum expon)))
(generate-output stream sign digits nsig
expon exponent-printer)
))
(defun complex-formatter (stream v formatter-fn args)
(declare (type complex v)
(ftype (function (t real &rest t) string) formatter-fn))
(apply 'format stream "#C(~A ~A)"
(mapcar (lambda (part)
(declare (type real part))
(apply formatter-fn nil part args))
(list (realpart v)
(imagpart v))
)))
;; -----------------------------------------------------------------
(defmethod engineering-format (stream (v float)
&key
(nsig 3)
(exponent-printer 'default-exponent-printer))
(cond ((finitep v)
(check-nsig nsig)
(format-with-computed-predigits stream v nsig 1000 3
exponent-printer))
(t (if (plusp v)
"+Infinity"
"-Infinity"))
))
(defmethod engineering-format (stream (v rational) &rest args &key &allow-other-keys)
(apply 'engineering-format stream (float v 1d0) args))
(defmethod engineering-format (stream (v complex) &rest args &key (nsig 3) &allow-other-keys)
(check-nsig nsig)
(complex-formatter stream v 'engineering-format args))
(defmethod engineering-format (stream v &key &allow-other-keys)
(declare (ignore v))
(format stream "NaN"))
;; -----------------------------------------------------------------
(defmethod scientific-format (stream (v float)
&key
(nsig 3)
(exponent-printer 'default-exponent-printer))
(cond ((finitep v)
(check-nsig nsig)
(format-with-computed-predigits stream v nsig 10 1
exponent-printer))
(t (if (plusp v)
"+Infinity"
"-Infinity"))
))
(defmethod scientific-format (stream (v rational) &rest args &key &allow-other-keys)
(apply 'scientific-format stream (float v 1d0) args))
(defmethod scientific-format (stream (v complex) &rest args &key (nsig 3) &allow-other-keys)
(check-nsig nsig)
(complex-formatter stream v 'scientific-format args))
(defmethod scientific-format (stream v &key &allow-other-keys)
(declare (ignore v))
(format stream "NaN"))
;; -----------------------------------------------------------------
;; -----------------------------------------------------------------
(defun default-exponent-printer (expon &optional (stream *standard-output*))
(declare (type fixnum expon))
(unless (zerop expon)
(princ #\E stream)
(princ expon stream)
))
(defun lower-case-e-exponent-printer (expon &optional (stream *standard-output*))
(declare (type fixnum expon))
(unless (zerop expon)
(princ #\e stream)
(princ expon stream)))
(defun always-signed-exponent-printer (expon &optional (stream *standard-output*))
(declare (type fixnum expon))
(unless (zerop expon)
(princ #\E stream)
(when (plusp expon)
(princ #\+ stream))
(princ expon stream)))
(defun paren-style-exponent-printer (expon &optional (stream *standard-output*))
(declare (type fixnum expon))
(unless (zerop expon)
(princ #\( stream)
(princ expon stream)
(princ #\) stream)))
(defun mathematica-style-exponent-printer (expon &optional (stream *standard-output*))
(declare (type fixnum expon))
(unless (zerop expon)
(princ "*^" stream)
(princ expon stream)))
;; -----------------------------------------------------------------
#|
(um:bind*
((x -1.23d-4)
(ndigits 8))
(engineering-format nil x :nsig ndigits))
|#