1
#!/usr/bin/env gsi-script
;; gsc-nbody.scm
;; NBODY exercise. Gambit Scheme version
;; compile: $ gsc -exe gsc-nbody
;; run: $ ./gsc-nbody <n-iterations>
;;
(declare (standard-bindings)
#;(extended-bindings)
(not safe)
#;(not interrupts-enabled)
)
(define PI (* 4.0 (atan 1)))
(define SOLAR-MASS (* 4 PI PI))
(define DAYS-PER-YEAR 365.24)
(define-syntax dpe (syntax-rules () ((_ v) (fl* v DAYS-PER-YEAR))))
(define-syntax slr (syntax-rules () ((_ m) (fl* m SOLAR-MASS))))
(define-syntax dotimes (syntax-rules () ((_ n body ...) (do ((i_ n (fx- i_ 1))) ((fxzero? i_)) body ...))))
(define-syntax subn! (syntax-rules () ((_ what withs) (set! what (fl- what withs)))))
(define-syntax addn! (syntax-rules () ((_ what withs) (set! what (fl+ what withs)))))
(define (take-energy n)
(define bodies
(f64vector
SOLAR-MASS 0.0 0.0 0.0 0.0 0.0 0.0 ;; SUN
(slr 9.54791938424326609e-04) 4.84143144246472090e+00 -1.16032004402742839e+00 -1.03622044471123109e-01 (dpe 1.66007664274403694e-03) (dpe 7.69901118419740425e-03) (dpe -6.90460016972063023e-05)
(slr 2.85885980666130812e-04) 8.34336671824457987e+00 4.12479856412430479e+00 -4.03523417114321381e-01 (dpe -2.76742510726862411e-03) (dpe 4.99852801234917238e-03) (dpe 2.30417297573763929e-05)
(slr 4.36624404335156298e-05) 1.28943695621391310e+01 -1.51111514016986312e+01 -2.23307578892655734e-01 (dpe 2.96460137564761618e-03) (dpe 2.37847173959480950e-03) (dpe -2.96589568540237556e-05)
(slr 5.15138902046611451e-05) 1.53796971148509165e+01 -2.59193146099879641e+01 1.79258772950371181e-01 (dpe 2.68067772490389322e-03) (dpe 1.62824170038242295e-03) (dpe -9.51592254519715870e-05)
))
(define body-size 7)
(define bodies-len (f64vector-length bodies)) ;; (2)
(define dt 0.01)
(define-syntax let-attr
(syntax-rules () ((_ ((a ix offset) ...) body ...)
(let ((a (f64vector-ref bodies (fx+ ix offset))) ...)
body ...))))
(define (subnv! ix offset val deltaval) (f64vector-set! bodies (fx+ ix offset) (fl- val deltaval)))
(define (addnv! ix offset val deltaval) (f64vector-set! bodies (fx+ ix offset) (fl+ val deltaval)))
(define-syntax ^2 (syntax-rules () ((_ x) (flsquare x))))
(define-syntax sum^2 (syntax-rules () ((_ x y z) (fl+ (^2 x) (^2 y) (^2 z)))))
(define (dsqt v) (fl* v (flsqrt v)))
(define E 0.0)
(dotimes
n
(do ((i 0 (fx+ i body-size)))
((fx= i bodies-len))
(do ((j (fx+ i body-size) (fx+ j body-size)))
((fx= j bodies-len))
(let-attr
((massi i 0) (Xi i 1) (Yi i 2) (Zi i 3)
(massj j 0) (Xj j 1) (Yj j 2) (Zj j 3)
(VXi i 4) (VYi i 5) (VZi i 6)
(VXj j 4) (VYj j 5) (VZj j 6))
(let* ((dx (fl- Xi Xj))
(dy (fl- Yi Yj))
(dz (fl- Zi Zj))
(mag (fl/ dt (dsqt (sum^2 dx dy dz)))))
(subnv! i 4 VXi (fl* dx massj mag))
(subnv! i 5 VYi (fl* dy massj mag))
(subnv! i 6 VZi (fl* dz massj mag))
(addnv! j 4 VXj (fl* dx massi mag))
(addnv! j 5 VYj (fl* dy massi mag))
(addnv! j 6 VZj (fl* dz massi mag))
))))
(do ((i 0 (fx+ i body-size)))
((fx= i bodies-len))
(let-attr
((X i 1) (Y i 2) (Z i 3)
(VX i 4) (VY i 5) (VZ i 6))
(addnv! i 1 X (fl* dt VX))
(addnv! i 2 Y (fl* dt VY))
(addnv! i 3 Z (fl* dt VZ)))))
(do ((i 0 (fx+ i body-size)))
((fx= i bodies-len) E)
(let-attr
((massi i 0)
(VX i 4)
(VY i 5)
(VZ i 6))
(addn! E (fl* 0.5 massi (sum^2 VX VY VZ)))
(do ((j (fx+ i body-size) (fx+ body-size j)))
((fx= j bodies-len))
(let-attr
((Xi i 1) (Yi i 2) (Zi i 3)
(Xj j 1) (Yj j 2) (Zj j 3) (massj j 0))
(subn! E (fl/ (fl* massi massj)
(flsqrt (sum^2
(fl- Xi Xj)
(fl- Yi Yj)
(fl- Zi Zj))))))
)))
)
(define (main . args)
(cond ((pair? args)
(let* ((n (string->number (car args)))
(st (current-second))
(e (take-energy n))
(et (current-second)))
(print e " per " (- et st) " sec on " n "\n")))
(else (print "Usage: gsc-nbody <n-iterations>~%"))))
For immediate assistance, please email our customer support: [email protected]