Program 5 - Correlation Dimension

Program Code

Off[General::"spell1",General::"spell"];

BeginPackage["CorrDim`"];

Unprotect[CorrDim];

CorrDim::usage = "CorrDim[points, r] calculates the approximate correlation integral of a list of 
	points in N dimensions for a specific value of r.  CorrDim examines the list of points and 
	chooses the appropriate (1- or N-Dimensional) algorithm.";	

Begin["`Private`"];

cd1 = Compile[{{t, _Real, 1}, {r, _Real}}, Module[{k = 0, i, j, len},
	len = Length[t];
	For[i = 1, i <= len, i++,
		For[j = i + 1, j <= len, j++,
			If[Abs[t[[i]] - t[[j]]] < r, k++];
		]
	];
	Return[If[k!=0,N[Log[2 k / (len^2 - len)] / Log[r]], 0.0]]
]];
cdn = Compile[{{t, _Real, 2}, {r, _Real}}, Module[{k = 0, i, j, len, d, r2},
	len = Length[t];
	r2 = r^2;
	For[i = 1, i <= len, i++,
		For[j = i + 1, j <= len, j++,
			d = t[[i]] - t[[j]];
			If[d.d < r2, k++];
		]
	];
	Return[If[k != 0, N[Log[2 k / (len^2 - len)] / Log[r]], 0.0]]
]];
			
RealQ[r_] := NumericQ[r] && Im[r] == 0;


CorrDim::badlist = "List of points (instead of `1`) expected at position 1 of CorrDim.";
CorrDim::badr = "Value of r in position 2 must be a positive real value (not equal to 1) instead of `1`.";

CorrDim[points_, rval_] := Module[{result},
	If[!RealQ[rval] || rval <= 0 || rval == 1, Message[CorrDim:badr, rval]; result = $Failed,
		If[VectorQ[points, RealQ],
			result = cd1[N[points], rval],
			If[MatrixQ[points, RealQ] && Depth[points] == 3,
				result = cdn[N[points], rval],
				Message[CorrDim::badlist, points]; result = $Failed
			]
		]
	];
	Return[result] /; result =!= $Failed
];

End[];

SetAttributes[CorrDim, ReadProtected];
Protect[CorrDim];

EndPackage[];

On[General::"spell1",General::"spell"];

Examples

Create the function governing the system.

f[y_] := 4  y^3 - 3  y;

Compute the first 400 iterations of the system.

p1 = NestList[f, 0.3, 399];

Look at last 100 iterations using the Take command.  The -100 means take the last 100 elements.

p = Take[p1, -100];

In the following table the the first entry is the value of r, the second is the "approximate correlation integral for that value of r".

t1 = Table[{r, CorrDim[p, r]}, {r, 0.01, 0.1, 0.01}]

plot1 = 
      ListPlot[t1, PlotRange -> {{0, .1}, {.6, 1}}, 
        Ticks -> {Range[0, .1, .05], Range[0, 1, .2]}, 
        AxesOrigin -> {0, .7}, PlotStyle -> PointSize[0.02], 
        DefaultFont -> {"Times", 12}]

Clear[x, y1];

y1[x_] = Fit[t1, {1, x}, x]

It appears that the correlation dimension is close to 1. The least square line has very small slope. 
Let us see how well it fits the obtained data.

plot2 = Plot[y1[x], {x, 0, .1},  PlotRange -> {{0, .1}, {0, 1}}];

Show[plot1, plot2];

Define the iteration functions for the 2-dimensional system.

f[x_, y_] := 3.8  x (1 - x - y);  g[x_, y_] := 3.8  x y

x = .4; y = .3;

q1 = Table[{u = x, x = f[x, y], y = g[u, y], m = {x, y}}, {400}];

q = q1[[Range[300, 399], 4]];

Now approximate the correlation integrals, using the same syntax as in the 1-dimensional example.

t2 = Table[{r, CorrDim[q, r]}, {r, 0.1, 0.11, 0.001}]

plot3 = 
      ListPlot[t2, PlotStyle -> PointSize[0.02], 
        PlotRange -> {{.1, .11}, {0, 2}}];

Clear[x, y2];

y2[x_] = Fit[t2, {1, x}, x]

plot4 = Plot[y2[x], {x, .1, .11}, PlotRange -> {{.1, .11}, {0, 2}}];

Show[plot3, plot4];