(***********************************************************************

                    Mathematica-Compatible Notebook

This notebook can be used on any computer system with Mathematica 4.0,
MathReader 4.0, or any compatible application. The data for the notebook 
starts with the line containing stars above.

To get the notebook into a Mathematica-compatible application, do one of 
the following:

* Save the data starting with the line of stars above into a file
  with a name ending in .nb, then open the file inside the application;

* Copy the data starting with the line of stars above to the
  clipboard, then use the Paste menu command inside the application.

Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode.  Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).

NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing the 
word CacheID, otherwise Mathematica-compatible applications may try to 
use invalid cache data.

For more information on notebooks and Mathematica-compatible 
applications, contact Wolfram Research:
  web: http://www.wolfram.com
  email: info@wolfram.com
  phone: +1-217-398-0700 (U.S.)

Notebook reader applications are available free of charge from 
Wolfram Research.
***********************************************************************)

(*CacheID: 232*)


(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[     14467,        414]*)
(*NotebookOutlinePosition[     15120,        437]*)
(*  CellTagsIndexPosition[     15076,        433]*)
(*WindowFrame->Normal*)



Notebook[{

Cell[CellGroupData[{
Cell["Remapping the image", "Section"],

Cell["\<\
The image as taken by the camera was intended to be looked at on a regular \
rectangular grid.  However, it proves quite fun to play tricks with that.\
\>", "Text"],

Cell[BoxData[
    \(\(data\  = \ ImageData[ball];\)\)], "Input"],

Cell["\<\
Unless you want to map the colors (that makes another interesting set of \
experiments), you need to separate the colors from the data.\
\>", "Text"],

Cell[BoxData[
    \(GetColors[data_]\  := \ 
      Map[Apply[RGBColor, #] &, data/256. , \ {2}]\)], "Input"],

Cell["\<\
This function gets the(x,y) coordinates which each pixel represents when you \
Show the image, assuming that x and y are between xmin and xmax and ymin and \
ymax, respectively.\
\>", "Text"],

Cell[BoxData[
    \(GetPoints[
        data_, {xmin_, \ xmax_}, \ {ymin_, \ 
          ymax_}]\  := \ \[IndentingNewLine]Module[{dx, 
          dy}, \[IndentingNewLine]{dx, dy}\  = \ 
          Take[Dimensions[data], 
            2]; \[IndentingNewLine]dx\  = \ \((xmax\  - \ xmin)\)/\((dx\  - 
                1)\)\ ; \ 
        dy\  = \ \((ymax\  - \ ymin)\)/\((dy\  - \ 
                1)\); \[IndentingNewLine]Developer`ToPackedArray[
          N[Table[{x, y}, {x, xmin, xmax, dx}, {y, ymin, ymax, \ 
                dy}]]]]\)], "Input"],

Cell["\<\
To actually see the image, you need to recombine the points and colors into \
\"pixels\"\
\>", "Text"],

Cell[BoxData[
    \(Pixelize[colors_, \ 
        points_]\  := \ \[IndentingNewLine]MapThread[{#1, \ 
            Point[#2]} &, \ {colors, \ points}, \ 2]\)], "Input"],

Cell["\<\
Theo:  Is there a way to make the default of point a little square instead of \
a circle?  Would look better here for small images, though once things get \
distorted, circles are probably better.\
\>", "Text"],

Cell[BoxData[
    \(\(Show[
        Graphics[
          Pixelize[GetColors[data], \ GetPoints[data, {0, 1}, {0, 1}]]], \ 
        AspectRatio \[Rule] Automatic];\)\)], "Input"],

Cell["\<\
An alternative way to do this is to consider each pixel a square with four \
corners:\
\>", "Text"],

Cell[BoxData[
    \(GetCorners[
        data_, {xmin_, \ xmax_}, \ {ymin_, \ 
          ymax_}]\  := \ \[IndentingNewLine]Module[{dx, \ 
          dy}, \[IndentingNewLine]{dx, dy}\  = \ 
          Take[Dimensions[data], 
            2]; \[IndentingNewLine]dx\  = \ \((xmax\  - \ xmin)\)/dx; \ 
        dy\  = \ \((ymax\  - \ ymin)\)/
            dy; \[IndentingNewLine]Developer`ToPackedArray[
          N[Table[{x, y}, {x, xmin, xmax, dx}, {y, ymin, ymax, \ 
                dy}]]]]\)], "Input"],

Cell["\<\
Then to produce the image, you make polygons out of each of the squares.\
\>", "Text"],

Cell[BoxData[
    \(Polygonalize[colors_, \ corners_]\  := \ 
      Module[{p}, \[IndentingNewLine]p\  = \ 
          Partition[
            corners, {2, 2}, {1, 1}]; \[IndentingNewLine]MapThread[{#1, \ 
              Polygon[Join[#2[\([1]\)], \ 
                  Reverse[#2[\([2]\)]]]]} &, {colors, \ p}, 2]]\)], "Input"],

Cell["\<\
Either way, if you have enough points, the graphs look basically the same.  \
\
\>", "Text"],

Cell[BoxData[
    \(\(Show[
        Graphics[
          Polygonalize[GetColors[data], \ 
            GetCorners[data, {0, 1}, {0, 1}]]], \ 
        AspectRatio \[Rule] Automatic];\)\)], "Input"],

Cell["After distortions, things will be a bit different.", "Text"],

Cell["\<\
An interesting example for a map is a discretization for a differential \
equation.  This is a first order discretization of Duffing's equation (for \
the motion of the end of a clamped beam).\
\>", "Text"],

Cell[BoxData[
    \(duffingmap[h_]\  := \ 
      Function[{v}, \ 
        Block[{x\  = \ v[\([1]\)], \ 
            y\  = \ v[\([2]\)]}, \[IndentingNewLine]x\  = \ 
            x\  + \ h\ y; \[IndentingNewLine]{x, 
            y\  + \ h\ \((x\  - \ x^3/25. )\)}]]\)], "Input"],

Cell["\<\
It is NOT Euler's method, since the updated value of x is used to compute the \
update for y.  This map has the property that it conserves to machine \
precision an approximate energy which is near (to order h) the actual energy \
of the system.  Euler's method would be defined by\
\>", "Text"],

Cell[BoxData[
    \(duffingeuler[h_]\  := \ 
      Function[{v}, \ 
        Block[{x\  = \ v[\([1]\)], \ 
            y\  = \ v[\([2]\)]}, \[IndentingNewLine]{x\  + \ h\ y, \ 
            y\  + \ h\ \((x\  - \ x^3/25. )\)}]]\)], "Input"],

Cell["\<\
It is upon iteration, that this map becomes interesting.  For example, by \
looking at a few initial conditions, we can get an idea of the action of the \
map\
\>", "Text"],

Cell[BoxData[{
    \(\(h\  = \  .1;\)\), "\[IndentingNewLine]", 
    \(\(Show[
        Block[{$DisplayFunction\  = \ Identity}, {\[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {2. , 0. }, 100], \ 
              PlotStyle \[Rule] 
                RGBColor[0, 0, 1]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {\(-2. \), 0. }, 100], 
              PlotStyle \[Rule] 
                RGBColor[0, 0, 1]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {0. , 1. }, 200], \ 
              PlotStyle \[Rule] 
                RGBColor[1, 0, 0]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {0. , 2. }, 400], \ 
              PlotStyle \[Rule] 
                RGBColor[1, 0, 0]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {0. ,  .05}, 
                1000]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {4.5, 0}, 100], \ 
              PlotStyle \[Rule] 
                RGBColor[0, 1, 0]], \[IndentingNewLine]ListPlot[
              NestList[duffingmap[h], {\(-4.5\), 0}, 100], \ 
              PlotStyle \[Rule] RGBColor[0, 1, 0]]}], \ 
        PlotRange \[Rule] {{\(-8\), 8}, {\(-4\), 4}}];\)\)}], "Input"],

Cell["\<\
The idea is to apply this map to all of the points used to represent the \
picture.  So that you can understand what is happening, lets consider a very \
small subset of the picture (the puppy's eye, for example)\
\>", "Text"],

Cell[BoxData[
    \(\(corners\  = \ 
        GetCorners[
          Take[data, {50, 55}, {50, 55}], \ {\(-1\), 1}, {\(-1\), 
            1}];\)\)], "Input"],

Cell["Now use the map (and Map)  to get to the new points", "Text"],

Cell[BoxData[
    \(\(newcorners\  = \ Map[duffingmap[ .1], \ corners, \ {2}];\)\)], "Input"],

Cell["\<\
This shows the original points in black and the new in blue, with a red line \
connecting corresponding points.\
\>", "Text"],

Cell[BoxData[
    \(\(Show[
        Block[{$DisplayFunction\  = \ Identity, \ p, \ 
            np}, \[IndentingNewLine]p\  = \ 
            Flatten[corners, \ 1]; \[IndentingNewLine]np\  = \ 
            Flatten[newcorners, 1]; \[IndentingNewLine]{ListPlot[p, \ 
              PlotStyle \[Rule] PointSize[0.02]], \ 
            ListPlot[np, \ 
              PlotStyle \[Rule] {PointSize[0.02], 
                  RGBColor[0, 0, 1]}], \[IndentingNewLine]Graphics[{RGBColor[
                  1, 0, 0], 
                MapThread[Line[{##}] &, {p, np}]}]}]];\)\)], "Input"],

Cell["\<\
Of course, we could iterate the map several times to get a greater effect:\
\>", "Text"],

Cell[BoxData[
    \(\(newcorners\  = \ 
        Nest[Map[duffingmap[ .1], #, {2}] &, \ corners, 5];\)\)], "Input"],

Cell[BoxData[
    \(\(Show[
        Block[{$DisplayFunction\  = \ Identity, \ p, \ 
            np}, \[IndentingNewLine]p\  = \ 
            Flatten[corners, \ 1]; \[IndentingNewLine]np\  = \ 
            Flatten[newcorners, 1]; \[IndentingNewLine]{ListPlot[p, \ 
              PlotStyle \[Rule] PointSize[0.02]], \ 
            ListPlot[np, \ 
              PlotStyle \[Rule] {PointSize[0.02], 
                  RGBColor[0, 0, 1]}], \[IndentingNewLine]Graphics[{RGBColor[
                  1, 0, 0], 
                MapThread[Line[{##}] &, {p, np}]}]}]];\)\)], "Input"],

Cell[BoxData[
    \(\(Show[
        Graphics[
          Polygonalize[GetColors[Take[data, {50, 55}, {50, 55}]], \ 
            newcorners]], \ %, \ AspectRatio \[Rule] Automatic];\)\)], "Input"],

Cell["\<\
Now catch a breath and then lets try it with the full picture:\
\>", "Text"],

Cell[BoxData[{
    \(\(corners\  = \ 
        GetCorners[
          data, {\(-1. \), 1. }, {\(-1. \), 
            1. }];\)\), "\[IndentingNewLine]", 
    \(\(colors\  = \ GetColors[data];\)\)}], "Input"],

Cell["Now map the points ...", "Text"],

Cell[BoxData[
    \(\(newcorners\  = \ 
        Nest[Map[duffingmap[ .1], #, {2}] &, \ corners, 5];\)\)], "Input"],

Cell["... and voila, the picture has been modified accordingly.  ", "Text"],

Cell[BoxData[
    \(\(Show[Graphics[Polygonalize[colors, newcorners]], \ 
        AspectRatio \[Rule] Automatic];\)\)], "Input"],

Cell["\<\
Now we are ready to try several iterations.  You can certainly do this with \
things as they are, but it can be slow.  This is a good example of the sort \
of command which can be sped up quite a bit by compiling it for machine \
numbers.\
\>", "Text"],

Cell["\<\
This defines a function which will take a given map and compile another \
function which will apply it to either the corners or the points related to \
the image.\
\>", "Text"],

Cell[BoxData[
    \(CompileIterations[map_]\  := \ 
      Compile[{{p, \ _Real, \ 3}, {n, \ _Integer}}, 
        Nest[Map[map, \ #, \ {2}] &, \ p, \ n]]\)], "Input"],

Cell[BoxData[
    \(cit\  = \ CompileIterations[duffingmap[ .1]]\)], "Input"],

Cell["\<\
This works just like a function which takes two arguments and it will apply \
the map to the points in the first argument the second argument number of \
times.  As a comparison, here are timings for 2 iterations\
\>", "Text"],

Cell[BoxData[
    \(Timing[\(cit[corners, 2];\)]\)], "Input"],

Cell[BoxData[
    \(Timing[\(Nest[Function[Map[duffingmap[ .1], \ #, \ {2}]], corners, 
          2];\)]\)], "Input"],

Cell["So the CompiledFunction is about 10 times faster.", "Text"],

Cell["\<\
However, if we do more than a few iterations, we start to get into trouble \
with the polygons!\
\>", "Text"],

Cell[BoxData[{
    \(\(newcorners\  = cit[\ corners, 100];\)\), "\[IndentingNewLine]", 
    \(\(Show[Graphics[Polygonalize[colors, newcorners]], \ 
        AspectRatio \[Rule] Automatic];\)\)}], "Input"],

Cell["So for bigger distortions, you are better off with points", "Text"],

Cell[BoxData[{
    \(\(points\  = \ 
        GetPoints[
          data, {\(-1\), 1}, {\(-1\), 1}];\)\), "\[IndentingNewLine]", 
    \(\(newpoints\  = \ cit[points, \ 100];\)\), "\[IndentingNewLine]", 
    \(\(Show[Graphics[Pixelize[colors, newpoints]], \ 
        AspectRatio \[Rule] Automatic];\)\)}], "Input"],

Cell["\<\
Near one of the stable equilibria, things don't get quite so distorted:\
\>", "Text"],

Cell[BoxData[{
    \(\(points\  = \ 
        GetPoints[data, {4, 6}, {\(-1\), 1}];\)\), "\[IndentingNewLine]", 
    \(\(Show[
        GraphicsArray[
          Block[{$DisplayFunction\  = \ Identity, \ pl}, 
            Table[pl\  = \ 
                Show[Graphics[Pixelize[colors, points]], \ 
                  AspectRatio \[Rule] Automatic, \ 
                  PlotRange \[Rule] {{3, 7}, {\(-2\), 2}}]; 
              points\  = \ cit[points, 100]; pl, {2}]]]];\)\)}], "Input"],

Cell["\<\
However, it is most interesting to contrast this with what Euler's method \
does...\
\>", "Text"],

Cell[BoxData[""], "Input"],

Cell[BoxData[{
    \(\(eit\  = \ 
        CompileIterations[duffingeuler[ .1]];\)\), "\[IndentingNewLine]", 
    \(\(points\  = \ 
        GetPoints[data, {4, 6}, {\(-1\), 1}];\)\), "\[IndentingNewLine]", 
    \(\(Show[
        GraphicsArray[
          Block[{$DisplayFunction\  = \ Identity, \ pl}, 
            Table[pl\  = \ 
                Show[Graphics[Pixelize[colors, points]], \ 
                  AspectRatio \[Rule] Automatic, \ 
                  PlotRange \[Rule] {{3, 7}, {\(-2\), 2}}]; 
              points\  = \ eit[points, 100]; pl, {2}]]]];\)\)}], "Input"],

Cell["Funny what can happen when you don't conserve energy!!!!", "Text"],

Cell["\<\
The other map, however, has some intriguing features of its own.  For example \
with a larger time step...\
\>", "Text"],

Cell[BoxData[
    \(\(cit\  = \ CompileIterations[duffingmap[ .7]];\)\)], "Input"],

Cell[BoxData[
    \(points\  = \ GetPoints[data, {\(-2\), 2}, {\(-2\), 2}]; 
    Do[Show[Graphics[Pixelize[colors, points]], \ 
        AspectRatio \[Rule] Automatic, \ 
        PlotRange \[Rule] {{\(-8\), 8}, {\(-4\), 4}}]; 
      points\  = \ cit[points, 1], {10}]\)], "Input"],

Cell["\<\
What you are seeing here is (part of ) the heteroclinic tangle at the \
crossing of the separatrices.\
\>", "Text"]
}, Open  ]]
},
FrontEndVersion->"4.0 for Microsoft Windows",
ScreenRectangle->{{0, 1024}, {0, 695}},
WindowSize->{919, 664},
WindowMargins->{{5, Automatic}, {Automatic, 0}}
]


(***********************************************************************
Cached data follows.  If you edit this Notebook file directly, not using
Mathematica, you must remove the line containing CacheID at the top of 
the file.  The cache data will then be recreated when you save this file 
from within Mathematica.
***********************************************************************)

(*CellTagsOutline
CellTagsIndex->{}
*)

(*CellTagsIndex
CellTagsIndex->{}
*)

(*NotebookFileOutline
Notebook[{

Cell[CellGroupData[{
Cell[1739, 51, 38, 0, 53, "Section"],
Cell[1780, 53, 174, 3, 33, "Text"],
Cell[1957, 58, 64, 1, 30, "Input"],
Cell[2024, 61, 159, 3, 33, "Text"],
Cell[2186, 66, 108, 2, 30, "Input"],
Cell[2297, 70, 201, 4, 33, "Text"],
Cell[2501, 76, 542, 11, 110, "Input"],
Cell[3046, 89, 112, 3, 33, "Text"],
Cell[3161, 94, 167, 3, 50, "Input"],
Cell[3331, 99, 220, 4, 52, "Text"],
Cell[3554, 105, 176, 4, 30, "Input"],
Cell[3733, 111, 109, 3, 33, "Text"],
Cell[3845, 116, 496, 10, 110, "Input"],
Cell[4344, 128, 96, 2, 33, "Text"],
Cell[4443, 132, 323, 6, 70, "Input"],
Cell[4769, 140, 102, 3, 33, "Text"],
Cell[4874, 145, 194, 5, 30, "Input"],
Cell[5071, 152, 66, 0, 33, "Text"],
Cell[5140, 154, 216, 4, 52, "Text"],
Cell[5359, 160, 276, 6, 70, "Input"],
Cell[5638, 168, 305, 5, 52, "Text"],
Cell[5946, 175, 237, 5, 50, "Input"],
Cell[6186, 182, 182, 4, 33, "Text"],
Cell[6371, 188, 1238, 23, 190, "Input"],
Cell[7612, 213, 236, 4, 52, "Text"],
Cell[7851, 219, 155, 4, 30, "Input"],
Cell[8009, 225, 67, 0, 33, "Text"],
Cell[8079, 227, 93, 1, 30, "Input"],
Cell[8175, 230, 135, 3, 33, "Text"],
Cell[8313, 235, 572, 11, 110, "Input"],
Cell[8888, 248, 98, 2, 33, "Text"],
Cell[8989, 252, 114, 2, 30, "Input"],
Cell[9106, 256, 572, 11, 110, "Input"],
Cell[9681, 269, 194, 4, 30, "Input"],
Cell[9878, 275, 86, 2, 33, "Text"],
Cell[9967, 279, 204, 5, 50, "Input"],
Cell[10174, 286, 38, 0, 33, "Text"],
Cell[10215, 288, 114, 2, 30, "Input"],
Cell[10332, 292, 75, 0, 33, "Text"],
Cell[10410, 294, 128, 2, 30, "Input"],
Cell[10541, 298, 262, 5, 52, "Text"],
Cell[10806, 305, 186, 4, 33, "Text"],
Cell[10995, 311, 165, 3, 30, "Input"],
Cell[11163, 316, 77, 1, 30, "Input"],
Cell[11243, 319, 236, 4, 52, "Text"],
Cell[11482, 325, 61, 1, 30, "Input"],
Cell[11546, 328, 117, 2, 30, "Input"],
Cell[11666, 332, 65, 0, 33, "Text"],
Cell[11734, 334, 119, 3, 33, "Text"],
Cell[11856, 339, 203, 3, 50, "Input"],
Cell[12062, 344, 73, 0, 33, "Text"],
Cell[12138, 346, 311, 6, 70, "Input"],
Cell[12452, 354, 95, 2, 33, "Text"],
Cell[12550, 358, 481, 10, 110, "Input"],
Cell[13034, 370, 107, 3, 33, "Text"],
Cell[13144, 375, 26, 0, 30, "Input"],
Cell[13173, 377, 575, 12, 130, "Input"],
Cell[13751, 391, 72, 0, 33, "Text"],
Cell[13826, 393, 130, 3, 33, "Text"],
Cell[13959, 398, 82, 1, 30, "Input"],
Cell[14044, 401, 279, 5, 70, "Input"],
Cell[14326, 408, 125, 3, 33, "Text"]
}, Open  ]]
}
]
*)




(***********************************************************************
End of Mathematica Notebook file.
***********************************************************************)