src/test/bench/shootout-fasta-redux.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) 2013-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 #![feature(slicing_syntax)]
  42
  43 use std::cmp::min;
  44 use std::io::{stdout, IoResult};
  45 use std::os;
  46 use std::slice::bytes::copy_memory;
  47
  48 const LINE_LEN: uint = 60;
  49 const LOOKUP_SIZE: uint = 4 * 1024;
  50 const LOOKUP_SCALE: f32 = (LOOKUP_SIZE - 1) as f32;
  51
  52 // Random number generator constants
  53 const IM: u32 = 139968;
  54 const IA: u32 = 3877;
  55 const IC: u32 = 29573;
  56
  57 const ALU: &'static str = "GGCCGGGCGCGGTGGCTCACGCCTGTAATCCCAGCACTTTG\
  58                             GGAGGCCGAGGCGGGCGGATCACCTGAGGTCAGGAGTTCGA\
  59                             GACCAGCCTGGCCAACATGGTGAAACCCCGTCTCTACTAAA\
  60                             AATACAAAAATTAGCCGGGCGTGGTGGCGCGCGCCTGTAAT\
  61                             CCCAGCTACTCGGGAGGCTGAGGCAGGAGAATCGCTTGAAC\
  62                             CCGGGAGGCGGAGGTTGCAGTGAGCCGAGATCGCGCCACTG\
  63                             CACTCCAGCCTGGGCGACAGAGCGAGACTCCGTCTCAAAAA";
  64
  65 const NULL_AMINO_ACID: AminoAcid = AminoAcid { c: ' ' as u8, p: 0.0 };
  66
  67 static IUB: [AminoAcid, ..15] = [
  68     AminoAcid { c: 'a' as u8, p: 0.27 },
  69     AminoAcid { c: 'c' as u8, p: 0.12 },
  70     AminoAcid { c: 'g' as u8, p: 0.12 },
  71     AminoAcid { c: 't' as u8, p: 0.27 },
  72     AminoAcid { c: 'B' as u8, p: 0.02 },
  73     AminoAcid { c: 'D' as u8, p: 0.02 },
  74     AminoAcid { c: 'H' as u8, p: 0.02 },
  75     AminoAcid { c: 'K' as u8, p: 0.02 },
  76     AminoAcid { c: 'M' as u8, p: 0.02 },
  77     AminoAcid { c: 'N' as u8, p: 0.02 },
  78     AminoAcid { c: 'R' as u8, p: 0.02 },
  79     AminoAcid { c: 'S' as u8, p: 0.02 },
  80     AminoAcid { c: 'V' as u8, p: 0.02 },
  81     AminoAcid { c: 'W' as u8, p: 0.02 },
  82     AminoAcid { c: 'Y' as u8, p: 0.02 },
  83 ];
  84
  85 static HOMO_SAPIENS: [AminoAcid, ..4] = [
  86     AminoAcid { c: 'a' as u8, p: 0.3029549426680 },
  87     AminoAcid { c: 'c' as u8, p: 0.1979883004921 },
  88     AminoAcid { c: 'g' as u8, p: 0.1975473066391 },
  89     AminoAcid { c: 't' as u8, p: 0.3015094502008 },
  90 ];
  91
  92 // FIXME: Use map().
  93 fn sum_and_scale(a: &'static [AminoAcid]) -> Vec<AminoAcid> {
  94     let mut result = Vec::new();
  95     let mut p = 0f32;
  96     for a_i in a.iter() {
  97         let mut a_i = *a_i;
  98         p += a_i.p;
  99         a_i.p = p * LOOKUP_SCALE;
 100         result.push(a_i);
 101     }
 102     let result_len = result.len();
 103     result[result_len - 1].p = LOOKUP_SCALE;
 104     result
 105 }
 106
 107 struct AminoAcid {
 108     c: u8,
 109     p: f32,
 110 }
 111
 112 struct RepeatFasta<'a, W:'a> {
 113     alu: &'static str,
 114     out: &'a mut W
 115 }
 116
 117 impl<'a, W: Writer> RepeatFasta<'a, W> {
 118     fn new(alu: &'static str, w: &'a mut W) -> RepeatFasta<'a, W> {
 119         RepeatFasta { alu: alu, out: w }
 120     }
 121
 122     fn make(&mut self, n: uint) -> IoResult<()> {
 123         let alu_len = self.alu.len();
 124         let mut buf = Vec::from_elem(alu_len + LINE_LEN, 0u8);
 125         let alu: &[u8] = self.alu.as_bytes();
 126
 127         copy_memory(buf.as_mut_slice(), alu);
 128         let buf_len = buf.len();
 129         copy_memory(buf[mut alu_len..buf_len],
 130                     alu[..LINE_LEN]);
 131
 132         let mut pos = 0;
 133         let mut bytes;
 134         let mut n = n;
 135         while n > 0 {
 136             bytes = min(LINE_LEN, n);
 137             try!(self.out.write(buf.slice(pos, pos + bytes)));
 138             try!(self.out.write_u8('\n' as u8));
 139             pos += bytes;
 140             if pos > alu_len {
 141                 pos -= alu_len;
 142             }
 143             n -= bytes;
 144         }
 145         Ok(())
 146     }
 147 }
 148
 149 fn make_lookup(a: &[AminoAcid]) -> [AminoAcid, ..LOOKUP_SIZE] {
 150     let mut lookup = [ NULL_AMINO_ACID, ..LOOKUP_SIZE ];
 151     let mut j = 0;
 152     for (i, slot) in lookup.iter_mut().enumerate() {
 153         while a[j].p < (i as f32) {
 154             j += 1;
 155         }
 156         *slot = a[j];
 157     }
 158     lookup
 159 }
 160
 161 struct RandomFasta<'a, W:'a> {
 162     seed: u32,
 163     lookup: [AminoAcid, ..LOOKUP_SIZE],
 164     out: &'a mut W,
 165 }
 166
 167 impl<'a, W: Writer> RandomFasta<'a, W> {
 168     fn new(w: &'a mut W, a: &[AminoAcid]) -> RandomFasta<'a, W> {
 169         RandomFasta {
 170             seed: 42,
 171             out: w,
 172             lookup: make_lookup(a),
 173         }
 174     }
 175
 176     fn rng(&mut self, max: f32) -> f32 {
 177         self.seed = (self.seed * IA + IC) % IM;
 178         max * (self.seed as f32) / (IM as f32)
 179     }
 180
 181     fn nextc(&mut self) -> u8 {
 182         let r = self.rng(1.0);
 183         for a in self.lookup.iter() {
 184             if a.p >= r {
 185                 return a.c;
 186             }
 187         }
 188         0
 189     }
 190
 191     fn make(&mut self, n: uint) -> IoResult<()> {
 192         let lines = n / LINE_LEN;
 193         let chars_left = n % LINE_LEN;
 194         let mut buf = [0, ..LINE_LEN + 1];
 195
 196         for _ in range(0, lines) {
 197             for i in range(0u, LINE_LEN) {
 198                 buf[i] = self.nextc();
 199             }
 200             buf[LINE_LEN] = '\n' as u8;
 201             try!(self.out.write(&buf));
 202         }
 203         for i in range(0u, chars_left) {
 204             buf[i] = self.nextc();
 205         }
 206         self.out.write(buf[..chars_left])
 207     }
 208 }
 209
 210 fn main() {
 211     let args = os::args();
 212     let args = args.as_slice();
 213     let n = if args.len() > 1 {
 214         from_str::<uint>(args[1].as_slice()).unwrap()
 215     } else {
 216         5
 217     };
 218
 219     let mut out = stdout();
 220
 221     out.write_line(">ONE Homo sapiens alu").unwrap();
 222     {
 223         let mut repeat = RepeatFasta::new(ALU, &mut out);
 224         repeat.make(n * 2).unwrap();
 225     }
 226
 227     out.write_line(">TWO IUB ambiguity codes").unwrap();
 228     let iub = sum_and_scale(&IUB);
 229     let mut random = RandomFasta::new(&mut out, iub.as_slice());
 230     random.make(n * 3).unwrap();
 231
 232     random.out.write_line(">THREE Homo sapiens frequency").unwrap();
 233     let homo_sapiens = sum_and_scale(&HOMO_SAPIENS);
 234     random.lookup = make_lookup(homo_sapiens.as_slice());
 235     random.make(n * 5).unwrap();
 236
 237     random.out.write_str("\n").unwrap();
 238 }