Not so smallpt rust port...

not-so-smallpt-min.rs

#![allow(non_snake_case)] use std::{collections::HashMap as HM,f64::consts::*,io
::Write,ops::*, sync::mpsc::channel, thread::spawn,time::Instant as I};enum S{U(
f64,usize,Vec<V3>),F}#[derive(Copy, Clone)]#[repr(C)]struct V3(f64,f64,f64);impl
Add for V3{type Output=Self; fn add(self,o:Self)->Self{V3(self.0+o.0,self.1+o.1,
self.2+o.2)}}impl Sub for V3{type Output=Self;fn sub(self,o:Self)->Self{V3(self.
0-o.0,self.1-o.1,self.2-o.2)}}impl Mul<f64>for V3{type Output=Self;fn mul(self,o
:f64)->Self{V3(self.0*o,self.1*o,self.2*o)}} impl Mul for V3{type Output=Self;fn
mul(self,o:Self) -> Self{V3(self.0*o.0,self.1*o.1,self.2*o.2)}} impl Rem for V3{
type Output=Self; fn rem(self,o:Self)->Self{V3(self.1*o.2-self.2*o.1,self.2*o.0-
self.0*o.2,self.0*o.1-self.1*o.0,)}} impl V3{ fn norm(self)->Self{self*(1./self.
dot(self).sqrt())}fn dot(self,o:V3)->f64{self.0*o.0+self.1*o.1+self.2*o.2}fn map
(self,f:impl Fn(f64)->f64) -> Self{V3(f(self.0), f(self.1),f(self.2))}}#[derive(
Copy,Clone)]struct Ray(V3,V3);enum R{D,S,R} struct Sph{R:f64,p:V3,e:V3,c:V3,r:R}
fn rng(x:&mut u32)->f64{*x=*x^*x<<13;*x=*x^*x>>17;*x=*x^*x<<5;*x as f64/std::u32
::MAX as f64}const END:f64=1e20; const E:f64=1e-4;const V3Z:V3=V3(0.,0.,0.);impl
Sph{fn intersect(&self,r:&Ray)->f64{let op=self.p-r.0;let b=op.dot(r.1);let det=
b*b-op.dot(op)+self.R*self.R;let det=if det<0.{return 0.}else{det.sqrt()};let t=
b-det;if t>E{t}else{let t=b+det;if t>E{t}else{0.}}}}const SPHERES:[Sph;9]=[Sph{R
:1e5,p:V3(1e5+1.,40.8,81.6),e:V3Z,c: V3(0.75,0.25,0.25),r:R::D},Sph{R:1e5,p:V3(-
1e5+99.,40.8,81.6),e:V3Z,c: V3(0.25,0.25,0.75),r: R::D},Sph{R:1e5,p:V3(50.,40.8,
1e5),e:V3Z,c:V3(0.75,0.75,0.75),r:R::D},Sph{R:1e5,p:V3(50.,40.8,-1e5+170.),e:V3Z
,c:V3Z,r:R::D},Sph{R:1e5,p: V3(50.,1e5,81.6),e:V3Z,c:V3(0.75,0.75,0.75),r:R::D},
Sph{R:1e5,p:V3(50.,-1e5+81.6,81.6),e:V3Z,c:V3(0.75,0.75,0.75),r:R::D},Sph{R:16.5
,p:V3(27.,16.5,47.),e:V3Z,c: V3(0.999,0.999,0.999),r: R::S},Sph{R:16.5,p:V3(73.,
16.5,78.),e:V3Z,c:V3(0.999,0.999,0.999),r: R::R},Sph{R:600.,p:V3(50.,681.6-0.27,
81.6),e:V3(12.,12.,12.),c:V3Z,r:R::D},]; fn clamp(x:f64)->f64{x.max(0.).min(1.)}
fn toint(x:f64)->u8{(clamp(x).powf(1./2.2)*255.+0.5)as u8}fn intersect(r:&Ray,t:
&mut f64,sph:&mut &Sph)->bool {let x=SPHERES.iter().fold((&SPHERES[0],END),|acc,
i|{let nt=i.intersect(r);if nt>0.&&nt<acc.1{(i,nt)}else{acc}});*sph=&x.0;*t=x.1;
*t<END}fn rad(r:&Ray,d:i32,Xi:&mut u32)->V3{let(mut t,mut obj)=(0.,&SPHERES[0]);
if d>100||!intersect(r,&mut t,&mut obj){return V3Z}let x=r.0+r.1*t;let n=(x-obj.
p).norm();let nl=if n.dot(r.1)<0.{n}else{n*-1.};let mut f=obj.c;let(d,p)=(d+1,f.
0.max(f.1).max(f.2));if d>5{if rng(Xi)<p{f=f*(1./p)} else{return obj.e}}if let R
::D=obj.r{let(r1,r2)=(2.*PI*rng(Xi),rng(Xi));let r2s=r2.sqrt();let w=nl;let u=((
if w.0.abs()>0.1{V3(0., 1., 0.)} else {V3(1.,0.,0.)})%w).norm();let v=w%u;return
obj.e+f*rad(&Ray(x,(u*r1.cos()*r2s+v*r1.sin()*r2s+w*(1.-r2).sqrt()).norm()),d,Xi
)} else if let R::S=obj.r{ return obj.e+f*rad(&Ray(x,r.1-n*2.*n.dot(r.1)),d,Xi)}
let reflRay=Ray(x,r.1-n*2.*n.dot(r.1));let into=n.dot(nl)>0.;let(nc,nt)=(1.,1.5)
;let (nnt,ddn)=(if into{nc/nt}else{nt/nc},r.1.dot(nl));let cos2t=1.-nnt*nnt*(1.-
ddn*ddn); if cos2t<0.{return obj.e+f*rad(&reflRay,d,Xi)}let tdir=(r.1*nnt-n*((if
into{1.}else{-1.})*(ddn*nnt+cos2t.sqrt()))).norm();let(a,b)=(nt-nc,nt+nc);let R0
=a*a/(b*b); let c=1.-if into{-ddn}else{tdir.dot(n)};let Re=R0+(1.-R0)*c.powi(5);
let(Tr,P)=(1.-Re,0.25+0.5*Re);let(RP,TP)=(Re/P,Tr/(1.-P));obj.e+f*(if d>2{if rng
(Xi)<P{rad(&reflRay,d,Xi)*RP}else{rad(&Ray(x,tdir),d,Xi)*TP}}else{rad(&reflRay,d
,Xi)*Re + rad(&Ray(x,tdir), d,Xi)*Tr})}fn main()->std::io::Result<()>{let(w,h)=(
1024,768); let mut args = std::env::args(); let ns = args.nth(1). unwrap_or("4".
to_string()).parse().unwrap_or(4);let ns=(if ns<4{4}else{ns})/4;let nt=args.next
(). unwrap_or("4".to_string()). parse(). unwrap_or(4); let cam=Ray(V3(50.,51.55,
295.6), V3(0.,-0.042612,-1.).norm());let cx=V3(w as f64*0.5135/h as f64, 0., 0.)
;let cy=(cx%cam.1).norm()*0.5135;let mut Ts=Vec::with_capacity(nt);let mut ts=(0
..nt).map(|x|(x,0.)).collect::<HM<_,_>>();let(stx,srx)=channel();for t in 0..nt{
let stx=stx.clone();Ts.push(spawn(move||{let(mut Xi,rs)=(t as u32+1,|Xi:&mut _|{
let r1=2.*rng(Xi);if r1<1.{r1.sqrt()-1.}else{1.-(2.-r1).sqrt()}});for y in(0..(h
+nt-1)/nt).map(|x|x*nt+t).take_while(|x| *x<h){ let mut c=Vec::with_capacity(w);
for x in 0..w{let C=[(0,0),(0,1),(1,1),(1,0)].iter().fold(V3Z,|acc,(sx,sy)|{acc+
(0..ns).fold(V3Z,|acc,_|{let(dx,dy)=(rs(&mut Xi),rs(&mut Xi));let d=cx*(((*sx as
f64+0.5+dx)/2.+ x as f64)/w as f64-0.5)+ cy*(((*sy as f64+0.5+dy)/2.+y as f64)/h
as f64-0.5)+cam.1;acc+rad(&Ray(cam.0+d*140.,d.norm()),0,&mut Xi)*(1./ns as f64)}
).map(clamp)});c.push(C*0.25)}stx.send((t,S::U(y as f64/(h-1) as f64,h-y-1,c))).
unwrap()}stx.send((t,S::F)).unwrap()}));ts.insert(t, 0.);} let(mut li,mut so,mut
ib)=(I::now(),0.,Vec::new());ib.resize(w*h,V3Z); while !ts.is_empty(){match srx.
recv().unwrap(){ (t, S::U(pct, sli, sl)) => {(&mut ib [sli * w .. (sli+1) * w]).
copy_from_slice(sl.as_slice());ts.entry(t).and_modify(|v|*v=pct);}(t,S::F)=> {so
+=1.; ts.remove(&t);}}let now=I::now();if now.duration_since(li).as_millis()>=50
{eprint!("\rRendering ({} spp) {:5.2}%",ns*4,(ts.iter().map(|(_,v)|*v).sum::<f64
>()+so)/nt as f64*100.); li=now}}eprintln!(); for i in Ts{i.join().unwrap ()}let
mut f=std::fs::File::create ("image.ppm")?;let mut fb=Vec::with_capacity(w*h*3);
f.write_all(format!("P6\n{} {}\n255\n", w, h).as_bytes())?; for i in ib{fb.push(
toint(i.0));fb.push(toint(i.1));fb.push(toint(i.2))}f.write_all(fb.as_slice())?;
Ok(())}

