Re: Compiling numerical iterations
- To: mathgroup at smc.vnet.net
- Subject: [mg129940] Re: Compiling numerical iterations
- From: Peter Klamser <klamser at googlemail.com>
- Date: Wed, 27 Feb 2013 03:06:38 -0500 (EST)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- Delivered-to: l-mathgroup@wolfram.com
- Delivered-to: mathgroup-newout@smc.vnet.net
- Delivered-to: mathgroup-newsend@smc.vnet.net
- References: <kgab2i$c7f$1@smc.vnet.net>
Hi Cornelius, your code is not slow because you use Mathematica but because you are not writing in functional and Mathematica style. Use /@, @@, pure Functions like #1^2& etc. to realize your algorithm. Plus @@ {a, b, c, d} is much faster then procedural programming. Kind regards from Peter 2013/2/26 firlefranz <cornelius.franz at gmx.net>: > Thanks a lot! To be honest, some of the commands Ray is using, I've never seen before. I stopped using mathematica before version 5 came out. > > Coming back to Peters statement of exporting the code from Mathematica to C. How this can be done starting from my or Ray's code? There is an automated C-code-gernerator implemented in Mathematica 9, am I right? > > > As a step forward to the first little peace of code, here is another try, which is not really optimized. The be concrete, I try to simulate the autocorrelation function of a random walk, which is doing a step at none equally distant time steps. This has to been done for a long random walk for many particles. Instead of doing an average over many random walks and calculate one autocorrelation function I want to simulate many correlation functions and make an average over them. Since the time steps are non equal, I wrote a sub-function, which creates a new time axis and taking the necessary value for the random walk from the first table. > > Here is what I come up with. It's running in a reasonable time for one particle, but for a real statistic ensemble, I have to do it over 1.000.000 particles for a long time. Optimizing this or (probably better) exporting it to C would hopefully help a lot. So how to export it? > > Clear["Global`*"] > SeedRandom[1234567890]; > zeitmax = 100;(* muss ein Integer sein *) > numteilchen = 1; > tauj = 1; > corr = Table[0, {i, 1, zeitmax/10}]; > > SucheIndex[zeitliste_, zeit_, maxindex_] := > Module[{i}, > For[i = 1, i <= maxindex, i++, > If[zeitliste[[i]] > zeit, Break[]]; > ]; > i - 1 > ]; > > For[j = 1, j <= numteilchen, j++, > (* Zeitachse generieren von 0 bis zeitmax *) > t = 0; > i = 1; > tabzeit = {}; > time = AbsoluteTiming[While[True, > tabzeit = Append[tabzeit, t]; > dt = -tauj*Log[1 - RandomReal[]]; > If[t > zeitmax, Break[]]; > t = t + dt; > i++; > ]; > ]; > Print[time]; > maxidx = i; > > (* Random Walk *) > time = AbsoluteTiming[ > tabwalk = Table[0, {i, 1, maxidx}]; > a = 0; > For[i = 1, i <= maxidx, i++, > tabwalk[[i]] = a; > If[RandomReal[{-1, 1}] > 0, a = a + 1, a = a - 1]; > ]; > ]; > Print[time]; > (*tabwalk=Table[Subscript[b, i],{i,1,maxidx}];*) > > (* Korrelationsfunktion berechnen *) > time = AbsoluteTiming[ > For[k = 1, k <= zeitmax/10, k++, > For[n = 1, n <= zeitmax/10*9, n++, > corr[[k]] = > corr[[k]] + > tabwalk[[SucheIndex[tabzeit, n - 1, maxidx]]]* > tabwalk[[SucheIndex[tabzeit, n + k - 2, maxidx]]]/ > zeitmax*10/9/numteilchen; > (*Print[corr//N];*) > ]; > ]; > ]; > Print[time]; > Print[corr // N]; > ]; > Table[{tabzeit[[i]], tabwalk[[i]]}, {i, 1, maxidx}] > corr // N >