src/test/bench/shootout-nbody.rs

   1 // The Computer Language Benchmarks Game
   2 // http://benchmarksgame.alioth.debian.org/
   3 //
   4 // contributed by the Rust Project Developers
   5
   6 // Copyright (c) 2011-2014 The Rust Project Developers
   7 //
   8 // All rights reserved.
   9 //
  10 // Redistribution and use in source and binary forms, with or without
  11 // modification, are permitted provided that the following conditions
  12 // are met:
  13 //
  14 // - Redistributions of source code must retain the above copyright
  15 //   notice, this list of conditions and the following disclaimer.
  16 //
  17 // - Redistributions in binary form must reproduce the above copyright
  18 //   notice, this list of conditions and the following disclaimer in
  19 //   the documentation and/or other materials provided with the
  20 //   distribution.
  21 //
  22 // - Neither the name of "The Computer Language Benchmarks Game" nor
  23 //   the name of "The Computer Language Shootout Benchmarks" nor the
  24 //   names of its contributors may be used to endorse or promote
  25 //   products derived from this software without specific prior
  26 //   written permission.
  27 //
  28 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  29 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  30 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
  31 // FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
  32 // COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
  33 // INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
  34 // (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
  35 // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
  36 // HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
  37 // STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  38 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
  39 // OF THE POSSIBILITY OF SUCH DAMAGE.
  40
  41 use std::num::Float;
  42
  43 const PI: f64 = 3.141592653589793;
  44 const SOLAR_MASS: f64 = 4.0 * PI * PI;
  45 const YEAR: f64 = 365.24;
  46 const N_BODIES: uint = 5;
  47
  48 static BODIES: [Planet;N_BODIES] = [
  49     // Sun
  50     Planet {
  51         x: 0.0, y: 0.0, z: 0.0,
  52         vx: 0.0, vy: 0.0, vz: 0.0,
  53         mass: SOLAR_MASS,
  54     },
  55     // Jupiter
  56     Planet {
  57         x: 4.84143144246472090e+00,
  58         y: -1.16032004402742839e+00,
  59         z: -1.03622044471123109e-01,
  60         vx: 1.66007664274403694e-03 * YEAR,
  61         vy: 7.69901118419740425e-03 * YEAR,
  62         vz: -6.90460016972063023e-05 * YEAR,
  63         mass: 9.54791938424326609e-04 * SOLAR_MASS,
  64     },
  65     // Saturn
  66     Planet {
  67         x: 8.34336671824457987e+00,
  68         y: 4.12479856412430479e+00,
  69         z: -4.03523417114321381e-01,
  70         vx: -2.76742510726862411e-03 * YEAR,
  71         vy: 4.99852801234917238e-03 * YEAR,
  72         vz: 2.30417297573763929e-05 * YEAR,
  73         mass: 2.85885980666130812e-04 * SOLAR_MASS,
  74     },
  75     // Uranus
  76     Planet {
  77         x: 1.28943695621391310e+01,
  78         y: -1.51111514016986312e+01,
  79         z: -2.23307578892655734e-01,
  80         vx: 2.96460137564761618e-03 * YEAR,
  81         vy: 2.37847173959480950e-03 * YEAR,
  82         vz: -2.96589568540237556e-05 * YEAR,
  83         mass: 4.36624404335156298e-05 * SOLAR_MASS,
  84     },
  85     // Neptune
  86     Planet {
  87         x: 1.53796971148509165e+01,
  88         y: -2.59193146099879641e+01,
  89         z: 1.79258772950371181e-01,
  90         vx: 2.68067772490389322e-03 * YEAR,
  91         vy: 1.62824170038242295e-03 * YEAR,
  92         vz: -9.51592254519715870e-05 * YEAR,
  93         mass: 5.15138902046611451e-05 * SOLAR_MASS,
  94     },
  95 ];
  96
  97 struct Planet {
  98     x: f64, y: f64, z: f64,
  99     vx: f64, vy: f64, vz: f64,
 100     mass: f64,
 101 }
 102
 103 impl Copy for Planet {}
 104
 105 fn advance(bodies: &mut [Planet;N_BODIES], dt: f64, steps: int) {
 106     for _ in range(0, steps) {
 107         let mut b_slice = bodies.as_mut_slice();
 108         loop {
 109             let bi = match shift_mut_ref(&mut b_slice) {
 110                 Some(bi) => bi,
 111                 None => break
 112             };
 113             for bj in b_slice.iter_mut() {
 114                 let dx = bi.x - bj.x;
 115                 let dy = bi.y - bj.y;
 116                 let dz = bi.z - bj.z;
 117
 118                 let d2 = dx * dx + dy * dy + dz * dz;
 119                 let mag = dt / (d2 * d2.sqrt());
 120
 121                 let massj_mag = bj.mass * mag;
 122                 bi.vx -= dx * massj_mag;
 123                 bi.vy -= dy * massj_mag;
 124                 bi.vz -= dz * massj_mag;
 125
 126                 let massi_mag = bi.mass * mag;
 127                 bj.vx += dx * massi_mag;
 128                 bj.vy += dy * massi_mag;
 129                 bj.vz += dz * massi_mag;
 130             }
 131             bi.x += dt * bi.vx;
 132             bi.y += dt * bi.vy;
 133             bi.z += dt * bi.vz;
 134         }
 135     }
 136 }
 137
 138 fn energy(bodies: &[Planet;N_BODIES]) -> f64 {
 139     let mut e = 0.0;
 140     let mut bodies = bodies.iter();
 141     loop {
 142         let bi = match bodies.next() {
 143             Some(bi) => bi,
 144             None => break
 145         };
 146         e += (bi.vx * bi.vx + bi.vy * bi.vy + bi.vz * bi.vz) * bi.mass / 2.0;
 147         for bj in bodies.clone() {
 148             let dx = bi.x - bj.x;
 149             let dy = bi.y - bj.y;
 150             let dz = bi.z - bj.z;
 151             let dist = (dx * dx + dy * dy + dz * dz).sqrt();
 152             e -= bi.mass * bj.mass / dist;
 153         }
 154     }
 155     e
 156 }
 157
 158 fn offset_momentum(bodies: &mut [Planet;N_BODIES]) {
 159     let mut px = 0.0;
 160     let mut py = 0.0;
 161     let mut pz = 0.0;
 162     for bi in bodies.iter() {
 163         px += bi.vx * bi.mass;
 164         py += bi.vy * bi.mass;
 165         pz += bi.vz * bi.mass;
 166     }
 167     let sun = &mut bodies[0];
 168     sun.vx = - px / SOLAR_MASS;
 169     sun.vy = - py / SOLAR_MASS;
 170     sun.vz = - pz / SOLAR_MASS;
 171 }
 172
 173 fn main() {
 174     let n = if std::os::getenv("RUST_BENCH").is_some() {
 175         5000000
 176     } else {
 177         std::os::args().as_slice().get(1)
 178             .and_then(|arg| arg.parse())
 179             .unwrap_or(1000)
 180     };
 181     let mut bodies = BODIES;
 182
 183     offset_momentum(&mut bodies);
 184     println!("{:.9}", energy(&bodies));
 185
 186     advance(&mut bodies, 0.01, n);
 187
 188     println!("{:.9}", energy(&bodies));
 189 }
 190
 191 /// Pop a mutable reference off the head of a slice, mutating the slice to no
 192 /// longer contain the mutable reference. This is a safe operation because the
 193 /// two mutable borrows are entirely disjoint.
 194 fn shift_mut_ref<'a, T>(r: &mut &'a mut [T]) -> Option<&'a mut T> {
 195     use std::mem;
 196     use std::raw::Repr;
 197
 198     if r.len() == 0 { return None }
 199     unsafe {
 200         let mut raw = r.repr();
 201         let ret = raw.data as *mut T;
 202         raw.data = raw.data.offset(1);
 203         raw.len -= 1;
 204         *r = mem::transmute(raw);
 205         Some(unsafe { &mut *ret })
 206     }
 207 }