Re: Mathematica 7 is now available
- To: mathgroup at smc.vnet.net
- Subject: [mg94326] Re: Mathematica 7 is now available
- From: Jon Harrop <jon at ffconsultancy.com>
- Date: Wed, 10 Dec 2008 04:49:31 -0500 (EST)
- References: <gg62p3$56g$1@smc.vnet.net> <gh8hkr$r1d$1@smc.vnet.net> <ghavsp$ock$1@smc.vnet.net> <200812081121.GAA15802@smc.vnet.net> <ghlmbp$jof$1@smc.vnet.net>
Daniel Lichtblau wrote: > Jon Harrop wrote: >> This is not specific to OCaml. Here are equivalent implementations of >> that program in C++, Java, OCaml, SML, Scheme, Lisp and Haskell and all >> well over 100,000x faster than the Mathematica I gave (and many are more >> concise as well!): >> >> http://www.ffconsultancy.com/languages/ray_tracer/ > > I do not understand this claim that the speed difference is 5 orders of > magnitude. Some experiments indicate 3, perhaps. The code I gave is over 5 orders of magnitude slower by my measurements. > Some modest improvements to the original Mathematica code, though, bring > it to maybe 2 orders of magnitude of the OCaml speed. That is a huge improvement. Thank you. > In the variant > below, I make sure everyone is a machine real number or machine integer, > whichever is appropriate. I use Compile on RaySphere. There are a few > other tweaks for general efficiency (some might apply to the OCaml > implementation, I'm not sure). I have a few questions about this. :-) > delta = Sqrt[$MachineEpsilon]; > inf = 100000.; This is nasty. Why is infinity not represented internally by the floating point number infinity in Mathematica? > neglight = Normalize[{1., 3., -2.}]; > nohit = {inf, {0., 0., 0.}}; > > Clear[RaySphere, Intersect, Create, RayTrace, Main]; > > RaySphere = With[{inf = inf}, Compile[ What is the purpose of "With[{inf = inf}, ..]" here? > {{c, _Real, 1}, {o, _Real, 1}, {d, _Real, 1}, {r, _Real}}, > Module[{v = c - o, b, disc, t}, > b = v.d; > disc = b^2 - v.v + r; > If[disc < 0., inf, > disc = Sqrt[disc]; > If[b > disc, b - disc, > t = b + disc; > If[t < 0., inf, t]]]]]]; > > Intersect[o_, d_, ln : {lambda_, n_}, List[c_, r_]] := Block[ > {lambda2 = RaySphere[c, o, d, r]}, > If[lambda2 >= lambda, ln, {lambda2, c}]] Defer calculation of the normal. Fair enough but it is only done once per pixel. > Intersect[o_, d_, ln : {lambda_, n_}, List[c_, r_, s_]] := Block[ > {lambda2 = RaySphere[c, o, d, r]}, > If[lambda2 >= lambda, ln, > Fold[Intersect[o, d, #1, #2] &, ln, s] > ] > ] No currying and a named subpattern, ok. > RayTrace[o_, d_, scene_] := Block[ > {lambda, n, g, p}, > {lambda, n} = Intersect[o, d, nohit, scene]; > If[lambda >= inf, 0., > n = Normalize[o + d*lambda - n]; > g = n.neglight; > If[g <= 0., 0., > {lambda, n} = > Intersect[o + lambda d + delta n, neglight, nohit, scene]; > If[lambda < inf, 0., g]]]] > > Create[level_, c_, r_] := > Block[{obj = List[c, r^2]}, > If[level == 1, obj, > Block[{a = 3*r/Sqrt[12], Aux}, > Aux[x1_, z1_] := Create[level - 1, c + {x1, a, z1}, 0.5 r]; > List[c, > 9*r^2, {obj, Aux[-a, -a], Aux[a, -a], Aux[-a, a], Aux[a, a]}]]]] > > prepareTable = Compile[{{n, _Integer}, {s, _Integer}}, > With[{nss = N[s]}, > Flatten[ > Table[Normalize[{(x + s1/nss)/n - .5, (y + s2/nss)/n - .5, > 1.}], {s1, 0, s - 1}, {s2, 0, s - 1}, {y, 0, n - 1}, {x, 0, > n - 1}], 1] > ]]; So you precompute a matrix of all the ray directions using a compiled function. > Main[level_, n_, ss_] := > Block[{scene = Create[level, {0., -1., 4.}, 1.], nss = N[ss]}, > scene = MapAll[Developer`ToPackedArray, scene]; > Total[Map[RayTrace[{0., 0., 0.}, #, scene] &, > prepareTable[n, ss], {3}], 1]/nss^2] > > In[22]:= Timing[ff512 = Main[9, 512, 4];] > Out[22]= {1542.27, Null} That is orders of magnitude faster than my version but I am not sure why. I suspect the main performance improvement comes from the use of packed arrays but I haven't been able to reproduce it in my own code yet. > At this point, it might make sense to use ParallelMap instead of Map, in > Main. Replacing Map with ParallelMap slows it down, even with Method -> "CoarsestGrained". How can I get a speedup from my 8 cores?! :-) Many thanks, -- Dr Jon D Harrop, Flying Frog Consultancy Ltd. http://www.ffconsultancy.com/?u
- Follow-Ups:
- Re: Mathematica 7 is now available
- From: Daniel Lichtblau <danl@wolfram.com>
- Re: Mathematica 7 is now available
- From: Daniel Lichtblau <danl@wolfram.com>
- Re: Mathematica 7 is now available
- References:
- Re: Mathematica 7 is now available
- From: Jon Harrop <jon@ffconsultancy.com>
- Re: Mathematica 7 is now available