(************** Content-type: application/mathematica **************
                     CreatedBy='Mathematica 5.1'

                    Mathematica-Compatible Notebook

This notebook can be used with any Mathematica-compatible
application, such as Mathematica, MathReader or Publicon. 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[     13323,        313]*)
(*NotebookOutlinePosition[     14023,        337]*)
(*  CellTagsIndexPosition[     13979,        333]*)
(*WindowFrame->Normal*)



Notebook[{

Cell[CellGroupData[{
Cell["Hermite denominators", "Title"],

Cell["\<\
Calculating the sequence of Hermite denominators for a rational number\
\>", "Subsubtitle"],

Cell[TextData[{
  "This notebook provides code that computes the sequence of Hermite \
denominators of a rational input between 0 and 1/2. The ",
  StyleBox["Mathematica",
    FontSlant->"Italic"],
  " function \"hermite\" takes as input a pair {p,q} of positive integers \
with 2p<q. It returns a list with 3 entries, each of which is again a list. \n\
\nThe first of these lists is the list of intermediate matrices generated by \
the Hermite algorithm as it reduces the input to a fraction of the form 1/k. \
The second is the list of fractions it generates during this. The third list \
is the list of Hermite denominators for the original input. These \
denominators are good denominators, about as good as those produced by the \
classical or the centered continued fraction algorithms. The list produced by \
the Hermite algorithm is not necessarily identical to that produced by either \
of these algorithms, though. As we shall see, there are examples in which the \
Hermite algorithm skips past some of the denominators generated by the \
classical algorithm, and examples in which it picks up denominators that are \
skipped by the centered continued fraction algorithm. \n\nIn more detail, the \
first object returned is a list of 2 by 2 matrices with rational entries in \
the first column and integer entries in the second column, each of \
determinant 1 or -1, constructed during the lattice reduction process. These \
matrices are each a basis of the same lattice, the one given by initial \
lattice basis {{1,0},{-p/q,1}}.\n\nThe standard Euclidean norm of a vector \
(x,y) is of course Sqrt[x^2+y^2], and a basis {u,v} for a two-dimensional \
lattice is a reduced basis (in two dimensions all versions of \"reduced\" \
coincide) if u.u\[LessEqual]v.v and Abs[u.v]\[LessEqual](1/2)Abs[u.u]. \n\nBy \
this measure, the initial basis, {{1,0},{-p/q,1}} is reduced. But when we \
change the inner product, redefining u.v to mean u[[1]]v[[1]]+c u[[2]]v[[2]], \
where 0<c<1, as c decreases, a reduced basis can cease to be reduced, and \
computations are required to effect the reduction of a basis with respect to \
the new inner product. Thus, when c=105/121, the norms of our first two \
vectors in the example below, with p=4 and q=11, are equal, and when c drops \
below 105/121, the formerly longer vector {-4/11,1} becomes shorter than the \
formerly shortest vector {1,0}. Thus, the previously reduced basis \
{{-4/11,1},{1,0}} for our lattice ceases to be reduced. Switching the two \
vectors restores it to the status of reduced. This is the first basis we \
construct, and the first entry in the list of bases returned. \n\nAs c drops \
below 7/363, this basis too ceases to be reduced. Now with u={-4/11,1} and \
v={1,0}, we need v'=u and u'=u+2 v to have a reduced basis. \n\nThe sequence \
of bases that result, as we slide c down toward zero, is the first entry \
returned by \"hermite\". \n\nThe second list returned is the list of positive \
rational numbers that effectively constitute the fraction that the hermite \
algorithm is currently working on. This fraction becomes steadily simpler, \
until finally it has numerator 1. At that point, just as with the standard \
continued fraction algorithm, we are pretty much done. \n\nThe entries of the \
second list can be extracted by inspection of the first column in the list of \
reduced lattice bases. Thus, in the third of our bases, the first entries of \
u and v are respectively 3/11 and -4/11, and Abs[(3/11)/(-4/11)]=3/4. \n\nThe \
third list is the sequence of integers that serve as denominators for Hermite \
approximations to the original target of approximation. These are also the \
second, integer, entries in the first vector of each of the successive \
reduced bases. Thus, in our example, the Hermite approximations to 4/11 are \
0/1,1/2, and 1/3, with finally 4/11 itself as the unlisted \"approximation\". \
The classical continued fraction expansion of 4/11 is 1/(2+1/(1+1/3)), so \
that the approximations computed by this means are the same: 1/2, 1/3, and \
finally, 4/11 itself. But the Hermite procedure does not always return the \
same list. Thus, the Hermite denominators for 47/129 are 1, 2, 3, and 11, \
while the classical continued fraction convergents are 1/2, 1/3, 3/8, and \
4/11. \n\nIt should be kept in mind that while the Hermite denominators give \
good approximations, and while the work of computing them is not all that \
taxing, it is still more trouble than getting the classical continued \
fraction expansion. The real reason we study this Hermite algorithm is to \
glean from it insight into how we might deal with higher-dimensional problems \
in simultaneous diophantine approximation. The Jacobi-Perron algorithm is, in \
effect, a direct translation to higher dimensions of the classical continued \
fraction algorithm, but in higher dimensions, its output is no longer \
best-possible. Indeed, it can be quite far from the mark, compared to what \
one gets with, say, the Lagarias algorithm. \n"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[
    \(hermite[{4, 11}]\)], "Input"],

Cell[BoxData[
    \({{{{\(-\(4\/11\)\), 1}, {1, 0}}, {{3\/11, 2}, {\(-\(4\/11\)\), 
            1}}, {{\(-\(1\/11\)\), 3}, {3\/11, 2}}}, {4\/11, 3\/4, 1\/3}, {1, 
        2, 3}}\)], "Output"]
}, Open  ]],

Cell[BoxData[
    \(nearest[u_] := Floor[u + 1/2]\)], "Input",
  InitializationCell->True],

Cell[BoxData[
    \(remnant[u_] := u - nearest[u]\)], "Input",
  InitializationCell->True],

Cell["\<\
These give the nearest integer and the offset from it, respectively, for a \
rational input u. \
\>", "Text"],

Cell[BoxData[
    \(\(\(centeredcfx[u_] := 
      Module[{qlist, v, a, p, q, r, s, mlist}, qlist = {1}; mlist = {}; 
        p = nearest[u]; q = 1; r = 1; s = 0; v = remnant[u]; 
        While[v \[NotEqual] 0, v = 1/v; a = nearest[v]; 
          v = remnant[
              v]; {{r, p}, {s, 
                q}} = {{r, p}, {s, q}} . {{0, 1}, {1, 
                  a}}; \[IndentingNewLine]AppendTo[qlist, q]; 
          AppendTo[mlist, {{r, p}, {s, q}}]]; qlist]\)\(\[IndentingNewLine]\)
    \)\)], "Input",
  InitializationCell->True],

Cell["\<\
This gives the centered continued fraction expansion of a rational number u, \
so as to provide another basis for comparison of the Hermite algorithm's \
outputs to those of familiar continued fraction algorithms. \
\>", "Text"],

Cell[BoxData[
    \(\(\(switchtrigger[m_] := 
      Module[{m11, m12, m21, 
          m22}, \[IndentingNewLine]{{m11, m12}, {m21, m22}} = 
          m; \[IndentingNewLine]If[Abs[m12] \[GreaterEqual] Abs[m22], 
          Return[0]]; \[IndentingNewLine]If[Abs[m11] \[LessEqual] Abs[m21], 
          Return[\(-1\)]]; \[IndentingNewLine]\((m11^2 - m21^2)\)/\((m22^2 - 
              m12^2)\)]\)\(\[IndentingNewLine]\)
    \)\)], "Input",
  InitializationCell->True],

Cell["\<\
And now we hit some technical code implementing one and another fragment of \
the Hermite algorithm in one dimension. \
\>", "Text"],

Cell[BoxData[
    \(addtrigger[m_] := 
      Module[{m11, m12, m21, 
          m22}, \[IndentingNewLine]{{m11, m12}, {m21, m22}} = 
          m; \[IndentingNewLine]If[Abs[m21 + m11] \[GreaterEqual] Abs[m21], 
          Return[0]]; \[IndentingNewLine]\((m21^2 - \((m11 + 
                    m21)\)^2)\)/\((\((m12 + m22)\)^2 - m22^2)\)]\)], "Input",
  InitializationCell->True],

Cell[BoxData[
    \(subtrigger[m_] := 
      Module[{m11, m12, m21, m22}, {{m11, m12}, {m21, m22}} = 
          m; \[IndentingNewLine]If[Abs[m11 - m21] \[GreaterEqual] Abs[m21], 
          Return[0]]; \[IndentingNewLine]\((m21^2 - \((m11 - 
                    m21)\)^2)\)/\((\((m12 - m22)\)^2 - m22^2)\)]\)], "Input",
  InitializationCell->True],

Cell[BoxData[
    \(move[m_] := 
      Module[{t1, t2, t3, t}, \[IndentingNewLine]t1 = switchtrigger[m]; 
        t2 = addtrigger[m]; t3 = subtrigger[m]; 
        t = Max[{t1, t2, t3}]; \[IndentingNewLine]If[t \[Equal] t1, 
          Return[{{0, 1}, {1, 0}} . m]]; \[IndentingNewLine]If[t \[Equal] t2, 
          Return[{{1, 0}, {1, 1}} . m]]; \[IndentingNewLine]If[t \[Equal] t3, 
          Return[{{1, 0}, {\(-1\), 1}} . m]]]\)], "Input",
  InitializationCell->True],

Cell[BoxData[
    \(hermite[u_] := 
      Module[{p, q, m, mlist, qlist, ratiolist, qq, oldq, ratio}, 
        If[u \[LessEqual] 0 || u \[GreaterEqual] 1/2 || 
            Head[u] \[NotEqual] Rational, 
          Print["\<Invalid input. Algorithm needs a rational number strictly \
between 0 and 1/2.\>"]; Return[False]]; \[IndentingNewLine]p = Numerator[u]; 
        q = Denominator[u]; \[IndentingNewLine]m = {{1, 0}, {\(-p\)/q, 1}}; 
        mlist = {m}; \[IndentingNewLine]qlist = {0}; qq = 0; ratiolist = {}; 
        oldq = 0; \[IndentingNewLine]While[
          m[\([2, 1]\)] \[NotEqual] 0, \[IndentingNewLine]m = 
            move[m]; \[IndentingNewLine]qq = 
            Abs[m[\([1, 2]\)]]; \[IndentingNewLine]If[qq > oldq, 
            ratio = \(-m[\([1, 1]\)]\)/
                m[\([2, 1]\)]; \[IndentingNewLine]AppendTo[ratiolist, 
              ratio]; \[IndentingNewLine]oldq = 
              qq; \[IndentingNewLine]AppendTo[mlist, 
              m]; \[IndentingNewLine]AppendTo[qlist, qq]]]; {Drop[mlist, 1], 
          ratiolist, Drop[qlist, 1]}]\)], "Input",
  InitializationCell->True],

Cell["\<\

The first list gives the successive matrices ``m'' that the hermite algorithm \
generates. The second gives the successive, and successively simpler, ratios \
it produces. The third gives the sequence of denominators q such that there \
exists an integer p for which p/q is a Hermite approximation to the input,  \
but with a `0' entry prepended. 

 \
\>", "Text"],

Cell[CellGroupData[{

Cell[BoxData[
    \(hermite[47/129]\)], "Input"],

Cell[BoxData[
    \({{{{\(-\(47\/129\)\), 1}, {1, 0}}, {{35\/129, 2}, {\(-\(47\/129\)\), 
            1}}, {{\(-\(4\/43\)\), 3}, {35\/129, 2}}, {{\(-\(1\/129\)\), 
            11}, {\(-\(4\/43\)\), 3}}}, {47\/129, 35\/47, 
        12\/35, \(-\(1\/12\)\)}, {1, 2, 3, 11}}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[BoxData[
    \(hermite[4831/29334]\)], "Input"],

Cell[BoxData[
    \({{{{\(-\(4831\/29334\)\), 1}, {1, 0}}, {{58\/4889, 
            6}, {\(-\(4831\/29334\)\), 1}}, {{41\/29334, 85}, {58\/4889, 
            6}}, {{10\/14667, \(-674\)}, {41\/29334, 85}}, {{1\/29334, 
            1433}, {10\/14667, \(-674\)}}}, {4831\/29334, 
        348\/4831, \(-\(41\/348\)\), \(-\(20\/41\)\), \(-\(1\/20\)\)}, {1, 6, 
        85, 674, 1433}}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[BoxData[
    \(centeredcfx[13457/97653]\)], "Input"],

Cell[BoxData[
    \({1, 7, 29, \(-283\), \(-820\), 2177, 7888, 17953, 97653}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[BoxData[
    \(\(hermite[13457/97653]\)[\([3]\)]\)], "Input"],

Cell[BoxData[
    \({1, 7, 29, 283, 537, 820, 1357, 2177, 7888, 17953}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[BoxData[
    \(ContinuedFraction[13457/97653]\)], "Input"],

Cell[BoxData[
    \({0, 7, 3, 1, 8, 1, 1, 1, 1, 1, 3, 2, 5}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[BoxData[
    \(Map[Denominator, 
      Map[FromContinuedFraction, 
        Table[Take[%, k], {k, 1, Length[%]}]]]\)], "Input"],

Cell[BoxData[
    \({1, 7, 22, 29, 254, 283, 537, 820, 1357, 2177, 7888, 17953, 
      97653}\)], "Output"]
}, Open  ]]
}, Open  ]]
},
FrontEndVersion->"5.1 for Microsoft Windows",
ScreenRectangle->{{0, 1024}, {0, 685}},
AutoGeneratedPackage->Automatic,
WindowSize->{495, 527},
WindowMargins->{{Automatic, 141}, {Automatic, 22}},
ShowSelection->True
]

(*******************************************************************
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[1776, 53, 37, 0, 95, "Title"],
Cell[1816, 55, 101, 2, 50, "Subsubtitle"],
Cell[1920, 59, 5040, 68, 1439, "Text"],

Cell[CellGroupData[{
Cell[6985, 131, 49, 1, 30, "Input"],
Cell[7037, 134, 191, 3, 76, "Output"]
}, Open  ]],
Cell[7243, 140, 90, 2, 30, "Input",
  InitializationCell->True],
Cell[7336, 144, 90, 2, 30, "Input",
  InitializationCell->True],
Cell[7429, 148, 119, 3, 52, "Text"],
Cell[7551, 153, 534, 11, 190, "Input",
  InitializationCell->True],
Cell[8088, 166, 238, 4, 71, "Text"],
Cell[8329, 172, 461, 9, 150, "Input",
  InitializationCell->True],
Cell[8793, 183, 142, 3, 52, "Text"],
Cell[8938, 188, 376, 7, 90, "Input",
  InitializationCell->True],
Cell[9317, 197, 346, 6, 110, "Input",
  InitializationCell->True],
Cell[9666, 205, 468, 8, 130, "Input",
  InitializationCell->True],
Cell[10137, 215, 1104, 19, 370, "Input",
  InitializationCell->True],
Cell[11244, 236, 375, 9, 166, "Text"],

Cell[CellGroupData[{
Cell[11644, 249, 48, 1, 30, "Input"],
Cell[11695, 252, 284, 4, 110, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[12016, 261, 52, 1, 30, "Input"],
Cell[12071, 264, 393, 6, 203, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[12501, 275, 57, 1, 30, "Input"],
Cell[12561, 278, 90, 1, 29, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[12688, 284, 66, 1, 30, "Input"],
Cell[12757, 287, 84, 1, 29, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[12878, 293, 63, 1, 30, "Input"],
Cell[12944, 296, 73, 1, 29, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[13054, 302, 131, 3, 70, "Input"],
Cell[13188, 307, 107, 2, 48, "Output"]
}, Open  ]]
}, Open  ]]
}
]
*)



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