(***********************************************************************
This file was generated automatically by the Mathematica front end.
It contains Initialization cells from a Notebook file, which typically
will have the same name as this file except ending in ".nb" instead of
".m".

This file is intended to be loaded into the Mathematica kernel using
the package loading commands Get or Needs.  Doing so is equivalent to
using the Evaluate Initialiation Cells menu command in the front end.

DO NOT EDIT THIS FILE.  This entire file is regenerated automatically 
each time the parent Notebook file is saved in the Mathematica front end.
Any changes you make to this file will be overwritten.
***********************************************************************)

\!\(BeginPackage["\<TriangleSequences`\>"]; \n\n
  T::"\<usage\>" = 
    "\<T[{a,b}] gives the result of the triangle  mapping on the point {a,b}.\
\>"; \nTriangle::usage\  = \ 
    "\<Triangle[{a,b}] gives the partition triangle in which the point {a,b} \
is located.\>"; \n
  TriangleSequence::"\<usage\>" = 
    "\<TriangleSequence[{a,b},n] returns a list containing the first n terms \
in the triangle sequence for the point {a,b}.\>"; \n
  CVectors::"\<usage\>" = 
    "\<CVectors[{a,b},n] gives you n vectors which     approach the plane \
x+ay+cz=0.\>"; \n
  PlotOrbit::"\<usage\>" = 
    "\<PlotOrbit[{a,b}] plots the orbit of the point {a,b} under the triangle \
iteration.\>"; \n
  NumberPoints::"\<usage\>" = 
    "\<NumberPoints is an option for PlotOrbit.  NumberPoints->True numbers \
the points in the orbit, beginning with 0. Default setting: \
NumberPoints->True.\>"; \n
  DisplayLines::"\<usage\>" = 
    "\<DisplayLines is an option for PlotOrbit.  DisplayLines->True draws \
lines between the points in the orbit. Default setting: \
DisplayLines->False.\>"; \n
  PartitionColor::"\<usage\>" = 
    "\<PartitionColor is an option for PlotOrbit. PartitionColor specifies \
the color used to display the partition lines of the triangle. Default \
setting: PartitionColor->RGBColor[1,.5,0].\>"; \n
  OrbitColor::"\<usage\>" = 
    "\<OrbitColor is an option for PlotOrbit.  OrbitColor specifies the color \
in which lines between orbit points will be displayed. Default setting: \
OrbitColor->RGBColor[0,.25,1].\>"; \n
  OrbitPoints::"\<usage\>" = 
    "\<OrbitPoints is an option for PlotOrbit.  OrbitPoints->n displays the \
first n points in the orbit.  Default setting: OrbitPoints->30.\>"; \n
  PartitionLines::"\<usage\>" = 
    "\<PartitionLines is an option for PlotOrbit.  PartitionLines->n displays \
the  first n-1 partition triangles.  Default setting: PartitionLines->10.\>"; 
  \nFindPeriodicPoint::"\<usage\>" = 
    "\<FindPeriodicPoint[{a1,a2,...}] finds a point with periodic triangle \
sequence {a1,a2,...,a1,a2,...}.  FindPeriodicPoint[{b1,b2,...}, {a1,a2,...}] \
finds the point with triangle sequence {b1,b2,...,a1,a2,...,a1,a2,...}.  \
NOTE: This function sometimes returns approximate rather than exact \
values.\>"; \n
  Spread::"\<usage\>" = 
    "\<Spread[{a,b}] begins with a region surrounding the point {a,b}, and \
plots the effect of the iterative triangle mapping on points inside this \
region.\>"; \n
  FindTerminatingPoint::usage\  = \ 
    "\<FindTerminatingPoint[{a1,a2,..,an},f[x]] finds a point {x,y} on the \
curve y = f[x] which has the terminating triangle sequence {a1,a2,..an}.\>"; 
  \n\nBegin["\<`Private`\>"]; \nK[a_, b_] := Floor[Simplify[\(1 - a\)\/b]]; \n
  \nTriangle[{a_, b_}] := \ K[a, b]; \n\n
  T[{a_, b_}] := FullSimplify[{b\/a, \(1 - a - K[a, b]\ b\)\/a}]; \n
  TriangleSequence[{a_, b_}, n_:  50] := 
    Module[{x, y, seq, c}, x = a; y = b; seq = {}; c = 0; 
      While[N[y] \[NotEqual] 0 && c < n, {
          seq = Append[seq, K[x, y]]; {x, y} = T[{x, y}]; c = c + 1}]; seq]; 
  \n\nCVectors[{aa_, bb_}, n_:  50] := 
    Module[{lst, c2, c1, c0, i, ak, new}, lst = {}; 
      ak = TriangleSequence[{aa, bb}, n]; c2 = {1, 0, 0}; c1 = {0, 1, 0}; 
      c0 = {0, 0, 1}; 
      For[i = 1, i \[LessEqual] Length[ak], 
        \(i++\), {
          new = c2 - c1 - ak\[LeftDoubleBracket]i\[RightDoubleBracket]\ c0, 
          lst = Append[lst, new], c2 = c1, c1 = c0, c0 = new}]; lst]; \n\n
  Options[PlotOrbit] = {NumberPoints \[Rule] True, 
      DisplayLines \[Rule] False, PartitionColor \[Rule] RGBColor[1,  .5, 0], 
      OrbitColor \[Rule] RGBColor[0,  .25, 1], OrbitPoints \[Rule] 30, 
      PartitionLines \[Rule] 10}; \n\n
  PlotOrbit[{a_, b_}, opts___] := 
    PlotOrb[{a, b}, 
      \({NumberPoints, DisplayLines, PartitionColor, OrbitColor, OrbitPoints, 
            PartitionLines} /. {opts}\) /. Options[PlotOrbit]]; \n\n
  PlotOrb[{a_, b_}, {numpts_, lines_, tcol_, ocol_, pts_, lns_}] := 
    Module[{len, triangle, points, orbit}, 
      len = Length[TriangleSequence[{a, b}, pts]]; 
      triangle = {tcol, Line[{{1, 0}, {0, 0}, {1, 1}}], 
          Table[Line[{{1, 0}, {1\/i, 1\/i}}], {i, 1, lns}]}; 
      points = NestList[T, {a, b}, len]; orbit = {ocol, Line[points]}; 
      plot = Table[
          Point[points\[LeftDoubleBracket]i\[RightDoubleBracket]], {i, 1, 
            len}]; labels = 
        Table[Text[i - 1, 
            points\[LeftDoubleBracket]i
                \[RightDoubleBracket] + { .018,  .008}], {i, 1, len}]; 
      display = {triangle, plot}; 
      If[lines, display = Append[display, orbit]]; 
      If[numpts, display = Append[display, labels]]; 
      Show[Graphics[display, Axes \[Rule] True, 
          PlotRange \[Rule] {{0, 1}, {0, 1}}, AspectRatio \[Rule] 1]]]; \n\n
  Points[x1_, x2_, y1_, y2_] := 
    Flatten[Table[{x1 + 1\/50\ \((x2 - x1)\)\ i, 
          y1 + 1\/50\ \((y2 - y1)\)\ j}, {i, 1, 10}, {j, 1, 10}], 1]; \n\n
  Spread[{a_, b_}, r_:  10] := 
    Module[{p, xch, ych, g}, 
      \((xch =  .002; ych =  .002; 
        p = Points[a - xch, a + xch, b - ych, b + ych]; \n\t\t\tg\  = \ {}; 
        g\  = \ Append[g, 
            Graphics[{RGBColor[1, 0,  .2], 
                Table[Point[p\[LeftDoubleBracket]i\[RightDoubleBracket]], {i, 
                    1, Length[p]}]}]]; 
        For[n = 0, n < r, 
          n = n + 1, {
            p = Table[
                T[p\[LeftDoubleBracket]i\[RightDoubleBracket]], {i, 1, 
                  Length[p]}]; 
            g\  = \ Append[g, 
                Graphics[{RGBColor[1, 0,  .2], 
                    Table[Point[
                        p\[LeftDoubleBracket]i\[RightDoubleBracket]], {i, 1, 
                        Length[p]}]}]]}]; \n\t\t\t
        ShowAnimation[g, Axes \[Rule] True, 
          PlotRange \[Rule] {{0, 1}, {0, 1}}, AspectRatio \[Rule] 1])\)]; \n\n
  P[k_] := {{0, 0, 1}, {1, 0, \(-1\)}, {0, 1, \(-k\)}}; \n\n
  InTriangle[{a_, b_}] := 
    If[Im[N[a]] \[NotEqual] 0 || Im[N[b]] \[NotEqual] 0, False, 
      If[N[a] > 1 || N[b] > N[a] || N[b] < 0, False, True]]; \n\n
  AlmostInTriangle[{a_, b_}] := \ InTriangle[N[{a, b}]]; \n\n
  FindPeriodicPoint[per_] := 
    Module[{m, sol}, m = IdentityMatrix[3]; 
      Do[m = m . P[per\[LeftDoubleBracket]i\[RightDoubleBracket]], {i, 1, 
          Length[per]}]; 
      sol = Solve[{1, x, y} . m == l\ {1, x, y}, {x, y}, l]; 
      el = Table[{x, y} /. sol\[LeftDoubleBracket]i\[RightDoubleBracket], {i, 
            1, Length[sol]}]; 
      ans\  = \ Simplify[Flatten[Select[el, InTriangle]]]; \n\t
      If[ans\  != \ {}, \ ans, \n\t\t\t
        Flatten[Select[Re[N[el]], AlmostInTriangle]]]]; \n\n
  S[{x_, y_}, k_] := {1\/\(1 + k\ x + y\), x\/\(1 + k\ x + y\)}; \n\n
  FindPeriodicPoint[beg_, per_] := 
    Module[{pt}, pt = FindPeriodicPoint[per]; 
      Simplify[Fold[S, pt, 
          Table[beg\[LeftDoubleBracket]\(-i\)\[RightDoubleBracket], {i, 1, 
              Length[beg]}]]]]; \n\n
  InvPoint[x_, y_, klist_] := 
    Fold[S, {x, y}, Table[klist[\([\(-i\)]\)], {i, 1, Length[klist]}]]; \n\n
  PartitionPts[klist_] := {InvPoint[0, 0, klist], InvPoint[1, 0, klist]}; \n\n
  FindTerminatingPoint[klist_, f_] := \n\t\t
    Module[{p, m, inLine, sol, pts}, 
      \((p = PartitionPts[klist]; 
        m = \((p[\([1, 2]\)] - p[\([2, 2]\)])\)/
            \((p[\([1, 1]\)] - p[\([2, 1]\)])\); 
        inLine[{x_, y_}] := 
          Evaluate[
            \((p[\([1, 1]\)] <= N[x] <= p[\([2, 1]\)])\) || 
              \((p[\([1, 1]\)] >= N[x] >= p[\([2, 1]\)])\)]; 
        sol = Solve[{
              Simplify[y - p[\([1, 2]\)] == m\ \((x - p[\([1, 1]\)])\)], 
              y == f}, {x, y}]; \n\t\t\t
        pts\  = \ Table[{x, y} /. sol[\([i]\)], {i, 1, Length[sol]}]; \n\t\t\t
        Flatten[Select[pts, inLine]])\)]; \n\nEnd[]\n
  EndPackage[]\)