not-so-smallpt.rs

// Based on smallpt by Kevin Beason: http://www.kevinbeason.com/smallpt/
#![allow(non_snake_case, non_camel_case_types)]
use std::{
    collections::HashMap as Map, f64::consts::*, io::Write, ops::*, sync::mpsc::channel,
    thread::spawn, time::Instant,
};
enum State {
    Update(f64, usize, Vec<V3>),
    Finished,
}
#[derive(Copy, Clone)]
#[repr(C)]
struct V3(f64, f64, f64);
impl Add for V3 {
    type Output = Self;
    fn add(self, o: Self) -> Self { V3(self.0 + o.0, self.1 + o.1, self.2 + o.2) }
}
impl Sub for V3 {
    type Output = Self;
    fn sub(self, o: Self) -> Self { V3(self.0 - o.0, self.1 - o.1, self.2 - o.2) }
}
impl Mul<f64> for V3 {
    type Output = Self;
    fn mul(self, o: f64) -> Self { V3(self.0 * o, self.1 * o, self.2 * o) }
}
impl Mul for V3 {
    type Output = Self;
    fn mul(self, o: Self) -> Self { V3(self.0 * o.0, self.1 * o.1, self.2 * o.2) }
}
impl Rem for V3 {
    type Output = Self;
    fn rem(self, o: Self) -> Self {
        V3(
            self.1 * o.2 - self.2 * o.1,
            self.2 * o.0 - self.0 * o.2,
            self.0 * o.1 - self.1 * o.0,
        )
    }
}
impl V3 {
    fn norm(self) -> Self { self * (1. / self.dot(self).sqrt()) }
    fn dot(self, o: V3) -> f64 { self.0 * o.0 + self.1 * o.1 + self.2 * o.2 }
    fn map(self, f: impl Fn(f64) -> f64) -> Self { V3(f(self.0), f(self.1), f(self.2)) }
}
#[derive(Copy, Clone)]
struct Ray(V3, V3);
enum Refl_t {
    DIFF,
    SPEC,
    REFR,
}
struct Sphere {
    rad: f64,
    p: V3,
    e: V3,
    c: V3,
    refl: Refl_t,
}
fn rng(x: &mut u32) -> f64 {
    // 32bit xorshift
    *x = *x ^ *x << 13;
    *x = *x ^ *x >> 17;
    *x = *x ^ *x << 5;
    *x as f64 / std::u32::MAX as f64
}
const END: f64 = 1e20;
const EPS: f64 = 1e-4;
const V3Z: V3 = V3(0., 0., 0.);
impl Sphere {
    fn intersect(&self, r: &Ray) -> f64 {
        let op = self.p - r.0;
        let b = op.dot(r.1);
        let det = b * b - op.dot(op) + self.rad * self.rad;
        let det = if det < 0. { return 0. } else { det.sqrt() };
        let t = b - det;
        if t > EPS {
            t
        } else {
            let t = b + det;
            if t > EPS {
                t
            } else {
                0.
            }
        }
    }
}
const WHITE: V3 = V3(0.75,  0.75,  0.75);
const BLUE: V3  = V3(0.25,  0.25,  0.75);
const RED: V3   = V3(0.75,  0.25,  0.25);
const LIGHT: V3 = V3(0.999, 0.999, 0.999);
const BRIGHT: V3= V3(12.,   12.,   12.);
const SPHERES: [Sphere; 9] = [
    Sphere { rad: 1e5,  p: V3(1e5 + 1., 40.8, 81.6),    e: V3Z,    c: RED,   refl: Refl_t::DIFF, },
    Sphere { rad: 1e5,  p: V3(-1e5 + 99., 40.8, 81.6),  e: V3Z,    c: BLUE,  refl: Refl_t::DIFF, },
    Sphere { rad: 1e5,  p: V3(50., 40.8, 1e5),          e: V3Z,    c: WHITE, refl: Refl_t::DIFF, },
    Sphere { rad: 1e5,  p: V3(50., 40.8, -1e5 + 170.),  e: V3Z,    c: V3Z,   refl: Refl_t::DIFF, },
    Sphere { rad: 1e5,  p: V3(50., 1e5, 81.6),          e: V3Z,    c: WHITE, refl: Refl_t::DIFF, },
    Sphere { rad: 1e5,  p: V3(50., -1e5 + 81.6, 81.6),  e: V3Z,    c: WHITE, refl: Refl_t::DIFF, },
    Sphere { rad: 16.5, p: V3(27., 16.5, 47.),          e: V3Z,    c: LIGHT, refl: Refl_t::SPEC, },
    Sphere { rad: 16.5, p: V3(73., 16.5, 78.),          e: V3Z,    c: LIGHT, refl: Refl_t::REFR, },
    Sphere { rad: 600., p: V3(50., 681.6 - 0.27, 81.6), e: BRIGHT, c: V3Z,   refl: Refl_t::DIFF, },
];
fn clamp(x: f64) -> f64 { x.max(0.0).min(1.0) }
fn toint(x: f64) -> u8 { (clamp(x).powf(1. / 2.2) * 255. + 0.5) as u8 }
fn intersect(r: &Ray, t: &mut f64, sph: &mut &Sphere) -> bool {
    let x = SPHERES.iter().fold((&SPHERES[0], END), |acc, i| {
        let nt = i.intersect(r);
        if nt > 0. && nt < acc.1 { (i, nt) } else { acc } });
    *sph = &x.0;
    *t = x.1;
    *t < END
}
fn radiance(r: &Ray, depth: i32, Xi: &mut u32) -> V3 {
    let (mut t, mut obj) = (0., &SPHERES[0]);
    if depth > 100 || !intersect(r, &mut t, &mut obj) {
        return V3Z;
    }
    let x = r.0 + r.1 * t;
    let n = (x - obj.p).norm();
    let nl = if n.dot(r.1) < 0. { n } else { n * -1. };
    let mut f = obj.c;
    let (depth, p) = (depth + 1, f.0.max(f.1).max(f.2));
    if depth > 5 { if rng(Xi) < p { f = f * (1. / p) } else { return obj.e } }
    if let Refl_t::DIFF = obj.refl {
        let (r1, r2) = (2. * PI * rng(Xi), rng(Xi));
        let r2s = r2.sqrt();
        let w = nl;
        let u = ((if w.0.abs() > 0.1 {
            V3(0., 1., 0.)
        } else {
            V3(1., 0., 0.)
        }) % w)
            .norm();
        let v = w % u;
        return obj.e + f * radiance(&Ray(x, (u * r1.cos() * r2s
                                           + v * r1.sin() * r2s
                                           + w * (1. - r2).sqrt()).norm()), depth, Xi)
    } else if let Refl_t::SPEC = obj.refl {
        return obj.e + f * radiance(&Ray(x, r.1 - n * 2. * n.dot(r.1)), depth, Xi)
    }
    let reflRay = Ray(x, r.1 - n * 2. * n.dot(r.1));
    let into = n.dot(nl) > 0.;
    let (nc, nt) = (1., 1.5);
    let (nnt, ddn) = (if into { nc / nt } else { nt / nc }, r.1.dot(nl));
    let cos2t = 1. - nnt * nnt * (1. - ddn * ddn);
    if cos2t < 0. { return obj.e + f * radiance(&reflRay, depth, Xi); }
    let tdir = (r.1 * nnt - n * ((if into { 1. } else { -1. }) * (ddn * nnt + cos2t.sqrt()))).norm();
    let (a, b) = (nt - nc, nt + nc);
    let R0 = a * a / (b * b);
    let c = 1. - if into { -ddn } else { tdir.dot(n) };
    let Re = R0 + (1. - R0) * c.powi(5);
    let (Tr, P) = (1. - Re, 0.25 + 0.5 * Re);
    let (RP, TP) = (Re / P, Tr / (1. - P));
    obj.e + f * (if depth > 2 { if rng(Xi) < P { radiance(&reflRay, depth, Xi) * RP }
                                else {           radiance(&Ray(x, tdir), depth, Xi) * TP } }
                 else { radiance(&reflRay, depth, Xi) * Re + radiance(&Ray(x, tdir), depth, Xi) * Tr })
}
fn main() -> std::io::Result<()> {
    let (w, h) = (1024, 768);
    let mut args = std::env::args();
    let num_samps = args.nth(1).unwrap_or("4".to_string()).parse().unwrap_or(4);
    let num_samps = (if num_samps < 4 { 4 } else { num_samps }) / 4;
    let num_threads = args.next().unwrap_or("4".to_string()).parse().unwrap_or(4);
    let cam = Ray(V3(50., 51.55 /*XXX 52*/, 295.6), V3(0., -0.042612, -1.).norm());
    let cx = V3(w as f64 * 0.5135 / h as f64, 0., 0.);
    let cy = (cx % cam.1).norm() * 0.5135;
    let mut threads = Vec::with_capacity(num_threads);
    let mut thread_states = (0..num_threads).map(|x| (x, 0.)).collect::<Map<_, _>>();
    let (status_tx, status_rx) = channel();
    for t in 0..num_threads {
        let status_tx = status_tx.clone();
        threads.push(spawn(move || {
            let (mut Xi, rs) = (t as u32 + 1, |Xi: &mut _| {
                let r1 = 2. * rng(Xi);
                if r1 < 1. { r1.sqrt() - 1.  } else { 1. - (2. - r1).sqrt()}});

            for y in (0..(h + num_threads - 1) / num_threads) // render only h/num_threads rows...
                .map(|x| x * num_threads + t)                 // render ony every num_threads'th row
                .take_while(|x| *x < h)                       // fix up larger then allowed...
            {
                let mut row = Vec::with_capacity(w);
                for x in 0..w {
                    row.push(
                        [(0, 0), (0, 1), (1, 1), (1, 0)]
                            .iter()
                            .fold(V3Z, |acc, (sx, sy)| {
                                acc + (0..num_samps)
                                    .fold(V3Z, |acc, _| {
                                        let (dx, dy) = (rs(&mut Xi), rs(&mut Xi));
                                        let d = cx * (((*sx as f64 + 0.5 + dx) / 2. + x as f64) / w as f64 - 0.5)
                                              + cy * (((*sy as f64 + 0.5 + dy) / 2. + y as f64) / h as f64 - 0.5)
                                              + cam.1;
                                        acc + radiance(&Ray(cam.0 + d * 140., d.norm()), 0, &mut Xi)
                                                  * (1. / num_samps as f64)
                                    }).map(clamp) // clamp here, to not blow out large-valued samples.
                            }) * 0.25);
                }
                status_tx.send((t, State::Update(y as f64 / (h - 1) as f64, h - y - 1, row))).unwrap();
            }
            status_tx.send((t, State::Finished)).unwrap();
        }));
        thread_states.insert(t, 0.);
    }

    let (mut last_instant, mut sum_offset, mut imgbuf) = (Instant::now(), 0., Vec::new());
    imgbuf.resize(w * h, V3Z);
    while !thread_states.is_empty() {
        match status_rx.recv().unwrap() {
            (t, State::Update(pct, sli, sl)) => {
                (&mut imgbuf[sli * w..(sli + 1) * w]).copy_from_slice(sl.as_slice());
                thread_states.entry(t).and_modify(|v| *v = pct);
            }
            (t, State::Finished) => {
                sum_offset += 1.;
                thread_states.remove(&t);
            }
        }

        let now = Instant::now();
        if now.duration_since(last_instant).as_millis() >= 50 {
            eprint!("\rRendering ({} spp) {:5.2}%", num_samps * 4,
                (thread_states.iter().map(|(_, v)| *v).sum::<f64>() + sum_offset) / num_threads as f64 * 100.);
            last_instant = now;
        }
    }
    eprintln!();

    for i in threads { i.join().unwrap() }
    let mut f = std::fs::File::create("image.ppm")?;
    let mut fb = Vec::with_capacity(w * h * 3);
    f.write_all(format!("P6\n{} {}\n255\n", w, h).as_bytes())?;
    for i in imgbuf { fb.push(toint(i.0)); fb.push(toint(i.1)); fb.push(toint(i.2)) }
    f.write_all(fb.as_slice())?;
    Ok(())
}