Hastebin

#!/usr/bin/env scheme-script
;; chezfl-colnbody.scm
;; NBODY exercise optimized in way of CL version - records splitted into field's columns.
;; compile: scheme> (parameterize ((optimize-level 3)) (compile-program "chezfl-colnbody.scm"))
;;     run: $ ./chezfl-colnbody.scm <n-iterations>

(import (chezscheme))

(module
 (take-energy)

 (define PI (* 4.0 (atan 1)))
 (define SOLAR-MASS (* 4 PI PI))
 (define DAYS-PER-YEAR 365.24)
 (define-syntax dpe (syntax-rules () ((_ v) (* v DAYS-PER-YEAR))))
 (define-syntax slr (syntax-rules () ((_ m) (* m SOLAR-MASS))))

 (define-syntax dotimes (syntax-rules () ((_ n body ...) (do ((i_ n (fx1- i_))) ((fxzero? i_)) body ...))))

 (define (take-energy n)
   (define bodies
     (flvector
      SOLAR-MASS 0.0 0.0 0.0 0.0 0.0 0.0
      (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-syntax body-size (identifier-syntax 7))
   (define bodies-len (flvector-length bodies))    ;; (2)
   (define n-body (fx/ bodies-len body-size))
   ;;(define-syntax n-body (identifier-syntax 5))
   (define-syntax dt (identifier-syntax 0.01))

   (define-syntax mkcol (syntax-rules () ((_ ix) (do ((i ix (+ i body-size))
                                                      (j 0 (+ j 1))
                                                      (column (make-flvector n-body)))
                                                     ((<= bodies-len i) column)
                                                   (flvector-set! column j (flvector-ref bodies i))))))
   ;; we rearrange bodies data into field columns
   (define-values (mass-column x-column y-column z-column veloce-x veloce-y veloce-z)
     (values (mkcol 0) (mkcol 1) (mkcol 2) (mkcol 3) (mkcol 4) (mkcol 5) (mkcol 6)))

   (define-syntax let-thval
     (syntax-rules () ((_ ((a ix col) ...) body ...)
                       (let-syntax ((a (identifier-syntax  ;; (3)
                                        (_ (flvector-ref col ix))
                                        ((set! _ val) (flvector-set! col ix val)))) ...)
                         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-syntax ^2 (syntax-rules () ((_ x) (fl* x x))))
   (define-syntax sum^2 (syntax-rules () ((_ x y z) (fl+ (^2 x) (^2 y) (^2 z)))))
   (define-syntax dsqt (syntax-rules () ((_ x) (let ((v x)) (fl* v (flsqrt v))))))

   (dotimes
    n
    (do ((i 0 (fx1+ i)))
        ((fx= i n-body))
      (do ((j (fx1+ i) (fx1+ j)))
          ((fx= j n-body))
        (let-thval
         ((massi i mass-column) (Xi i x-column) (Yi i y-column) (Zi i z-column)
          (massj j mass-column) (Xj j x-column) (Yj j y-column) (Zj j z-column)
          (VXi i veloce-x) (VYi i veloce-y) (VZi i veloce-z)
          (VXj j veloce-x) (VYj j veloce-y) (VZj j veloce-z))
         (let* ((dx (fl- Xi Xj))
                (dy (fl- Yi Yj))
                (dz (fl- Zi Zj))
                (mag (fl/ dt (dsqt (sum^2 dx dy dz)))))
           (subn! VXi (fl* dx massj mag))
           (subn! VYi (fl* dy massj mag))
           (subn! VZi (fl* dz massj mag))
           (addn! VXj (fl* dx massi mag))
           (addn! VYj (fl* dy massi mag))
           (addn! VZj (fl* dz massi mag))
         ))))
    (do ((i 0 (fx1+ i)))
        ((fx= i n-body))
      (let-thval
       ((X i x-column) (Y i y-column) (Z i z-column)
        (VX i veloce-x) (VY i veloce-y) (VZ i veloce-z))
       (addn! X (fl* dt VX))
       (addn! Y (fl* dt VY))
       (addn! Z (fl* dt VZ)))))

   (do ((i 0 (fx1+ i)) (E 0.0))
       ((fx= i n-body) E)
     (let-thval
      ((massi i mass-column)
       (VX i veloce-x)
       (VY i veloce-y)
       (VZ i veloce-z))
      (addn! E (fl* 0.5 massi (sum^2 VX VY VZ)))
      (do ((j (fx1+ i) (fx1+ j)))
          ((fx= j n-body))
        (let-thval
         ((Xi i x-column) (Yi i y-column) (Zi i z-column)
          (Xj j x-column) (Yj j y-column) (Zj j z-column) (massj j mass-column))
         (subn! E (fl/ (fl* massi massj)
                       (flsqrt (sum^2 (fl- Xi Xj)
                                      (fl- Yi Yj)
                                      (fl- Zi Zj)))))))))
   )
 ) ;; end-of-module

(define (tf t)
  (+ (time-second t) (* 1e-9 (time-nanosecond t))))

(cond ((pair? (cdr (command-line)))
       (let* ((n (string->number (cadr (command-line))))
              (st (current-time 'time-monotonic))
              (e (take-energy n))
              (et (current-time 'time-monotonic)))
         (printf "~,9f  per ~,3f sec on ~a ~%" e (tf (time-difference et st)) n)))
      (else (printf "Usage: ~a <n-iterations>~%" (car (command-line)))))


;;; i5 1035G1 1GHz  50000000 now ~ 5.4s vs chezfl-nbody 6.4s