MathGroup Archive 2004

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Plotting a function of an interpolated function


DrBob wrote:

>I'll take the last problem first -- the Plot:
>
>Plot[W[\[Phi]], {t, 0, 100/m}]
>
>It can't be evaluated until Phi is computed, but regardless of what Phi turns out to be, the plot is determined anyway:
>
>Clear[m, mp, G, P, \[Rho], \[Phi], W,
>   Inflaton]
>V = Function[x, (1/2)*m^2*x^2];
>\[Rho] = Function[x, D[x, t]^2/2 + V[x]];
>H = Function[x, Sqrt[((8*Pi*G)/3)*
>       \[Rho][x]]];
>P = Function[x, D[x, t]^2/2 - V[x]];
>Inflaton = Derivative[2][\[Phi]][t] +
>     3*H[\[Phi][t]]*Derivative[1][\[Phi]][t] +
>     Derivative[1][V][\[Phi][t]];
>W = Function[x, P[x]/\[Rho][x]];
>W[\[Phi]]
>
>-1
>
>W is constant, so the Plot should show -1 across the board (if it works at all).
>
>Now for the problem with Phi. Your definition of Phi isn't the best, but if you must use it, the Plot needs Evaluate:
>
>m = 10^16;
>mp = 1.2211*10^19;
>G = mp^(-2);
>Clear[\[Phi]]
>\[Phi] = \[Phi][t] /. First[NDSolve[
>      {Inflaton == 0, Derivative[1][\[Phi]][
>         0] == 0, \[Phi][0] == 1}, \[Phi][t],
>      {t, 0, 100/m}, MaxSteps -> 10^6]];
>Plot[Evaluate[W[\[Phi]]], {t, 0, 100/m}]
>
>But I don't see the expected constant -1. Instead, I see a Sin wave, I think due to noise.
>
>(W simplified to -1 when its argument was a symbol, but when the argument is numeric, the fraction P[x]/rho[x] appears to be vulnerable to numerical problems.)
>
>Define Phi this way, and I get a much more reasonable result:
>
>Clear[\[Phi]]
>\[Phi] = \[Phi] /. First[NDSolve[
>      {Inflaton == 0, Derivative[1][\[Phi]][
>         0] == 0, \[Phi][0] == 1}, \[Phi],
>      {t, 0, 100/m}, MaxSteps -> 10^6]];
>Plot[W[\[Phi][t]], {t, 0, 100/m}]
>
>The curve wavers, but it's -1 to machine precision, as expected.
>  
>
Sure, thanks for the help. I ended up rewriting the function like so:

mGUT = 10^16;
mp = 1.2211*10^19;
G = mp^(-2);

P = Function[x, D[x, t]^2/2 -
     V[x]];
W = Function[x, P[x]/\[Rho][x]];

parameters["m"] =
   {m -> mGUT};

plotW[phi0_, tmax_,
   opts___] :=
  Plot[Evaluate[
    W[Phi[phi0, tmax]] /.
     parameters["m"]],
   {t, 0, tmax}, PlotRange ->
    All, Frame -> True,
   PlotPoints -> 2500,
   FrameLabel -> {"t",
     "W(t)",
     "Pressure coefficient",
     ""}, opts]

line[5] = {Hue[0.3]};
plotW[10^19, 10000/mGUT,
   PlotStyle -> {line[5]}]; 

And it plotted exactly what I was looking for, which was the reheat time 
after inflation when w diverges significantly from -1.

>Bobby
>  
>
Thanks again!

-- Adam Getchell


  • Prev by Date: Re: How input stacked characters with vertical bar
  • Next by Date: Re: GUIKit - ScrollPane Tables within Wizard
  • Previous by thread: Re: Plotting a function of an interpolated function
  • Next by thread: NonlinearRegress+NIntegrate