Re: Time, Inverse, Simulation simple dynamical system, speed issue
- To: mathgroup at smc.vnet.net
- Subject: [mg77421] Re: Time, Inverse, Simulation simple dynamical system, speed issue
- From: Jean-Marc Gulliet <jeanmarc.gulliet at gmail.com>
- Date: Fri, 8 Jun 2007 05:32:44 -0400 (EDT)
- Organization: The Open University, Milton Keynes, UK
- References: <f48f6c$e1p$1@smc.vnet.net>
kristoph 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*)
> ]
Hi Kristoph,
Some of the characters have been messed up (I guess they stand for "less
than or equal to" followed by some values). You should post again your
code being sure that the incriminated character are translated either in
input form <= or in full form LessEqual[].
Regards,
Jean-Marc