Re: Time, Inverse, Simulation simple dynamical system, speed issue
- To: mathgroup at smc.vnet.net
- Subject: [mg77428] Re: Time, Inverse, Simulation simple dynamical system, speed issue
- From: kristoph <kristophs.post at web.de>
- Date: Fri, 8 Jun 2007 05:36:22 -0400 (EDT)
- References: <f48f6c$e1p$1@smc.vnet.net>
On 7 Jun., 10:21, kristoph <kristophs.p... 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 o= > f 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*) > ] Sorry. Here is the right code. Copy plain text did not work. Thanks again for help. << 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}}; Div = {{2, 1}, {2, 4}, {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<=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<=1000, t++,(*begin loop for the dynamical system consiting of t = 1,...,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]<=1/3, divi = Div[[1, All]], If[1/3 < x[t]<= 2/3, divi = Div[[2, All]], divi = Div[[3, All]]]]; Wel = N[Inverse[IdentityMatrix[2] - ThetaMatrix.Lam.Rho, 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*) ]