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"}]