Reply to thread

Message: [QUOTE="Xah Lee, post: 3692612"] alright, here's my improved code, pasted near the bottom. let me say a few things about Jon's code. If we rate that piece of mathematica code on the level of: Beginner Mathematica programer, Intermediate, Advanced, where Beginner is someone who just learned tried to program Mathematica no more than 6 months, then that piece of code is Beginner level. Here's some basic analysis and explanation. The program has these main functions: â€¢ RaySphere â€¢ Intersect â€¢ RayTrace â€¢ Create â€¢ Main The Main calls Create then feed it to RayTrace. Create calls itself recursively, and basically returns a long list of a repeating element, each of the element differ in their parameter. RayTrace calls Intersect 2 times. Intersect has 2 forms, one of them calls itself recursively. Both forms calls RaySphere once. So, the core loop is with the Intersect function and RaySphere. Some 99.99% of time are spent there. ------------------ I didn't realize until after a hour, that if Jon simply give numerical arguments to Main and Create, the result timing by a factor of 0.3 of original. What a incredible sloppiness! and he intended this to show Mathematica speed with this code? The Main[] function calls Create. The create has 3 parameters: level, c, and r. The level is a integer for the recursive level of raytracing . The c is a vector for sphere center i presume. The r is radius of the sphere. His input has c and r as integers, and this in Mathematica means computation with exact arithmetics (and automatic kicks into infinite precision if necessary). Changing c and r to float immediately reduced the timing to 0.3 of original. ------------------ now, back to the core loop. The RaySphere function contain codes that does symbolic computation by calling Im, which is the imaginary part of a complex number!! and if so, it returns the symbol Infinity! The possible result of Infinity is significant because it is used in Intersect to do a numerical comparison in a If statement. So, here in these deep loops, Mathematica's symbolic computation is used for numerical purposes! So, first optimization at the superficial code form level is to get rid of this symbolic computation. Instead of checking whethere his â€œdisc = Sqrt[b^2 - v.v + r^2]â€ has imaginary part, one simply check whether the argument to sqrt is negative. after getting rid of the symbolic computation, i made the RaySphere function to be a Compiled function. I stopped my optimization at this step. The above are some _fundamental_ things any dummy who claims to code Mathematica for speed should know. Jon has written a time series Mathematica package that he's selling commercially. So, either he got very sloppy with this Mathematica code, or he intentionally made it look bad, or that his Mathematica skill is truely beginner level. Yet he dares to talk bullshit in this thread. Besides the above basic things, there are several aspects that his code can improve in speed. For example, he used pattern matching to do core loops. e.g. Intersect[o_, d_][{lambda_, n_}, Bound[c_, r_, s_]] any Mathematica expert knows that this is something you don't want to do if it is used in a core loop. Instead of pattern matching, one can change the form to Function and it'll speed up. Also, he used â€œBlockâ€, which is designed for local variables and the scope is dynamic scope. However the local vars used in this are local constants. A proper code would use â€œWithâ€ instead. (in lisp, this is various let, let*. Lispers here can imagine how lousy the code is now.) Here's a improved code. The timing of this code is about 0.2 of the original. Also, optimization is purely based on code doodling. That is, i do not know what his code is doing, i do not have experience in writing a ray tracer. All i did is eyeballing his code flow, and improved the form. norm=Function[#/Sqrt@(Plus@@(#^2))]; delta=Sqrt[$MachineEpsilon]; myInfinity=10000.; Clear[RaySphere]; RaySphere = Compile[{o1, o2, o3, d1, d2, d3, c1, c2, c3, r}, Block[{v = {c1 - o1, c2 - o2, c3 - o3}, b = d1*(c1 - o1) + d2*(c2 - o2) + d3*(c3 - o3), discriminant = -(c1 - o1)^2 - (c2 - o2)^2 + (d1*(c1 - o1) + d2*(c2 - o2) + d3*(c3 - o3))^2 - (c3 - o3)^2 + r^2, disc, t1, t2}, If[discriminant < 0., myInfinity, disc = Sqrt[discriminant]; If[(t1 = b - disc) > 0., t1, If[(t2 = b + disc) <= 0., myInfinity, t2]]]]]; Remove[Intersect]; Intersect[{o1_,o2_,o3_},{d1_,d2_,d3_}][{lambda_,n_},Sphere [{c1_,c2_,c3_},r_]]:= Block[{lambda2=RaySphere[o1,o2,o3,d1,d2,d3,c1,c2,c3,r]}, If[lambda2â‰¥lambda,{lambda,n},{lambda2, norm[{o1,o2,o3}+lambda2 *{d1,d2,d3}-{c1,c2,c3}]}]] Intersect[{o1_,o2_,o3_},{d1_,d2_,d3_}][{lambda_,n_}, Bound[{c1_,c2_,c3_},r_,s_]]:= Block[{lambda2=RaySphere[o1,o2,o3,d1,d2,d3,c1,c2,c3,r]}, If[lambda2â‰¥lambda,{lambda,n}, Fold[Intersect[{o1,o2,o3},{d1,d2,d3}],{lambda,n},s]]] Clear[neglight,nohit] neglight=N@norm[{1,3,-2}]; nohit={myInfinity,{0.,0.,0.}}; Clear[RayTrace]; RayTrace[o_,d_,scene_]:= Block[{lambda,n,g,p},{lambda,n}=Intersect[o,d][nohit,scene]; If[lambda\[Equal]myInfinity,0,g=n.neglight; If[gâ‰¤0, 0,{lambda,n}=Intersect[o+lambda d+delta n,neglight] [nohit,scene]; If[lambda<myInfinity,0,g]]]] Clear[Create]; Create[level_,c_,r_]:= Block[{obj=Sphere[c,r]}, If[level\[Equal]1,obj, Block[{a=3*r/Sqrt[12],Aux}, Aux[x1_,z1_]:=Create[level-1,c+{x1,a,z1},0.5 r]; Bound[c,3 r,{obj,Aux[-a,-a],Aux[a,-a],Aux[-a,a],Aux[a,a]}] ] ]] Main[level_,n_,ss_]:= With[{scene=Create[level,{0.,-1.,4.},1.]}, Table[Sum[ RayTrace[{0,0,0}, N@norm[{(x+s/ss/ss)/n-1/2,(y+Mod[s,ss]/ss)/ n-1/2,1}],scene],{s,0, ss^2-1}]/ss^2,{y,0,n-1},{x,0,n-1}]] Timing[Export["image.pgm",Graphics@Raster@Main[2,100,4.]]] Note to those who have Mathematica. Mathematica 6 has Normalize, but that's not in Mathematica 4, so i cooked up above. Also, Mathematica 6 has AbsoluteTiming, which is intended to be equivalent if you use stop watch to measure timing. Mathematica 4 has only Timing, which measures CPU time. My speed improvement is based on Timing. But the same factor will shown when using Mathematica 6 too. I'm pretty sure further speed up by 0.5 factor of above's timing is possible. Within 2 more hours of coding. Jon wrote: Â«The Mathematica code is 700,000x slower so a 50% improvement will be uninteresting. Can you make my Mathematica code five orders of magnitude faster or not?Â» If anyone pay me $300, i can try to make it whatever the level of F# or OCaml's speed is as cited in Jon's website. ( [URL]http://www.ffconsultancy.com/languages/ray_tracer/index.html[/URL] ). Please write out or write to me the detail exactly what speed is required in some precise terms. If i agree to do it, spec satisfaction is guaranteed or your money back. PS Thanks Thomas M Hermann. It was fun. Xah âˆ‘ [URL]http://xahlee.org/[/URL] â˜„ [/QUOTE]

Verification