source: project/release/3/mathh/trunk/mathh.scm @ 8609

Last change on this file since 8609 was 8609, checked in by Kon Lovett, 12 years ago

Fix for MSVC

File size: 4.6 KB
Line 
1;;;; mathh.scm
2;;;; Provides access to ISO C math functions in <math.h>
3;;;; that are not defined by Chicken
4;;;; Public domain
5
6;; Issues
7;;
8;; - Windows Visual C++ 2005 deprecates
9;; hypot, j0, j1, jn, y0, y1, yn
10;; & renames
11;; _hypot, _j0, etc.
12;;
13;; - Windows does not provide
14;; log2, log1p, lgamma, tgamma. scalbn
15;;
16;; - 'gamma' is deprecated in favor of 'tgamma' but not available
17;; yet on some platforms, so we use 'gamma' for now.
18;;
19;; - Solaris log2 in sunmath.h
20
21;;; Scheme-level declarations
22(declare
23        (usual-integrations)
24        (no-procedure-checks-for-usual-bindings)
25        (unused
26          unimplemented-error)
27        (bound-to-procedure
28                bessel-j0 bessel-j1 bessel-jn
29                bessel-y0 bessel-y1 bessel-yn
30                cosh sinh tanh
31                hypot
32                gamma tgamma lgamma
33                log10 log2 log1p
34                fpmod modf ldexp scalbn frexp)
35        (export
36                bessel-j0 bessel-j1 bessel-jn
37                bessel-y0 bessel-y1 bessel-yn
38                cosh sinh tanh
39                hypot
40                gamma tgamma lgamma
41                log10 log2 make-log/base log1p
42                fpmod modf ldexp scalbn frexp))
43
44;;; Include proper C header
45#>
46#include <math.h>
47#include <float.h>
48<#
49
50;;; Unimplemented Support
51
52;from posixwin.scm
53
54(define (unimplemented-error name)
55  (error name (##core#immutable '"this function is not available on this platform")) )
56
57(define-macro (define-unimplemented name)
58  `(define (,name . _) (unimplemented-error ',name)) )
59
60(define-macro (lambda-unimplemented name)
61  `(lambda _ (unimplemented-error ',name)) )
62
63;;; Bessel functions of the 1st kind
64(define bessel-j0 (foreign-lambda double j0 double))
65(define bessel-j1 (foreign-lambda double j1 double))
66(define bessel-jn (foreign-lambda double jn int double))
67
68;;; Bessel functions of the 2nd kind
69(define bessel-y0 (foreign-lambda double y0 double))
70(define bessel-y1 (foreign-lambda double y1 double))
71(define bessel-yn (foreign-lambda double yn int double))
72
73;;; Hyperbolic functions
74(define cosh (foreign-lambda double cosh double))
75(define sinh (foreign-lambda double sinh double))
76(define tanh (foreign-lambda double tanh double))
77
78;;; Euclidean distance function
79(define hypot (foreign-lambda double hypot double double))
80
81;;; Gamma function
82(define gamma
83  (cond-expand
84    [windows
85      (lambda-unimplemented 'gamma)]
86    [linux
87      (foreign-lambda double "tgamma" double)]
88    [macosx
89      (foreign-lambda double "tgamma" double)]
90    [else
91      (foreign-lambda double "gamma" double)]))
92
93(define tgamma gamma)
94
95;;; Ln Abs Gamma function
96(define lgamma
97  (cond-expand
98    [windows
99      (lambda-unimplemented 'lgamma)]
100    [else
101      (foreign-lambda double lgamma double)]))
102
103;;; Base-10 logarithm
104(define log10 (foreign-lambda double log10 double))
105
106;;; Base-2 logarithm
107(define log2
108  (cond-expand
109    [(or linux macosx bsd)
110      (foreign-lambda double "log2" double)]
111    [else
112      (foreign-lambda* double ((double x))
113       "#ifndef M_LN2\n"
114       "#define #define M_LN2 0.693147180559945309417232121458176568\n"
115       "#endif\n"
116       "return (log(x) / M_LN2);")]))
117
118;; Log function for base n
119(define (make-log/base b)
120        (when (fixnum? b)
121                (set! b (exact->inexact b)))
122        (unless (flonum? b)
123                (error 'make-log/base "invalid flonum" b))
124        (cond
125                [(= 2.0 b)
126                        log2]
127                [(= 10.0 b)
128                        log10]
129                [else
130                        (let ([lnb (log b)]
131                                          [logf
132                                            (foreign-lambda* double ((double x) (double lnb))
133                                              "return (log(x) / lnb);")])
134                                (lambda (n)
135                                        (logf n lnb)))]) )
136
137;;; Natural logarithm of 1+x accurate for very small x
138(define log1p
139  (cond-expand
140    [windows ; potentially inaccurate but ...
141      (foreign-lambda* double ((double x))
142        "return (log(1.0 + x));")]
143    [else
144                  (foreign-lambda double log1p double)]))
145
146;;; Flonum remainder
147(define fpmod (foreign-lambda double fmod double double))
148
149;;; Return integer, fraction (as multiple values) of a flonum
150(define modf (foreign-primitive ((double x))
151        "double ipart;
152         double result  = modf(x, &ipart);
153         C_word* values = C_alloc(2 * C_SIZEOF_FLONUM);
154         C_word value1  = C_flonum(&values, ipart);
155         C_word value2  = C_flonum(&values, result);
156         C_values(4, C_SCHEME_UNDEFINED, C_k, value1, value2);
157"))
158
159;;; Compute x * 2**n
160(define ldexp (foreign-lambda double ldexp double integer))
161
162;;; Efficiently compute x * 2**n
163(define scalbn
164  (cond-expand
165    [windows ; not efficient but ...
166      (foreign-lambda* double ((double x) (integer n))
167        "return (ldexp(x, n));")]
168    [else
169      (foreign-lambda double scalbn double integer)]) )
170
171;;; Return mantissa, exponent (as multiple values) of a flonum
172(define frexp (foreign-primitive ((double x))
173        "int exp;
174         double result = frexp(x, &exp);
175         C_word* values = C_alloc(C_SIZEOF_FLONUM);
176         C_word value1 = C_flonum(&values, result);
177         C_word value2 = C_fix(exp);
178         C_values(4, C_SCHEME_UNDEFINED, C_k, value1, value2);
179"))
Note: See TracBrowser for help on using the repository browser.