Re: Time, Inverse, Simulation simple dynamical system, speed issue

*To*: mathgroup at smc.vnet.net*Subject*: [mg77438] Re: [mg77380] Time, Inverse, Simulation simple dynamical system, speed issue*From*: DrMajorBob <drmajorbob at bigfoot.com>*Date*: Fri, 8 Jun 2007 05:41:34 -0400 (EDT)*References*: <21255603.1181208515635.JavaMail.root@m35>*Reply-to*: drmajorbob at bigfoot.com

1) Some of that is gibberish (= E2 = 89 =, for instance); copy as plain text, THEN paste into e-mail. 2) Since you're setting RelWel={} every iteration AND Clearing it ever y iteration, why append anything to it or use it at all? You're actually appending to Value[t]. But see point (3). 3) AppendTo is a slow way to gather results over hundreds or thousands of iterations. Look into Sow and Reap, instead. For instance, AppendTo[ThetaTable, ThetaMatrix]; would become Sow[ThetaMatrix,"thetaTable"] AppendTo[RelWel, Rho[[1, 1]]*Wel[[1]]/Total[Rho.Wel]] would become Sow[Rho[[1, 1]]*Wel[[1]]/Total[Rho.Wel],"Value"] and Do[AppendTo[Value[t], RelWel[[t]]], {t, 1, 1000}] would just go away. You'd start the loop with something like results=Reap at Timing... or Timing@Reap@For... and "results" would include a tagged list for each Sow in the loop. You'd no longer use ThetaTable, RelWel, or Value[t]. It's possible WRI has greatly improved AppendTo at some point and I didn't hear about it, but it used to be, at least, well known that appending elements to a long list was prohibitively slow. I suspect AppendTo is using far more time than your calculations. 4) N[Inverse[IdentityMatrix[2] - ThetaMatrix.LamRho, Method -> DivisionFreeRowReduction].ThetaMatrix.divi, 1000] is NOT faster than the same thing without N. But using Wel in a later statement IS faster (probably), since Wel is approximate. Otherwise, the rational results are likely to grow huge numerators and denominators. It's possible you'd do just as well replacing N with Simplify, but probably not. 5) Do[AppendTo[Value[t], RelWel[[t]]], {t, 1, 1000}] is MUCH slower than the equivalent Value[t] = Join[Value[t],RelWel] (* if RelWel has 1000 elements *) or Value[t] = Join[Value[t],Take[RelWel,1000]] (* if RelWel has more than 1000 elements, but you only want 1000 *) But, in view of (3) above, RealWel and Value should go away entirely, and so does this Do loop. 6) If the theta function will be called thousands of times, it's much more efficient to define it this way: theta[i_, k_] := theta[i, k] = Lam[[i, k]]*Rho[[i, i]]* Wel[[i]]/Sum[Lam[[j, k]]*Rho[[j, j]]*Wel[[j]], {j, 1, 2}]; That doesn't actually work, since Wel changes during the run, but you can do this instead: Clear[a, b, temp] Wel = {a, b}; Lam = {{1/2, 2/3}, {1/2, 1/3}}; Rho = {{1/2, 0}, {0, 648/1000}}; temp[i_, k_] := Lam[[i, k]]*Rho[[i, i]]* Wel[[i]]/Sum[Lam[[j, k]]*Rho[[j, j]]*Wel[[j]], {j, 1, 2}]; ThetaMatrix := Evaluate@Table[Simplify@temp[i, k], {i, 1, 2}, {k, 1, 2}]; ThetaMatrix {{(125 a)/(125 a + 162 b), (125 a)/(125 a + 81 b)}, {(162 b)/( 125 a + 162 b), (81 b)/(125 a + 81 b)}} ThetaMatrix stores those results, and it always uses the current value of a and b. The theta function is not needed and neither is Clear[Wel], Clear[ThetaMatrix], or ThetaMatrix = Table[theta[i, k], {i, 1, 2}, {k,1, 2}]. (temp isn't needed again, either... unless you want to use a new Lam or Rho.) Replace every mention of Wel with {a,b}. That is, {a,b} = {5,5} instead of Wel = {5,5}, a instead of Wel[[1]], Rho.{a,b} rather than Rho.Wel, and of course {a,b} = N[Inverse[IdentityMatrix[2] - ThetaMatrix.LamRho, Method -> DivisionFreeRowReduction].ThetaMatrix.divi, 1000]; 7) There's a missing semicolon after ListPlot, before Clear. (Clear wasn't needed anyway.) In version 6 you'd need Print@ListPlot to see the charts . The ListPlot would be awfully dull, anyway... doesn't RealWel have only two points? You're supposedly appending to it each iteration, but you're also reinitializing it to {} each iteration, so isn't it just... ListPlot[Rho[[1, 1]]*Wel[[1]]/Total[Rho.Wel], PlotStyle -> {Hue[.8]}, PlotRange -> {0, 1}, ImageSize -> 350] I hope some of that helps. Bobby On Thu, 07 Jun 2007 03:03:33 -0500, kristoph <kristophs.post at web.de> wrote: > Attached you find a simulation of 5 dynamical systems each consisting > of 1000 Periods. For each period t = 1,...,1000 an inverse of a 2x2m > matrix needs to be computed. > It takes about 6 seconds to simulate the 5 dynamical systems. I would > be grateful for any hint in getting the simulation done much quicker. > It seems that calculating the inverse > > Wel = N[Inverse[ > IdentityMatrix[ > 2] - ThetaMatrix.LamRho, Method -> \ > DivisionFreeRowReduction].ThetaMatrix.divi, 1000] > > takes most of the time, therefore I tried to solve the problem > numerically (see above). This is a lot faster then using just > Inverse[...]. I tried LinearSolve, but it is slower then the above. > Decreasing the precision can result in a significant error, depending > on the parameter constellation. > > As you might have guessed I have to simulate not only 5 dynamical > systems for different parameter constellations. I would also like to > increase the order of the matrix that needs to be inverted. > > Here is the code. Thanks in advance. Kristoph > > > << Statistics`ContinuousDistributions`;(*packages and starting values > needed for the simulation*) > << Graphics`MultipleListPlot`; > > Wel = {5, 5}; > Lam = {{1/2, 2/3}, {1/2, 1/3}}; > Rho = {{1/2, 0}, {0, 648/1000}}; > RelDiv = {{2/3, 1/3}, {1/3, 2/3}, {1, 0}}; > theta[i_, k_] := Lam[[i, k]]*Rho[[i, i]]*Wel[[i]]/Sum[Lam[[j, > k]]*Rho[[j, j]]*Wel[[j]], {j, 1, 2}]; > > Timing[For[n = 1, n =E2=89=A4 5, n++,(*begin loop for the 5 simulations > of the > dynamical systeme*) > Wel = {5, 5};(*starting values*) > RelWel = {};(*needed for the plots, see below*) > ThetaTable = {};(*values needed for recursive calculations, saves > time*) > > For[t = 0, t =E2=89=A4 1000, t++,(*begin loop for the dynamical > system > consiting of t = 1000 Periods*) > Clear[ThetaMatrix]; > ThetaMatrix = Table[theta[i, k], {i, 1, 2}, {k, 1, 2}]; > AppendTo[ThetaTable, ThetaMatrix]; > AppendTo[RelWel, Rho[[1, 1]]*Wel[[1]]/Total[Rho.Wel]]; > Clear[Wel]; > x[t] = Random[];(*random pertubations are drawn from the above > matrix*) > If[x[t] =E2=89=A4 1/3, divi = Div[[1, All]], If[1/3 < x[t] =E2=89= > =A4 2/3, > divi = Div[[2, All]], divi = Div[[3, All]]]]; > > Wel = N[Inverse[IdentityMatrix[2] - > ThetaMatrix.LamRho, > Method -> DivisionFreeRowReduction].ThetaMatrix.divi, > 1000];(*the inverse needed for t + 1*) > ];(*end loop dynamical system*) > > Do[AppendTo[Value[t], RelWel[[t]]], {t, 1, 1000}]; > ListPlot[RelWel, PlotStyle -> {Hue[.8]}, PlotRange -> {0, 1}, > ImageSize -> 350] > Clear[RelWel, ThetaTable]; > > ] (*end loop simulations*) > ] > > > -- DrMajorBob at bigfoot.com

**Re: NDSolve, Do loop, and Plot**

**Re: Time, Inverse, Simulation simple dynamical system, speed issue**

**Re: Time, Inverse, Simulation simple dynamical system, speed issue**

**Re: Time, Inverse, Simulation simple dynamical system, speed issue**