diff --git a/benchmarks/nbody.el b/benchmarks/nbody.el new file mode 100644 index 0000000..7b1b2eb --- /dev/null +++ b/benchmarks/nbody.el @@ -0,0 +1,142 @@ +;; -*- lexical-binding: t; -*- + +;; Copyright (C) 2019 Free Software Foundation, Inc. + +;; This file is part of GNU Emacs. + +;; GNU Emacs is free software: you can redistribute it and/or modify +;; it under the terms of the GNU General Public License as published by +;; the Free Software Foundation, either version 3 of the License, or +;; (at your option) any later version. + +;; GNU Emacs is distributed in the hope that it will be useful, +;; but WITHOUT ANY WARRANTY; without even the implied warranty of +;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +;; GNU General Public License for more details. + +;; You should have received a copy of the GNU General Public License +;; along with GNU Emacs. If not, see . + +;;; Commentary: + +;; Adapted to elisp from CL version from: +;; https://benchmarksgame-team.pages.debian.net/benchmarksgame/program/nbody-sbcl-2.html + +(require 'cl-lib) + +(defconst elb-days-per-year 365.24) +(defconst elb-solar-mass (* 4 pi pi)) + +(cl-defstruct (elb-body (:type vector) + (:conc-name nil) + (:constructor make-elb-body (x y z vx vy vz mass))) + x y z + vx vy vz + mass) + +(defvar elb-jupiter + (make-elb-body 4.84143144246472090 + -1.16032004402742839 + -1.03622044471123109e-1 + (* 1.66007664274403694e-3 elb-days-per-year) + (* 7.69901118419740425e-3 elb-days-per-year) + (* -6.90460016972063023e-5 elb-days-per-year) + (* 9.54791938424326609e-4 elb-solar-mass))) + +(defvar elb-saturn + (make-elb-body 8.34336671824457987 + 4.12479856412430479 + -4.03523417114321381e-1 + (* -2.76742510726862411e-3 elb-days-per-year) + (* 4.99852801234917238e-3 elb-days-per-year) + (* 2.30417297573763929e-5 elb-days-per-year) + (* 2.85885980666130812e-4 elb-solar-mass))) + +(defvar elb-uranus + (make-elb-body 1.28943695621391310e1 + -1.51111514016986312e1 + -2.23307578892655734e-1 + (* 2.96460137564761618e-03 elb-days-per-year) + (* 2.37847173959480950e-03 elb-days-per-year) + (* -2.96589568540237556e-05 elb-days-per-year) + (* 4.36624404335156298e-05 elb-solar-mass))) + +(defvar elb-neptune + (make-elb-body 1.53796971148509165e+01 + -2.59193146099879641e+01 + 1.79258772950371181e-01 + (* 2.68067772490389322e-03 elb-days-per-year) + (* 1.62824170038242295e-03 elb-days-per-year) + (* -9.51592254519715870e-05 elb-days-per-year) + (* 5.15138902046611451e-05 elb-solar-mass))) + +(defvar elb-sun + (make-elb-body 0.0 0.0 0.0 0.0 0.0 0.0 elb-solar-mass)) + +(defun elb-applyforces (a b dt) + (let* ((dx (- (x a) (x b))) + (dy (- (y a) (y b))) + (dz (- (z a) (z b))) + (distance (sqrt (+ (* dx dx) (* dy dy) (* dz dz)))) + (mag (/ dt (* distance distance distance))) + (dxmag (* dx mag)) + (dymag (* dy mag)) + (dzmag (* dz mag))) + (cl-decf (vx a) (* dxmag (mass b))) + (cl-decf (vy a) (* dymag (mass b))) + (cl-decf (vz a) (* dzmag (mass b))) + (cl-incf (vx b) (* dxmag (mass a))) + (cl-incf (vy b) (* dymag (mass a))) + (cl-incf (vz b) (* dzmag (mass a)))) + nil) + +(defun elb-advance (system dt) + (cl-loop for (a . rest) on system + do (dolist (b rest) + (elb-applyforces a b dt))) + (dolist (b system) + (cl-incf (x b) (* dt (vx b))) + (cl-incf (y b) (* dt (vy b))) + (cl-incf (z b) (* dt (vz b)))) + nil) + +(defun elb-energy (system) + (let ((e 0.0)) + (cl-loop for (a . rest) on system do + (cl-incf e (* 0.5 + (mass a) + (+ (* (vx a) (vx a)) + (* (vy a) (vy a)) + (* (vz a) (vz a))))) + (dolist (b rest) + (let* ((dx (- (x a) (x b))) + (dy (- (y a) (y b))) + (dz (- (z a) (z b))) + (dist (sqrt (+ (* dx dx) (* dy dy) (* dz dz))))) + (cl-decf e (/ (* (mass a) (mass b)) dist))))) + e)) + +(defun elb-offset-momentum (system) + (let ((px 0.0) + (py 0.0) + (pz 0.0)) + (dolist (p system) + (cl-incf px (* (vx p) (mass p))) + (cl-incf py (* (vy p) (mass p))) + (cl-incf pz (* (vz p) (mass p)))) + (setf (vx (car system)) (/ (- px) elb-solar-mass) + (vy (car system)) (/ (- py) elb-solar-mass) + (vz (car system)) (/ (- pz) elb-solar-mass)) + nil)) + +(defun elb-nbody (n) + (let ((system (list elb-sun elb-jupiter elb-saturn elb-uranus elb-neptune))) + (elb-offset-momentum system) + (dotimes (_ n) + (elb-advance system 0.01)) + (elb-energy system))) + +(defun elb-nbody-entry () + (elb-nbody 300000)) + +(provide 'elb-nbody)