Mathematica 9 is now available
Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2004
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2004

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

Search the Archive

Re: Re: the circle map

  • To: mathgroup at smc.vnet.net
  • Subject: [mg52305] Re: [mg52267] Re: the circle map
  • From: DrBob <drbob at bigfoot.com>
  • Date: Sun, 21 Nov 2004 07:23:50 -0500 (EST)
  • References: <cneuml$qbc$1@smc.vnet.net> <200411200841.DAA08749@smc.vnet.net>
  • Reply-to: drbob at bigfoot.com
  • Sender: owner-wri-mathgroup at wolfram.com

The mapping {x,y} -> {-a x - y, x} (without Mod or FractionalPart) has determinant one and maps an ellipse into itself for any starting point:

Clear[p];
a = 0.41209;
m = {{-a, -1}, {1, 0}};
p[n_] := p[n] = m.p[n - 1]
p[0] = {0.7, .65};
data = Table[p@n, {n, 0, 10000}];
short = Take[data, 10];
ListPlot[data, PlotRange -> All]
Show[Graphics@{Text @@@ Transpose@{Range@Length@short, short}}, Axes -> False]

The same thing occurs with Fractional Part, but there are false starts (transient behavior) until the series reaches a point whose corresponding ellipse falls entirely within the square -1<x<1, -1<y<1. FractionalPart "clips" the eventual ellipse to that square. For instance:

Clear[p];
a = 0.41209;
m = {{-a, -1}, {1, 0}};
p[n_] := p[n] = FractionalPart[m.p[n - 1]]
p[0] = {0.7, .65};
data = Table[p@n, {n, 0, 10000}];
short = Take[data, 10];
ListPlot[data, PlotRange -> All]
Show[Graphics@{Text @@@ Transpose@{Range@Length@short, short}}, Axes -> False]

Mod[#,1]& clips to a different square, 0 < x < 1, 0 < y < 1, and it adds numerical instability -- hence the fractal nature. 1139 points isn't enough to see any fractal behavior at all, but the next two points leave the original ellipse behind. (They're at edges of the unit square.) Things get crazy after that.

Clear[p];
a = 0.41209;
m = {{-a, -1}, {1, 0}};
p[n_] := p[n] = Mod[m.p[n - 1], 1]
p[0] = {0.7, .65};
data = Table[p@n, {n, 0, 10000}];
limited = Take[data, 1139];
firstFractal = Take[data, 1141];
ListPlot[data, PlotRange -> All]
ListPlot[limited, PlotRange -> All]
ListPlot[firstFractal, PlotRange -> All]

Points outside the initial ellipse are mapped into OTHER ellipses, but they keep escaping to start new ones.

When and if numerical instability occurs depends on p[0] as well as a. Changing a to 0.41208 destroys it, at least in the first 175000 entries (I checked). But change p[0] slightly, and it's back:

Clear[p];
a = 0.41208;
m = {{-a, -1}, {1, 0}};
p[n_] := p[n] = Mod[m.p[n - 1], 1]
p[0] = {0.71, .65};
data = Table[p@n, {n, 0, 10000}];
ListPlot[data, PlotRange -> All]

Fitting data to an ellipse isn't too hard:

Clear[a, b, c, p];
m = {{-0.41209, -1}, {1, 0}};
p[n_] := p[n] = m . p[n - 1]
p[0] = {0.7, 0.65};
data = Table[p[n], {n, 0,
      1000}];
onePtError[a_, b_, c_][
    {x_, y_}] :=
   (a*x^2 + b*x*y + c*y^2 - 1)^2
error[a_, b_, c_] :=
   Total[onePtError[a, b, c] /@
     data]
NMinimize[error[a, b, c],
   {a, b, c}]
b/a /. Last[%]
{3.0907663843655494*^-28,
   {a -> 0.9090901239676206,
    b -> 0.3746269491858168,
    c -> 0.9090901239676206}}
0.41209000000000007

b/a is the upper-left entry of the m matrix (the original a0), and c == a, so we can do it this way:

onePtError[a_, b_][{x_, y_}] := (x^2 + b*x*y + y^2 - a)^2
error[a_, b_] := Total[onePtError[a, b] /@ data]
NMinimize[error[a, b], {a, b}]
{3.8495156975525818*^-28, {a -> 1.1000009499999996, b -> 0.41209}}

We can solve for the ellipse's RHS this way, too:

x^2 + b*x*y + y^2 == a /. Thread[{x, y} -> p[0]] /. b -> 0.41209
1.10000095 == a

When FractionalPart is used to define p[n], NMinimize still works perfectly if transient points at the beginning of the series are omitted from the error definition.

Bobby

On Sat, 20 Nov 2004 03:41:36 -0500 (EST), Peter Valko <p-valko at tamu.edu> wrote:

> Roger,
> Can you tell me why is it that in the following code of yours:
>
> Clear[x, y, n];
> a0 = 0.41209;
> x[n_] := x[n] = Mod[-a0*x[n - 1] - y[n - 1], 1];
> y[n_] := y[n] = Mod[x[n - 1], 1] ;
> x[0] = 0.7;
> y[0] = .65;
> a = Table[{x[n], y[n]}, {n, 0, 10000}];
> ListPlot[a, PlotRange -> All] ;
>
> we get a fractal-like pic, but changing to a0 = 0.41208 we do not?
> (And why is that replacing Mod[-,1] by FractionalPart[-] in the above
> code will not give the same phenomenon?
>
> Peter
>
>
> Roger Bagula <tftn at earthlink.net> wrote in message news:<cneuml$qbc$1 at smc.vnet.net>...
>> I've done a lot of searches on chaos
>>  and Mathematica and have never seem this.
>> It is sensative chaos , in both the angle based a0 and the
>> initial starting point.
>> The circle was used by Chua as a starting point in his lectures on Chaos.
>>
>> Clear[x,y,a,b,s,g,a0]
>> (* circle map: from  Chaos in Digital Filters ,Chua,Lin,
>>   IEEE transactions on Circuits and Systems,vol 35 no 6 June 1988*)
>>   (* very sensitive to intial conditions*)
>> a0=Cos[Pi/6]/2;
>> x[n_]:=x[n]=Mod[-a0*x[n-1]-y[n-1],1]
>> y[n_]:=y[n]=Mod[x[n-1],1]
>> x[0]=0.7;y[0]=.65;
>> a=Table[{x[n],y[n]},{n,0, 10000}];
>> ListPlot[a, PlotRange->All]
>>
>> Respectfully, Roger L. Bagula
>>
>> tftn at earthlink.net, 11759Waterhill Road, Lakeside,Ca 92040-2905,tel: 619-5610814 :
>> alternative email: rlbtftn at netscape.net
>> URL :  http://home.earthlink.net/~tftn
>
>
>
>



-- 
DrBob at bigfoot.com
www.eclecticdreams.net


  • Prev by Date: Re: the circle map
  • Next by Date: limits
  • Previous by thread: Re: the circle map
  • Next by thread: Re: the circle map