aboutsummaryrefslogtreecommitdiff
path: root/src/volume-bodies.sc
blob: f4ee2a3079ee87e719119a3ff3314b6772401842 (plain)
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
(sc-include "macros/macros")

(pre-include "stdio.h")
(pre-include "stdlib.h")
(pre-include "extent-sampling.h")
(pre-include "oracles.h")

(pre-define EDGE 1.0)
(pre-define RTOL 0.1)
(pre-define ELLIPSOID-AXES-START 0.5)
(pre-define ELLIPSOID-AXES-END 1.0)

(sc-define-syntax* (with-cube-oracle edge body ...)
  (let ((params (sc-gensym)))
    `(begin
       (declare ,params (struct-variable cube-params
                                         (edge ,edge)
                                         (center NULL)))
       (declare oracle (struct-variable extent-oracle-t
                                        (oracle cube-extent-oracle)
                                        (params (address-of ,params))))
       ,@body)))

(sc-define-syntax* (with-ellipsoid-oracle axes-start axes-end dimension body ...)
  (let ((params (sc-gensym))
        (axes (sc-gensym)))
    `(with-vector ,axes ,dimension
       (for-i i ,dimension
         (gsl-vector-set ,axes i (+ ,axes-start (* i (/ (- ,axes-end ,axes-start)
                                                        (- ,dimension 1))))))
       (declare ,params (struct-variable ellipsoid-params
                                         (axes ,axes)))
       (declare oracle (struct-variable extent-oracle-t
                                        (oracle ellipsoid-extent-oracle)
                                        (params (address-of ,params))))
       ,@body)))

(sc-define-syntax* (experiment-dimension dimension true-volume filename)
  `(begin
     (with-rstats samples
       (with-rstats stats
         (volume (address-of oracle)
                 (,true-volume dimension (struct-get oracle params))
                 rng
                 dimension
                 RTOL
                 stats)
         (gsl-rstat-add (gsl-rstat-n stats)
                        samples))
       (with-data-file fp ,filename "a"
         (fprintf fp "%d\t%g\n" dimension (gsl-rstat-mean samples))))))

(sc-define-syntax* (cube-experiment)
  (let* ((filename "volume-cube.dat"))
    `(begin
       (with-data-file fp ,filename "w"
         (fprintf fp "dimension\tsamples\n"))
       (for-i-step dimension 10 100 10
         (with-cube-oracle EDGE (experiment-dimension dimension cube-true-volume ,filename))))))

(sc-define-syntax* (ellipsoid-experiment)
  (let* ((filename "volume-ellipsoid.dat"))
    `(begin
       (with-data-file fp ,filename "w"
         (fprintf fp "dimension\tsamples\n"))
       (for-i-step dimension 10 100 10
         (with-ellipsoid-oracle ELLIPSOID-AXES-START ELLIPSOID-AXES-END dimension
                                (experiment-dimension dimension ellipsoid-true-volume ,filename))))))

(define (main) (int)
  (with-rng rng
    (cube-experiment)
    (ellipsoid-experiment)))