Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

making a fast simulation faster? Mathm. 6.0

• To: mathgroup at smc.vnet.net
• Subject: [mg79259] making a fast simulation faster? Mathm. 6.0
• From: kristoph <kristophs.post at web.de>
• Date: Sun, 22 Jul 2007 04:24:30 -0400 (EDT)

Dear All,

I know it is much to ask but I'm trying desperately to make the code
below faster. Due to important help from this group for parts of the
system it is already fast but it is not yet fast enough. I'm using an
intel 2 GHz. machine which needs about 6 seconds in the below
configuration.
I'm simulating a stochastic dynamical system which I simulate 5 times
for each configuration. There exist a total of 499^2 configurations.
You might imaging that time is a big issue. For each simulation I
extract data which I process into average date of the 5 simulations.
This is needed for further processing.


I'm very grateful for any suggestions to the simulation much faster.
If you have any questions regarding the code please do not hesitate to
email me.

Thanks a lot,
Kristoph


divi = {{3, 1}, {3, 1}, {2, 0}};
inve = 2;
wert = Dimensions[divi][[2]];
rela = Table[
   divi[[i, k]]/Total[divi[[i]]], {i, 1, Length[divi]}, {k, 1,
    wert}];
lambdaGOP = Table[Mean[rela[[All, i]]], {i, 1, wert}];
lambdaEQ = Table[1/wert, {i, 1, wert}];
lambda = Transpose[{lambdaGOP, lambdaEQ}];

Clear[a, b, theta, preCalc, iterate, meanadd, meantransform, varadd,
  vartransform];
ClearAll[listmean, listvar, meantransform, data, intervall];
rho[h_, g_] := {{h/500, 0}, {0, g/500}};
relativesver[a_, b_] := (rho[h, g][[1, 1]] a)/
   Total[rho[h, g].{a, b}];
iterate[{a_, b_}] :=
  N[preCalc[a, b].divi[[Ceiling[3 RandomReal[]]]]];
rw[w_] := relativesver @@@ lo[[w]];
listmean[g_] := listmean[g] = Join[listmean[g - 1], bigloop[[1]]];
listvar[g_] := listvar[g] = Join[listvar[g - 1], bigloop[[2]]];
simuls = 5;
periods = 800;
band = 0.8;

Timing[For[g = 1, g <= 3, g++, Print["Lauf-", g];
  bigloop =
   Reap[For[h = g - Floor[g*band], h <= 499, h++, verm = {a, b};
      theta[a_, b_] =
       Simplify[
        Table[(lambda[[k, i]]*rho[h, g][[i, i]]*verm[[i]])/
          Sum[lambda[[k, j]]*rho[h, g][[j, j]]*verm[[j]], {j, 1,
            inve}], {i, 1, inve}, {k, 1, wert}]];
      preCalc[a_, b_] =
       Inverse[IdentityMatrix[inve] -
           theta[a, b].lambda.rho[h, g]].theta[a, b] // Simplify;
      lo =
       Reap[For[n = 1, n <= simuls,
          n++' Sow[
             Delete[NestList[iterate, {5, 5}, periods],
              1]];];(*end For*)][[-1, 1]];(*end Reap*)
      meanreihe = 1/simuls*Sum[rw[x], {x, 1, simuls}];
      totalstandardabw =
       Sqrt[1/simuls*Sum[Variance[rw[x]], {x, 1, simuls}]];
      Sow[meanreihe[[-1]], "EndwertMeanReihe"];
      Sow[totalstandardabw, "totalstandardabw"];
      If[rho[h, g][[1, 1]] >= rho[h, g][[2, 2]] &&
        totalstandardabw < 0.01, h = 500];
      Clear[meanreihe, totalstandardabw];], {"EndwertMeanReihe",
      "totalstandardabw"}][[-1]];
  bigloop[[2, 1]] =
   PadLeft[bigloop[[2, 1]],
    Length[bigloop[[2, 1]]] + (g - Floor[g*band] - 1)];
  bigloop[[1, 1]] =
   PadLeft[bigloop[[1, 1]],
    Length[bigloop[[1, 1]]] + (g - Floor[g*band] - 1)];
  If[g == 1, listmean[g] = bigloop[[1]]; listvar[g] = bigloop[[2]],
   listmean[g]; listvar[g]];
  If[g >= 3, listmean[g - 2] =.; listvar[g - 2] =.];
  Clear[a, b, theta, preCalc];]]