Tracking Differential Equations

*To*: mathgroup at smc.vnet.net*Subject*: [mg128227] Tracking Differential Equations*From*: williamduhe at hotmail.com*Date*: Thu, 27 Sep 2012 03:06:49 -0400 (EDT)*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com*Delivered-to*: l-mathgroup@wolfram.com*Delivered-to*: mathgroup-newout@smc.vnet.net*Delivered-to*: mathgroup-newsend@smc.vnet.net

Bellow is a program intended to track three properties as a function of time via differential equations. a[t] T[t] and Rm[t]. a[t] roughly equals the volume of the universe T[t] is the temperature of the universe and Rm[t] is the the matter density of the universe. ( these are all approximations!!! LOL ) What I want to do is find a way to track how a[t] has changed over the course of one of the cycles demonstrated by the syntax. In other words I would like to note by what amount a[t] has changed from the first point where T[t] = Rm[t] to the second point where T[t] = Rm[t]. I would like to plot the change in a[t] as a function of my R value. Then perhaps a 3D plot of a[t] as a function of alpha and R. Thanks ahead of time for any help. Best wishes, Duhe curveColor[n_Integer] := Hue[Mod[2/3 - 2 + Sqrt[5]*(n - 1), 1], 0.6, 0.6]; Manipulate[ Module[{M, g, ti, tf, m, RB, V, i, rR, pR, pM, r, p, imgSize, aspRatio, s1, plt}, M = 100; g = 1; ti = 0;(*initial time*)tf = 12000;(*final plot time*)m = 1; RB[t] = (m*T[t]^3*g)/(2*\[Pi])^(3/2)*(m/T[t])^(3/2)*E^(-m/T[t]); V = 5*10^-12*M^4; i = Sqrt[(g*m^4 - V + g/(2*\[Pi])^(3/2)*E^-1)/(3*M^2)]*m; rR = g*T[t]^4;(*Energy Density-Radiation*)pR = g*T[t]^4;(*Pressure- Radiation*)pM = 0;(*Pressure-Phi Particle*) r = rR + Rm[t] - V;(*Overall Energy Density*) p = pR + V;(*Overall Pressure*)imgSize = 150; s1 = NDSolve[{a''[t] == -a[t]*1/(6*M^2) (2*g*T[t]^4 + 2 V + Rm[t]), T'[t] == (R*Rm[t])/(4* T[t]^3) + (Sqrt[3/2]*\[Alpha]^2)/(4*Sqrt[T[t]]* m^(11/2))*((RB[t])^2 - Rm[t]^2) - (a'[t]*T[t])/a[t], Rm'[t] == -R* Rm[t] + (Sqrt[3/2]*\[Alpha]^2*T[t]^(5/2))/ m^(11/2)*((RB[t])^2 - Rm[t]^2) - 3 a'[t]/a[t]*Rm[t], T[0] == m, Rm[0] == g/(2*\[Pi])^(3/2)*E^-1, a'[0] == i, a[0] == 1}, {a[t], T[t], Rm[t]}, {t, ti, tf}, Method -> "BDF"][[1]]; Row[{Show[ plt = Plot[ Evaluate[{Tooltip[Rm[t]/RB[t], "Rm[t]"], Tooltip[20 T[t], "T[t]"], Tooltip[a[t], "a[t]"]} /. s1], {t, ti, tf}, Frame -> True, Axes -> False, ImageSize -> imgSize, PlotRange -> {-1, 22}, AspectRatio -> 2.5], FrameTicks -> {{Automatic, (FrameTicks /. AbsoluteOptions[plt, FrameTicks])[[2]] /. {y_, yl_?NumericQ, r__} :> {y, yl/20, r}}, {(FrameTicks /. AbsoluteOptions[plt, FrameTicks])[[1]] /. {x_, xl_?NumericQ, r__} :> {x, xl/1000, r}, Automatic}}, FrameLabel -> {{Row[{Style["Rm[t]", Bold, curveColor[1]], " & ", Style["a[t]", Bold, curveColor[3]]}], Style["T[t]", Bold, curveColor[2]]}, {"Time/1000", None}}], ParametricPlot[Evaluate[{T[t], a[t]} /. s1], {t, ti, tf}, Frame -> True, FrameLabel -> {"T[t]", "a[t]"}, Axes -> False, AspectRatio -> 2, ImageSize -> imgSize, PlotRange -> {{0.14, 0.5}, {2, 7}}, ColorFunction -> (ColorData["Rainbow"][#3] &)]}]], {{\[Alpha], 0, "|\[Alpha]|"}, 0, 5, Appearance -> "Labeled"}, {{R, 0}, 0, 0.001, Appearance -> "Labeled"}]