In the past I had done an approximate match of sunspot number plots to a
"beat" of three cycles. I have posted several results to my own egroups.
Here is a picture of the sunspot activity 1610 to 1994:
http://profile.imeem.com/GUmj0c/photo/b-mFxDqCxu/
http://www.google.com/codesearch?q=fractal+dimension+lang%3Amathematica&hl=en&btnG=Search+Code
http://www.internationalmathematicasymposium.org/IMS99/paper25/ims99paper25.nb
Using this correlation dimension method from the Italian paper I get
very close to a Cantor dimension
for the sunspot data.
Working with it the older data before 1700 seems suspect?
Correlation dimension:
s1=0.6376053551998503
This method appears to be a very good one for one dimensional or below
times series.
My own count based estimate gives ( not very good because I use a Floor[]
to limit the m\number of counts needed):
s2= 0.3781094484568289
s1+s2~ 1
On the calculation of the "beat" function: a Weierstrass like approach
in products like:( for fundamental frequency w0 and dimensiuon s)
amp=175.1
w0=(10.9259 +10.9231)/2; (* the average frequency from the neural net
solution*)
s=0.6376053551998503;
amp2 = Max[Table[g[x] - g[0], {x, 1, Length[b]}]]
p0[x_] = x*Product[Sin[w0^(s*n)*x]^2/w0^((s - 2)*n), {n, -2, 2}]
g[x_] = Fit[b, {1, p0[x]}, x]
g2 = Plot[amp*(g[x] - g[0])/amp2, {x, 1, Length[b]}, PlotRange -> All]
g3 = ListPlot[b, PlotStyle -> {Hue[1]}, PlotJoined -> True]
Show[{g2, g3}]
Gives an approximation with the observed beat structure.
It isn't really a very good result.
I hope someone comes up with a better answer.
Mathematica:
Directory[]
raw = Flatten[ReadList["yearrg.dat", Number, RecordLists -> True]]
Dimensions[raw]
Clear[b]
b = Table[raw[[n]], {n, 2, 772, 2}]/175.1
(*Form From : Testing Chaos and Fractal Properties in Economic Time Series
Maria I.Loffredo
Dipartimento di Matematica, Università di Siena, I - 53100 Siena, Italy
e - mail : loffredo@[EMAIL PROTECTED]
*)
P1[i_] := P1[i] = b[[i]]
P2[i_] := P2[i] = {b[[i]], b[[i + 1]]}
P3[i_] := P3[i] = {b[[i]], b[[i + 1]], b[[i + 2]]}
Hs[s_] := If[s > 0, 1, 0] (* Heaviside function *)
dist[P_, Q_] := N[Sqrt[Sum[(P[[k]] - Q[[k]])2, {k, 1, 2}] ]]
Corr[R_] := N[(Sum[Hs[R - dist[P3[i], P3[j]]], {i, 1, Length[b] - 2}, {j,
1,
i - 1}]
+ Sum[Hs[
R - dist[P3[i], P3[j]]], {i, 1, Length[b] - 2}, {j, i +
1, Length[b] - 2}])/(Length[b] - 2)2]
datacorr = Table[{R, Corr[R]}, {R, 0.125, 1, 0.0625}]
Needs["Graphics`Graphics`"]
ga = LogLogListPlot[datacorr]
LP = LogLogListPlot[datacorr, DisplayFunction -> Identity];
lp = Table[LP[[1, i, 1]], {i, 1, Length[datacorr]}];
fb[x_] = Fit[lp, {1, x}, x] // Simplify
0.04205939538432238+ 0.6376053551998503 x
gb = Plot[fb[x], {x, -1, 0}]
Show[{gb, ga}]
N[Log[2]/Log[3]]
0.6309297535714573`
Mathematica :
Directory[]
raw = Flatten[ReadList["yearrg.dat", Number, RecordLists -> True]]
Dimensions[raw]
b = Table[Floor[raw[[n]]], {n, 2, 772, 2}]
bm = Max[b]
c = Table[Count[b, n], {n, 1, bm}]
Max[c]
ListPlot[c]
d = Table[Count[c, n], {n, 0, Max[c]}]
e = Delete[Union[Table[If[Log[1 + c[[n]]] == 0, {}, N[{Log[n], Log[
1 + c[[n]]]}]], {n, 1.Length[c]}]], 1]
g0 = ListPlot[e, PlotJoined -> True]
h[x_] = Fit[e, {1, x}, x]
2.639117941421544- 0.3781094484568289 x
g1 = Plot[h[x], {x, 0, 5}]
Show[{g0, g1}]
This following program is a very crude matching of a Cantor cartoon
Biscuit (
Besicovitch -Ursell ) function with
the Sunspot data which the correlation dimension study says has a
dimension very near the Cantor one.
Even at that is is the best fitting I've gotten: better than several
days trying with Weierstrass product functions!
The major point is that it predicts a major cut off in solar activity in
about 30 years,
but a steady raise until then. If so that would be a good thing:
a mini-ice age is just what we need. Until then for the next
part of our century with solar activity going up and
greenhouse gas warming as well,
it is going to get hot.
The scale of solar warming is much slower of less magnitude
than greenhouse/ CO2 warming by a factor of about 10
I think.
Now, to the process that has a Cantor type of dynamics in the Sun?
My guess is that it is a Lorenz type " weather " circulation like that
that gives Jupiter
it's spots as well. The sunspot number is like the number of hurricanes
per year in earth weather terms. The sunspots like eyes of huuicanes
are slower / cooler than the storm they are central to.
My Cantor output is more like the solar flux output than the sunspot
number.
Mathematica:
Clear[f, g, h, k, ff, kk, ll]
f[x_] := x /; 0 <= x <= 1/3
f[x_] := 0 /; 1/3 < x <= 2/3
f[x_] := x /; 2/3 < x <= 1
ff[x_] = f[Mod[Abs[x], 1]];
s0 = Log[2]/Log[3];
kk[x_] = Sum[ff[3^k*x]/3^(s0*k), {k, 0, 20}];
ll[x_] = Sum[ff[3^k*(x + 1/2)]/3^(s0*k), {k, 0, 20}];
(* adjusting axes crudely*)
Floor[6000/(772/2)];
ga = Table[175.1*kk[(n - 1610)/10000]/2, {n, 4000 + 1610, 11000 + 1610,
15}];
g0 = ListPlot[ga, PlotJoined -> True]
Directory[];
raw = Flatten[ReadList["yearrg.dat", Number, RecordLists -> True]];
Dimensions[raw]
b = Table[raw[[n]], {n, 2, 772, 2}];
amp = Max[b]
c = Apply[Plus, b]/Length[b]
{772}
175.1`
34.19611398963731`
ticks = Map[{#, StringForm["`1`", Sequence @[EMAIL PROTECTED]
raw[[1 + 2*#]]]} &,
Range[1, Length[b], 100]];
g1 = ListPlot[b, PlotJoined -> True, AspectRatio -> 0.2,
Frame -> True, FrameTicks -> {ticks, Automatic, None,
None}, PlotStyle -> Green, FrameLabel -> "sunspot #", Axes -> None];
Show[{g0, g1}]
Here is a picture of the two curves together:
http://profile.imeem.com/GUmj0c/photo/uHsjrgaAz1/
In conclusion, this fractal approach to sunspot numbers
seems to have some virtues over a simple "beats" of cycles approach
or the Mathematica neural net approach.
|
