Time, Inverse, Simulation simple dynamical system, speed issue
- To: mathgroup at smc.vnet.net
- Subject: [mg77380] Time, Inverse, Simulation simple dynamical system, speed issue
- From: kristoph <kristophs.post at web.de>
- Date: Thu, 7 Jun 2007 04:03:33 -0400 (EDT)
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*)
]