| Author |
Comment/Response |
Stan
|
 06/30/12 12:02pm
In Response To 'Re: Re: Maximizing solutions of NDSolve'
---------
That is strange. I'm using Mathematica at university and it's a version 7. I can't get hold of it until Monday. However I downloaded a trial version of version 8 on my laptop and it indeed found the correct maximum in this example.
Nevertheless when I try to calculate what I actually need even v8 gives the wrong result:
a[t_] := (1 + 3 t/8)^(2/3)
q[q0_, t_] := q0 a[t]^-3
\[Omega][\[Kappa]02_, q0_, t_] :=
2 q[q0, t]^(1/2) Sqrt[\[Kappa]02/(2 a[t]^(1/2) q[q0, t]^(1/2)) +
Sin[t]^2]
s[(\[Kappa]02_)?NumericQ, (q0_)? NumericQ] :=
NDSolve[{y''[
x] + (4 q[q0, x]) (\[Kappa]02/(2 Sqrt[a[x] q[q0, x]]) +
Sin[x]^2) y[x] == 0, y'[0] == 0, y[0] == 0.5}, y, {x, 0, 200},
MaxSteps -> 500000, Method -> {Automatic}, MaxStepFraction -> 1/100]
n[\[Kappa]02_, q0_,
t_] := \[Omega][\[Kappa]02, q0, t]/
2 (Abs[y'[t] /. s[\[Kappa]02, q0]]^2/\[Omega][\[Kappa]02, q0,
t]^2 + Abs[y[t] /. s[\[Kappa]02, q0]]^2) - 1/2
mueff[\[Kappa]02_, q0_, t_] := 1/(2 t) Log[2 n[\[Kappa]02, q0, t]]
ptdata = Table[{\[Kappa]02,
First[Re[mueff[\[Kappa]02, (32 Pi)^2, 30.5 Pi]]]}, {\[Kappa]02, 0,
1.5, 0.001}];
int = Interpolation[ptdata]
Maximize[{int[x], 0 < x < 1.5}, x]
gives
{0.121158, {x -> 0.101868}}
when the actual global maximum is around x=0.2:
In[161]:= int[0.2]
Out[161]= 0.126081
URL: , |
|
