I am trying to find the times of maxima and minima in order to calculate the period of some oscillations. However for certain values of my system (essentially a driven simple harmonic oscillator with an additional term) the period jumps hugely to give ridiculous answers, some of order 10^8 times the expected results and I receive the errors InterpolatingFunction::dmval: and FindMinimum::nrnum:. It appears that the FindMaximum function isn't recognising peaks sometimes and is continuing infinitely. My code is given below:
kR = 5.822559;
kZ = 2.62326;
amax = 2
adelta = 0.2
a[Z_, R_, T_] := a0*Exp[-(((Z - T)/kZ)^2)/2]*Exp[- ((R/kR)^2)];
Tmin = -20;
Tmax = 200;
Zmin = 0;
Zmax = 10;
Rmin = -5;
Rmax = 5;
T0 = -20;
Ttest = 50;
tablePeriod1 =
Table[2*Abs[((T /.
Last[FindMinimum[
u1[0, 0, T] /.
NDSolve[{D[
u1[Z, R, T], {T, 2}]*(1 + 1.5*D[u1[Z, R, T], T]^2 +
0.5*D[v1[Z, R, T], T]^2) + u1[Z, R, T] -
D[v1[Z, R, T], T]*D[u1[Z, R, T], T]*
v1[Z, R, T] == -0.25*D[a[Z, R, T]^2, Z],
D[v1[Z, R, T], {T, 2}]*(1 + 1.5*D[v1[Z, R, T], T]^2 +
0.5*D[u1[Z, R, T], T]^2) + v1[Z, R, T] -
D[u1[Z, R, T], T]*D[v1[Z, R, T], T]*
u1[Z, R, T] == -0.25*D[a[Z, R, T]^2, R],
(D[u1[Z, R, T], T] /. T -> T0) == 0, v1[Z, R, T0] == 0,
u1[Z, R, T0] == 0, (D[v1[Z, R, T], T] /. T -> T0) ==
0}, {u1, v1}, {Z, Zmin, Zmax}, {R, Rmin, Rmax}, {T,
Tmin, Tmax}, MaxSteps -> Infinity,
Method -> "Automatic"][[1]], {T, Ttest},
WorkingPrecision -> 30]]) - (T /.
Last[
FindMaximum[
u1[0, 0, T] /.
NDSolve[{D[
u1[Z, R, T], {T, 2}]*(1 + 1.5*D[u1[Z, R, T], T]^2 +
0.5*D[v1[Z, R, T], T]^2) + u1[Z, R, T] -
D[v1[Z, R, T], T]*D[u1[Z, R, T], T]*
v1[Z, R, T] == -0.25*D[a[Z, R, T]^2, Z],
D[v1[Z, R, T], {T, 2}]*(1 + 1.5*D[v1[Z, R, T], T]^2 +
0.5*D[u1[Z, R, T], T]^2) + v1[Z, R, T] -
D[u1[Z, R, T], T]*D[v1[Z, R, T], T]*
u1[Z, R, T] == -0.25*D[a[Z, R, T]^2, R],
(D[u1[Z, R, T], T] /. T -> T0) == 0, v1[Z, R, T0] == 0,
u1[Z, R, T0] == 0, (D[v1[Z, R, T], T] /. T -> T0) ==
0}, {u1, v1}, {Z, Zmin, Zmax}, {R, Rmin, Rmax}, {T,
Tmin, Tmax}, MaxSteps -> Infinity,
Method -> "Automatic"][[1]], {T, Ttest},
WorkingPrecision -> 30]]))], {a0, 0.2, amax, adelta}]
And the results were:
{6.221796534875934665845156568, 6.279119955159784643716071509, \
6.318407692782765006990372270, 6.302787259988871206694852995, \
6.359180123726873636918626166,
1.48489160439112764037865586341*10^8, 6.491986737115338573852543549, \
46.069994156461475901853847230, 6.727526312055413488476783498, \
20.818040498618419289660167133}
I expected answers around 2Pi, but while some results stay close to this value, others jump around hugely. Does anybody have any ideas as to what I might be doing wrong?
Sorry if I've formatted this terribly, I can't seem to find a way to preview the post. Thank you for you time.
Attachments:
