(* Content-type: application/mathematica *)

(*** Wolfram Notebook File ***)
(* http://www.wolfram.com/nb *)

(* CreatedBy='Mathematica 6.0' *)

(*CacheID: 234*)
(* Internal cache information:
NotebookFileLineBreakTest
NotebookFileLineBreakTest
NotebookDataPosition[       145,          7]
NotebookDataLength[    173692,       4816]
NotebookOptionsPosition[    169728,       4691]
NotebookOutlinePosition[    170135,       4708]
CellTagsIndexPosition[    170092,       4705]
WindowFrame->Normal*)

(* Beginning of Notebook Content *)
Notebook[{

Cell[CellGroupData[{
Cell[TextData[{
 StyleBox["Calculation of the critical radius of a fissionable device\n",
  FontFamily->"Times New Roman"],
 StyleBox[" Comments on a note in the", "Subsubtitle",
  FontFamily->"Times New Roman"],
 StyleBox[" \"Los Alamos Primer\" ", "Subsubtitle",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["(by Robert Serber and Richard Rhodes) and an article on neutron \
multiplication by Rudolf Peierls. \n ", "Subsubtitle",
  FontFamily->"Times New Roman"],
 StyleBox[" ", "Subsubtitle",
  FontFamily->"Times New Roman",
  FontWeight->"Plain",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["A Mathematica document by S\[ODoubleDot]ren Molander\n  Change \
notes:\n  2005-09-09: Fixed a typo in Eq 1''. Minor changes in text and \
elaborated an integral.\n  2005-09-10: Fixed a minus sign in the Eigenvalue \
calculation of Peierls formula. The plot of the Lf/f Eigen function is now \
correct. Also compared Peierls approximation with a more accurate one and \
made color coded plots. \n  2005-09-11: Checked the solution of Peierls \
forumla for \[Beta]\[Rule]\[Alpha] which checks much better with the ", 
  "Subsubtitle",
  FontFamily->"Times New Roman",
  FontSize->12,
  FontWeight->"Plain",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Primer", "Subsubtitle",
  FontFamily->"Times New Roman",
  FontSize->12,
  FontWeight->"Plain",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" version. Rewrote some parts of the text\n  2005-09-14: Made \
numerical values modular for swift evaluation of the document for either U235 \
or PU239. \n  2005-09-17: Changed notation from total transport cross section \
to the more correct transport cross section. \n  2005-10-09: Fixed a factor \
of 1/2 missing in eq. 1''.  \n", "Subsubtitle",
  FontFamily->"Times New Roman",
  FontSize->12,
  FontWeight->"Plain",
  FontVariations->{"CompatibilityType"->0}]
}], "Title"],

Cell[CellGroupData[{

Cell["Introduction", "Section"],

Cell[TextData[{
 StyleBox["One can say that modern nuclear physics started with Werner \
Heisenberg who applied quantum mechanics to the atomic nucleus 1932, the same \
year that James Chadwick discovered the neutron. In 1934 Enrico Fermi \
discovered unknown radioactive elements when bombarding a substance with \
neutrons and in 1938  Strassman and Hahn performed the first controlled \
fission experiment. In 1939 Lise Meitner and Otto Frisch used Niels Bohr's \
liquid drop model to conclude that Uranium was transmuted into Barium under \
the release of large amounts of energy.  In 1939 an explosion of theoretical \
and experimental work was done on nuclear fission and the same year the \
notion of critical mass first saw the day, at least in published form. Leo \
Szilard and and Enrico Fermi were among the first to realize that nuclear \
fission could trigger a chain reaction, but they decided not to publish their \
findings to prevent Germany from using their results. It has been disputed if \
Werner Heisenberg tried delibrately to slow down the wartime German efforts \
towards a nuclear bomb, what is known is that the efforts were behind those \
of the allies. \n\nWhen I did nuclear physics in undergraduate school I was \
interested in the physics of neutron multiplication and how it gave rise to \
the notion of \"critical mass\". Naturally I never expected to find a \
realistic calculation using detailed models (which is rightly classified \
information), but was more interested in a \"back of the envelope\" argument \
that highlighted the basic physics. However, I had never found any \
information on the internet or in books on this subject, so I decided to \
compile this document. It is an attempt to fill in some of the details of the \
critical mass calculations based on the article by Rudolf Peierls published \
1939 [1]  and the notes in the book ",
  FontFamily->"Times New Roman"],
 StyleBox["\"The Los Alamos Primer\"",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["[2]. The ",
  FontFamily->"Times New Roman"],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" was a set of 5 lectures given by Robert Serber in April 1943 for \
the physicists at Los Alamos on the basic physics and the related engineering \
problems of building a nuclear bomb. Much of the material had been worked out \
at a small group of people at Berkely headed by Robert Oppenheimer and \
including names like Hans Bethe, Edward Teller and Richard Tolman. The Primer \
was originally declassified in 1965 and has now been published in an extended \
book form by Robert Serber and Richard Rhodes, where additional material has \
been included. \n\n Many of the physical properties of Uranium were not well \
known at the time when the document was written down and here I have tried to \
use the modern values of cross sections and densities (I am not sure about \
the value for the absolute geometric cross sections). The effect of an \
external shell (or \"tamper\", an outer shell which mirrors neutrons makes \
the effective critical mass even smaller)  has been disregarded, readers are \
referred to the ",
  FontFamily->"Times New Roman"],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" and the seminal book by Richard Rhodes ",
  FontFamily->"Times New Roman"],
 StyleBox["\"The Making of the Atom Bomb\"",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["[3] for more details. One of the first notes on the critical mass \
and neutron multiplication was published by Francis Perrin [4]; I have not \
read this paper so I can't comment on it here.\n \nFirst let us define some \
important physical properties of fissionable materials, the most central ones \
being the \"cross section\" (see the ",
  FontFamily->"Times New Roman"],
 StyleBox["Endnotes",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" in the ",
  FontFamily->"Times New Roman"],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["). The geometric cross section is the area taken up by a nucleus \
of radius ",
  FontFamily->"Times New Roman"],
 StyleBox["R",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". Assuming ",
  FontFamily->"Times New Roman"],
 StyleBox["R = ",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["10", 
    RowBox[{"-", "12"}]], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["cm one gets ",
  FontFamily->"Times New Roman"],
 StyleBox["\[Sigma]=\[Pi]",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["R", "2"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["\[TildeEqual] 3 ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   RowBox[{
    SuperscriptBox["10", 
     RowBox[{
      RowBox[{"-", "24"}], " "}]], " ", 
    SuperscriptBox["cm", 
     RowBox[{"-", "2"}]]}], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[".  Assume a thin foil of thickness d and area A with a number \
density of ",
  FontFamily->"Times New Roman"],
 StyleBox["n ",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["= ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[Rho]", " ", 
     FormBox[
      SubscriptBox["N", "A"],
      TraditionalForm]}], "M"], TraditionalForm]],
  FontSlant->"Italic"],
 "  (",
 Cell[BoxData[
  FormBox[
   SubscriptBox["N", "A"], TraditionalForm]],
  FontWeight->"Bold",
  FontSlant->"Italic"],
 " is Avogadro's number 6.022 ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox[
    RowBox[{"10", " "}], "23"], TraditionalForm]]],
 ") ",
 StyleBox["neutrons per unit volume and M is the mass in atomic units. For \
U235 n \[Tilde] ",
  FontFamily->"Times New Roman"],
 " 4.84 ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["10", 
    RowBox[{"22", " "}]], TraditionalForm]]],
 StyleBox["(see below). To get some idea on how sparse the atoms are \
scattered in the foil, let's calculate the average distance between atoms. \
This is simply 1/",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   RadicalBox["n", "3"], TraditionalForm]]],
 StyleBox[" \[Tilde] 2.7 ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["10", 
    RowBox[{"-", "8"}]], TraditionalForm]]],
 " or 27000 times the size of an atom and from the point of view of a small \
particle entering the foil at high speed, it sees mostly empty space. ",
 StyleBox["The total area filled with nuclei is thus",
  FontFamily->"Times New Roman"],
 StyleBox[" ",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["nAd\[Sigma] ",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" and a neutron passing through the foil has a chance of hitting a \
foil nucleus with probablity ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"area", " ", "of", " ", "nuceli"}], 
    RowBox[{"total", " ", "area"}]], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["= ",
  FontFamily->"Times New Roman"],
 StyleBox["nd\[Sigma]",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". In a similar way, the probability of fission is given as a cross \
section ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   StyleBox[
    SubscriptBox["\[Sigma]", "f"],
    FontFamily->"Courier"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[" and the fraction of collisions that give rise to fission is",
  FontFamily->"Times New Roman"],
 StyleBox[" ",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 Cell[BoxData[
  FormBox[
   StyleBox[
    FractionBox[
     SubscriptBox["\[Sigma]", "f"], "\[Sigma]"],
    FontFamily->"Courier"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". The number of fissions  in the foil is  ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   StyleBox[
    SubscriptBox["C", "f"],
    FontWeight->"Bold"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["=n",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox["\[Sigma]", TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["d",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   StyleBox[
    RowBox[{
     SubscriptBox["\[Sigma]", "f"], "/", "\[Sigma]"}],
    FontWeight->"Bold"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["=n",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "f"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["d.",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" ",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["The time it takes a neutron to pass through the foil is ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["d=v\[Tau]",
  FontFamily->"Times New Roman",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[", where v is the velocity which for 1 Mev neutrons can be \
approximated by ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["E=",
  FontFamily->"Times New Roman",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   FractionBox[
    SuperscriptBox["mv", "2"], "2"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["and is about 1.4 ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["10", "9"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[" cm/s. This means that the rate of fissions is ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["C", "f"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 "/\[Tau]",
 StyleBox[" and the inverse gives the time per fission, ",
  FontFamily->"Times New Roman"],
 StyleBox["\[Tau]=",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", 
    RowBox[{"n", " ", 
     SubscriptBox["\[Sigma]", "f"], "v"}]], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". ",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["\n",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" ",
  FontFamily->"Times New Roman"],
 StyleBox["\nThe cross sections here are given below for the relevant \
energies involved, around 1 Mev (the values here are given in the extended \
version of the ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["). Measurements made by several groups before 1939 (including \
Niels Bohr) showed that U235 had a significant fission cross section over a \
wide range of energies, which was important since only fast neutrons could \
make a significant contribution to a possible chain reaction in a bomb \
(unlike a nuclear fission reactor, where the opposite is true). Early \
experiments showed that the fission cross section of Uranium 238 drops \
sharply at energies < 1 Mev so natural Uranium has to be enriched to extract \
the isotope 235. Plutonium is a more efficient medium (higher fission cross \
section), it was discovered in 1942 by Glenn Seaborg and the isotope PU239 \
was used in the Trinity explosion 1945. \n\n",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Technical notes", "Subsubsection",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["\nThe ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Mathematica",
  FontFamily->"Times New Roman",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" document can evaluated using values for U235 or for PU 239, see \
below. The first time the notebook is evaluated there will be some minor \
warning messages mostly related to the naming of variables (they will \
disappear the 2:nd time).  The evaluation of the integral in eq. 12 takes a \
long time, please be patient. Many of the physical properties of PU239 are \
classified, so the numerical results for this case should be taken with a \
grain of salt, especially the total cross section. ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Mathematica",
  FontFamily->"Times New Roman",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" consistently outputs numerical values to 15 decimals in the \
html-document. This is an artefact and the real accuracy is in the order of \
1-2 decimal places. ",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}]
}], "Text",
 CellFrame->{{0, 0}, {0.5, 0}}],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{
   RowBox[{"MaterialIndex", " ", "=", " ", "1"}], ";"}], "  ", 
  RowBox[{"(*", " ", 
   RowBox[{
    RowBox[{"Use", " ", "1", " ", "for", " ", "U235"}], ",", " ", 
    RowBox[{"2", " ", "for", " ", "PU239"}]}], " ", "*)"}], 
  "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"Clear", "[", 
    RowBox[{
    "sf", ",", "st", ",", "ny", ",", "L", ",", "tf", ",", "Eq", ",", "En", 
     ",", "v", ",", "rho", ",", "nd", ",", "Na", ",", "m", ",", "b"}], "]"}], 
   ";"}], " "}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"rhoU235", " ", "=", " ", "18.9"}], ";"}], 
  " "}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"rhoPU239", " ", "=", " ", "19.84"}], " ", ";"}], " ", 
  RowBox[{"(*", " ", 
   RowBox[{"alpha", " ", "phase"}], " ", "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"MU235", " ", "=", " ", "235"}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"MPU239", " ", "=", " ", "239"}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"sfU235", " ", "=", " ", 
    RowBox[{"1.22", " ", 
     RowBox[{"10", "^", 
      RowBox[{"-", "24"}]}]}]}], ";"}], " ", 
  RowBox[{"(*", " ", 
   RowBox[{"1.5", " ", "in", " ", "primer"}], " ", 
   "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"sfPU239", " ", "=", " ", 
   RowBox[{"1.73", " ", 
    RowBox[{"10", "^", 
     RowBox[{"-", "24"}]}]}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"stU235", " ", "=", " ", 
   RowBox[{"4.0", " ", 
    RowBox[{"10", "^", 
     RowBox[{"-", "24"}]}]}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"stPU239", " ", "=", " ", 
   RowBox[{"8.5", " ", 
    RowBox[{"10", "^", 
     RowBox[{"-", "24"}]}]}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"nyU235", " ", "=", " ", "2.52"}], ";"}], " ", 
  RowBox[{"(*", " ", 
   RowBox[{"2.2", " ", "in", " ", "primer"}], " ", 
   "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"nyPU239", " ", "=", " ", "2.95"}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"Na", " ", "=", " ", 
   RowBox[{"6.022", " ", 
    RowBox[{"10", "^", "23"}]}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"rhoList", " ", "=", " ", 
   RowBox[{"{", 
    RowBox[{"rhoU235", ",", "rhoPU239"}], "}"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"Mlist", " ", "=", " ", 
   RowBox[{"{", 
    RowBox[{"MU235", ",", "MPU239"}], "}"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"sfList", " ", "=", " ", 
   RowBox[{"{", 
    RowBox[{"sfU235", ",", "sfPU239"}], "}"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"stList", " ", "=", " ", 
   RowBox[{"{", 
    RowBox[{"stU235", ",", "stPU239"}], "}"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"nyList", " ", "=", " ", 
    RowBox[{"{", 
     RowBox[{"nyU235", ",", "nyPU239"}], "}"}]}], ";"}], 
  "\[IndentingNewLine]", "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"m", " ", "=", " ", 
   RowBox[{"1.6759", " ", 
    RowBox[{"10", "^", 
     RowBox[{"-", "27"}]}]}]}], " ", ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"En", " ", "=", " ", 
   RowBox[{
    RowBox[{"10", "^", "6"}], " ", "1.602177", " ", 
    RowBox[{"10", "^", 
     RowBox[{"-", "19"}]}]}]}], " ", ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"v", " ", "=", " ", 
   RowBox[{
    RowBox[{"Sqrt", "[", 
     RowBox[{"2", " ", 
      RowBox[{"En", "/", "m"}]}], "]"}], "*", "100"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"rho", " ", "=", " ", 
   RowBox[{"rhoList", "[", 
    RowBox[{"[", "MaterialIndex", "]"}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"nd", " ", "=", " ", 
    RowBox[{
     RowBox[{
      RowBox[{"rhoList", "[", 
       RowBox[{"[", "MaterialIndex", "]"}], "]"}], "/", 
      RowBox[{"Mlist", "[", 
       RowBox[{"[", "MaterialIndex", "]"}], "]"}]}], " ", "Na"}]}], ";"}], 
  "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"st", " ", "=", " ", 
   RowBox[{"stList", "[", 
    RowBox[{"[", "MaterialIndex", "]"}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"sf", " ", "=", " ", 
   RowBox[{"sfList", "[", 
    RowBox[{"[", "MaterialIndex", "]"}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"ny", " ", "=", " ", 
   RowBox[{"nyList", "[", 
    RowBox[{"[", "MaterialIndex", "]"}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"tf", " ", "=", " ", 
   RowBox[{"1", "/", 
    RowBox[{"(", 
     RowBox[{"nd", " ", "sf", " ", "v"}], ")"}]}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{" ", 
  RowBox[{
   RowBox[{"L", " ", "=", " ", 
    RowBox[{"1", "/", 
     RowBox[{"(", 
      RowBox[{"nd", " ", "st"}], ")"}]}]}], ";"}], "\[IndentingNewLine]", 
  "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
 RowBox[{"TableForm", "[", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{
    "sf", ",", "st", ",", "ny", ",", "rho", ",", "v", ",", "tf", ",", "L"}], 
    "}"}], ",", "\[IndentingNewLine]", 
   RowBox[{"TableHeadings", " ", "\[Rule]", " ", 
    RowBox[{"{", 
     RowBox[{
      RowBox[{"{", 
       RowBox[{
       "\"\<Fission cross section \!\(\*SubscriptBox[\(\[Sigma]\), \(f\)]\) [\
\!\(\*SuperscriptBox[\(cm\), \(2\)]\)]\>\"", ",", " ", "\[IndentingNewLine]", 
        "\"\<Total cross section \!\(\*SubscriptBox[\(\[Sigma]\), \(\(t\)\(\\\
\ \)\)]\)[\!\(\*SuperscriptBox[\(cm\), \(2\)]\)]\>\"", ",", 
        "\[IndentingNewLine]", 
        "\"\<Average number of neutrons per fission \[Nu]\>\"", ",", 
        "\[IndentingNewLine]", 
        "\"\<Density [g/\!\(\*SuperscriptBox[\(cm\), \(3\)]\)]\>\"", ",", 
        "\[IndentingNewLine]", "\"\<Velocity of 1 Mev Neutron [cm/s]\>\"", 
        ",", "\[IndentingNewLine]", "\"\<Time between fissions [s]\>\"", ",", 
        "\[IndentingNewLine]", 
        "\"\<Mean free path of U 235 Neutron [cm]\>\""}], "}"}], ",", 
      "None"}], "}"}]}], ",", 
   RowBox[{"TableSpacing", "\[Rule]", " ", 
    RowBox[{"{", 
     RowBox[{"1", ",", "10"}], "}"}]}]}], "]"}]}], "Input",
 CellChangeTimes->{3.38866298921875*^9}],

Cell[BoxData[
 TagBox[
  TagBox[GridBox[{
     {
      TagBox["\<\"Fission cross section \\!\\(\\*SubscriptBox[\\(\[Sigma]\\), \
\\(f\\)]\\) [\\!\\(\\*SuperscriptBox[\\(cm\\), \\(2\\)]\\)]\"\>",
       HoldForm], "1.22`*^-24"},
     {
      TagBox["\<\"Total cross section \\!\\(\\*SubscriptBox[\\(\[Sigma]\\), \
\\(\\(t\\)\\(\\\\ \\)\\)]\\)[\\!\\(\\*SuperscriptBox[\\(cm\\), \\(2\\)]\\)]\"\
\>",
       HoldForm], "4.0000000000000004`*^-24"},
     {
      TagBox["\<\"Average number of neutrons per fission \[Nu]\"\>",
       HoldForm], "2.52`"},
     {
      TagBox["\<\"Density [g/\\!\\(\\*SuperscriptBox[\\(cm\\), \\(3\\)]\\)]\"\
\>",
       HoldForm], "18.9`"},
     {
      TagBox["\<\"Velocity of 1 Mev Neutron [cm/s]\"\>",
       HoldForm], "1.3827580447248204`*^9"},
     {
      TagBox["\<\"Time between fissions [s]\"\>",
       HoldForm], "1.223937616012325`*^-8"},
     {
      TagBox["\<\"Mean free path of U 235 Neutron [cm]\"\>",
       HoldForm], "5.161849233586198`"}
    },
    GridBoxAlignment->{
     "Columns" -> {{Left}}, "ColumnsIndexed" -> {}, "Rows" -> {{Baseline}}, 
      "RowsIndexed" -> {}},
    GridBoxDividers->{
     "Columns" -> {False, {True}, False}, "ColumnsIndexed" -> {}, 
      "Rows" -> {{False}}, "RowsIndexed" -> {}},
    GridBoxSpacings->{"Columns" -> {
        Offset[0.27999999999999997`], {
         Offset[0.5599999999999999]}, 
        Offset[0.27999999999999997`]}, "ColumnsIndexed" -> {}, "Rows" -> {
        Offset[0.2], {
         Offset[0.4]}, 
        Offset[0.2]}, "RowsIndexed" -> {}}],
   OutputFormsDump`HeadedColumn],
  Function[BoxForm`e$, 
   TableForm[
   BoxForm`e$, 
    TableHeadings -> {{
      "Fission cross section \!\(\*SubscriptBox[\(\[Sigma]\), \(f\)]\) \
[\!\(\*SuperscriptBox[\(cm\), \(2\)]\)]", 
       "Total cross section \!\(\*SubscriptBox[\(\[Sigma]\), \(\(t\)\(\\ \
\)\)]\)[\!\(\*SuperscriptBox[\(cm\), \(2\)]\)]", 
       "Average number of neutrons per fission \[Nu]", 
       "Density [g/\!\(\*SuperscriptBox[\(cm\), \(3\)]\)]", 
       "Velocity of 1 Mev Neutron [cm/s]", "Time between fissions [s]", 
       "Mean free path of U 235 Neutron [cm]"}, None}, 
    TableSpacing -> {1, 10}]]]], "Output",
 CellChangeTimes->{3.388662976359375*^9, 3.38866306759375*^9, 
  3.407484848984375*^9, 3.407485487765625*^9, 3.439910638625*^9}]
}, Open  ]]
}, Open  ]],

Cell[CellGroupData[{

Cell["Critical mass", "Section"],

Cell[TextData[{
 "In kinetic gas theory the mean free path is defined as",
 StyleBox[" ",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox["\[ScriptL]=",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", 
    RowBox[{"n", " ", 
     SubscriptBox["\[Sigma]", "t"]}]], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["where ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["n",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" ",
  FontWeight->"Bold",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["= number of particles per volume and ",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "t"], TraditionalForm]]],
 StyleBox[" is the transport cross section. The transport cross section is \
defined as the elastic scattering ",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "el"], TraditionalForm]]],
 StyleBox["\[Times] 1-cos(\[Theta]) where \[Theta] is the scattering angle. \
If the neutrons were scattered homogenously in all directions ",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "el"], TraditionalForm]]],
 "= ",
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "t"], TraditionalForm]]],
 StyleBox[" but in practice it is smaller and has to determined by laboratory \
measurements. ",
  FontVariations->{"CompatibilityType"->0}],
 "An important relation is the one between the diffusion constant ",
 StyleBox["D",
  FontSlant->"Italic"],
 " and the mean free path, ",
 StyleBox["D=",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[ScriptL]", " ", "v"}], "3"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[" ",
  FontSlant->"Italic"],
 StyleBox["where v is the average velocity and the constant ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["D",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" is a part in Fick's law of diffusion which states that the \
current ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["J= -D\[Del]N(r)",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[". ",
  FontWeight->"Bold",
  FontVariations->{"CompatibilityType"->0}],
 " As can be seen above, this gives an approximate mean free path of Uranium \
235 of 5.1 cm.",
 StyleBox[" ",
  FontWeight->"Bold",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Since the mean free path turns out to be in the same order as the \
critical radius to be calculated, standard diffusion theory is not really \
applicable. This was the approach of P",
  FontVariations->{"CompatibilityType"->0}],
 "eierls",
 StyleBox[" in his paper of 1939 and he derives an integral equation for the \
number of neutrons ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["N",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" ",
  FontWeight->"Bold",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["as a function of time and position. The reasoning is as follows: \
The number of scattered and secondary neutrons from an element at ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["(x',y',z')",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" in time ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["dt ",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["is ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["\[Beta] vdt' N(x',y',z',t')dx'dy'dz' ",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["where",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" v",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" = neutron velocity and \[Beta] is the number of neutrons emerging \
after each fission. When these neutrons have travelled a distance ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["r",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" they are reduced by a factor of ",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   StyleBox[
    SuperscriptBox["e", 
     RowBox[{"-", "\[Alpha]r"}]],
    FontSlant->"Italic"], TraditionalForm]]],
 " (",
 StyleBox["1/\[Alpha]",
  FontSlant->"Italic"],
 " is the mean free path) and they will fill a shell with radius ",
 StyleBox["r",
  FontSlant->"Italic"],
 " and thickness ",
 StyleBox["vdt .",
  FontSlant->"Italic"],
 "The density of neutrons after a time",
 StyleBox[" t' = t-",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["r", "c"], TraditionalForm]],
  FontSlant->"Italic"],
 "is therefore ",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Beta]", 
    RowBox[{"4", " ", "\[Pi]", " ", 
     SuperscriptBox["r", "2"]}]], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["N(x',y',z',t') ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["e", 
    RowBox[{"-", "\[Alpha]r"}]], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["dx'dy'dz'. ",
  FontSlant->"Italic"],
 "Integrating over the interior of the sphere leads to the integral equation",
 StyleBox["\n\n",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  RowBox[{
   StyleBox[
    RowBox[{"N", 
     RowBox[{"(", 
      RowBox[{"x", ",", "y", ",", "z", ",", "t"}], ")"}]}],
    FontSlant->"Plain"], 
   StyleBox["=",
    FontWeight->"Bold",
    FontSlant->"Italic"], 
   StyleBox[
    RowBox[{
     StyleBox[
      FractionBox["\[Beta]", 
       RowBox[{"4", "\[Pi]"}]],
      FontWeight->"Bold"], 
     StyleBox[
      RowBox[{"\[Integral]", 
       RowBox[{"N", 
        RowBox[{"(", 
         RowBox[{
          RowBox[{"x", "'"}], ",", 
          RowBox[{"y", "'"}], ",", 
          RowBox[{"z", "'"}], ",", 
          RowBox[{"t", "-", 
           RowBox[{"R", "/", "c"}]}]}], ")"}], 
        FractionBox[
         SuperscriptBox["e", 
          RowBox[{"-", "\[Alpha]r"}]], 
         SuperscriptBox["r", "2"]], 
        RowBox[{"\[DifferentialD]", "V"}]}]}],
      FontWeight->"Plain"]}],
    FontSlant->"Plain"]}]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 "'",
 StyleBox[" \t",
  FontSlant->"Italic"],
 "\t\t\t\t\t\t\t\t\t\t\t(1)\n\nSeparating the variables N(x,y,z,t)=N(x,y,x)",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["e", 
    RowBox[{
     RowBox[{"-", "\[Lambda]t"}], " "}]], TraditionalForm]],
  FontWeight->"Bold"],
 "and setting \[Lambda]= 0 for criticality,\n\n",
 Cell[BoxData[
  RowBox[{
   StyleBox[
    RowBox[{"N", 
     RowBox[{"(", 
      RowBox[{"x", ",", "y", ",", "z"}], ")"}]}],
    FontSlant->"Plain"], 
   StyleBox["=",
    FontWeight->"Bold",
    FontSlant->"Italic"], 
   StyleBox[
    RowBox[{
     StyleBox[
      FractionBox["\[Beta]", 
       RowBox[{"4", "\[Pi]"}]],
      FontWeight->"Bold"], 
     RowBox[{"\[Integral]", 
      RowBox[{"N", 
       RowBox[{
        StyleBox["(",
         FontWeight->"Plain"], 
        StyleBox[
         RowBox[{
          RowBox[{"x", "'"}], ",", 
          RowBox[{"y", "'"}], ",", 
          RowBox[{"z", "'"}]}],
         FontWeight->"Plain"], 
        StyleBox[")",
         FontWeight->"Bold"]}], 
       StyleBox[
        FractionBox[
         SuperscriptBox["e", 
          RowBox[{"-", "\[Alpha]r"}]], 
         SuperscriptBox["r", "2"]],
        FontWeight->"Plain"], 
       StyleBox[
        RowBox[{"\[DifferentialD]", "V"}],
        FontWeight->"Plain"]}]}]}],
    FontSlant->"Plain"]}]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 "'\t\t\t\t\t\t\t\t\t\t\t\t\t\t(1')\n\n Here ",
 StyleBox["\[Alpha]",
  FontSlant->"Italic"],
 " is defined as",
 StyleBox[" \[Alpha]=n(",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "s"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["+",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "a"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["+",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "d"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[")",
  FontSlant->"Italic"],
 StyleBox[" ",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["- scattering, diffusion and absorption cross sections - which is \
the inverse of the path length, i.e. ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["\[ScriptL]=1/\[Alpha]",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[".  \[Beta] is defined as the number of neutrons emerging from the \
fission process and",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" \[Beta]=n(",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "s"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["+\[Nu]",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "d"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[")",
  FontSlant->"Italic"],
 ". ",
 StyleBox["For spontaneous fission - criticality - to occur \[Beta] must \
exceed \[Alpha]. P",
  FontVariations->{"CompatibilityType"->0}],
 "eierls",
 StyleBox[" evaluates this integral in two extreme cases: ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["1-",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", "\[Beta]"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["\[LessLess] 1 ",
  FontSlant->"Italic"],
 "and",
 StyleBox[" ",
  FontWeight->"Bold",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" 1-",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", "\[Beta]"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["\[Tilde]1",
  FontSlant->"Italic"],
 " . The former case is the one leading to standard diffusion theory, and the \
other leads to an Eigenvalue problem, briefly discussed below. The integrand \
is assumed to be slowly varying and is expanded in a Taylor series, retaining \
the 2:nd order terms. Changing to spherical coordinates leads to 3 angular \
and 2 exponential integrals. The three angular integrals all yield the same \
value, and the x-integral is\n\n",
 Cell[BoxData[
  RowBox[{
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "\[Pi]"], 
    RowBox[{"sin", 
     SuperscriptBox[
      RowBox[{"(", "\[Theta]", ")"}], "3"], 
     RowBox[{"\[DifferentialD]", "\[Theta]"}], 
     RowBox[{
      SubsuperscriptBox["\[Integral]", "0", 
       RowBox[{"2", "\[Pi]"}]], 
      RowBox[{"cos", 
       SuperscriptBox[
        RowBox[{"(", "\[CurlyPhi]", ")"}], "2"], 
       RowBox[{"\[DifferentialD]", "\[Phi]"}]}]}]}]}], "=", 
   FractionBox[
    RowBox[{"4", "\[Pi]"}], "3"]}]],
  FontFamily->"Times New Roman"],
 "\n\nand one of the exponential integrals is\n\n",
 Cell[BoxData[
  RowBox[{
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "\[Infinity]"], 
    StyleBox[
     RowBox[{
      SuperscriptBox["r", "2"], 
      SuperscriptBox["e", 
       RowBox[{
        RowBox[{"-", "\[Alpha]"}], "r"}]], 
      RowBox[{"\[DifferentialD]", "r"}]}],
     FontFamily->"Courier"]}], 
   StyleBox["=",
    FontFamily->"Courier"], 
   StyleBox[
    FractionBox["2", 
     SuperscriptBox["\[Alpha]", "3"]],
    FontFamily->"Courier"]}]],
  FontFamily->"Times New Roman"],
 "\n\nCollecting terms gives \n\nN(x,y,z)=\[Beta] ",
 Cell[BoxData[
  FormBox[
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "\[Infinity]"], 
    RowBox[{
     RowBox[{"(", 
      RowBox[{
       RowBox[{"N", "(", 
        RowBox[{"x", ",", "y", ",", "z"}], ")"}], " ", "+", 
       RowBox[{
        FractionBox[
         SuperscriptBox["r", "2"], "3"], 
        FractionBox["1", "2"], 
        RowBox[{"(", 
         RowBox[{
          FractionBox[
           RowBox[{
            SuperscriptBox["\[PartialD]", "2"], "N"}], 
           RowBox[{"\[PartialD]", 
            SuperscriptBox["x", "2"]}]], "+", 
          FractionBox[
           RowBox[{
            SuperscriptBox["\[PartialD]", 
             RowBox[{"2", " "}]], "N"}], 
           RowBox[{"\[PartialD]", 
            SuperscriptBox["y", "2"]}]], " ", "+", " ", 
          FractionBox[
           RowBox[{
            SuperscriptBox["\[PartialD]", 
             RowBox[{"2", " "}]], "N"}], 
           RowBox[{"\[PartialD]", 
            SuperscriptBox["z", "2"]}]]}], ")"}]}]}], ")"}], 
     SuperscriptBox["e", 
      RowBox[{"-", "\[Alpha]r"}]], " ", 
     RowBox[{"\[DifferentialD]", "r"}]}]}], TraditionalForm]]],
 " = \[Beta](N(x,y,z)",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Alpha]"], TraditionalForm]]],
 " + ",
 Cell[BoxData[
  FormBox[
   RowBox[{
    FractionBox["\[Beta]", 
     RowBox[{"3", 
      SuperscriptBox["\[Alpha]", "3"]}]], 
    SuperscriptBox["\[Del]", "2"]}], TraditionalForm]]],
 "N(x,y,z))\t",
 StyleBox[" \t\t\t\t\t",
  FontWeight->"Bold"],
 "(1'')\n\nThe limit is taken to infinity since the radius of the sphere is \
assumed to be much larger than the mean free path \[ScriptL]. Putting this \
together and yields a Helmholtz equation: \n\n",
 Cell[BoxData[
  StyleBox[
   RowBox[{
    RowBox[{
     RowBox[{
      RowBox[{
       SuperscriptBox["\[Del]", "2"], "N"}], 
      RowBox[{"(", "r", 
       StyleBox[")",
        FontFamily->"Courier"]}]}], "+", 
     RowBox[{"3", 
      SuperscriptBox["\[Alpha]", "2"], 
      RowBox[{"(", 
       RowBox[{"1", "-", 
        FractionBox["\[Alpha]", "\[Beta]"]}], ")"}], "N", 
      RowBox[{"(", "r", ")"}]}]}], "=", "0"}],
   FontWeight->"Plain"]],
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 "  \t\t\t\t\t\t\t\t\t\t\t\t\t\t\t(2)\n\nThe solution is a spherical Bessel \
function of order 0 or\n\nN(r) = ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"sin", "(", "kr", ")"}], "kr"], TraditionalForm]]],
 "\t\t\t\t\t\t\t\t\t\t\t\t    \t\t\t        \t\t\t ",
 StyleBox["(3)",
  FontVariations->{"CompatibilityType"->0}],
 "\n\nwith the wavenumber ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["k", "2"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["= 3 ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["\[Alpha]", "2"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["(1-",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", "\[Beta]"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[")",
  FontSlant->"Italic"],
 StyleBox[". ",
  FontWeight->"Bold"],
 StyleBox["In the ",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["Primer",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[", Serber simply states the following equation: \n\n",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FractionBox[
   RowBox[{"dN", 
    RowBox[{"(", "r", ")"}]}], "dt"]],
  FontFamily->"Times New Roman"],
 StyleBox["+ div J = ",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  RowBox[{
   FractionBox[
    RowBox[{"(", 
     RowBox[{"\[Nu]", "-", "1"}], ")"}], "\[Tau]"], "N", 
   RowBox[{"(", "r", ")"}]}]],
  FontFamily->"Times New Roman"],
 "\t\t\t\t\t\t\t\t\t\t\t\t\t\t\t\t",
 StyleBox["(4)\n\nI",
  FontVariations->{"CompatibilityType"->0}],
 "n the extended annotated version of the ",
 StyleBox["Primer",
  FontSlant->"Italic"],
 " there is a note to this chapter that derives the equation in full. \
Assuming that diffusion theory is valid Fick's law ",
 StyleBox["J= -D\[Del]N(r)",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 " gives the Laplace operator term. The source term on the r.h.s. is \
explained as follows: In each fission, \[Nu] number of neutrons produced, and \
1 neutron is always absorbed (the one that triggers a fission in each step in \
the chain reaction) which together with the time between fissions give the \
number of neutrons produced per unit time, ",
 StyleBox["\n\n",
  FontWeight->"Bold"],
 Cell[BoxData[
  StyleBox[
   RowBox[{
    FractionBox[
     RowBox[{"dN", 
      RowBox[{"(", "r", ")"}]}], "dt"], "=", 
    RowBox[{
     RowBox[{"D", 
      RowBox[{
       SuperscriptBox["\[Del]", "2"], "N"}]}], "+", 
     RowBox[{
      FractionBox[
       RowBox[{"(", 
        RowBox[{"\[Nu]", "-", "1"}], ")"}], "\[Tau]"], "N", 
      RowBox[{"(", "r", ")"}]}]}]}],
   FontFamily->"Courier"]],
  FontFamily->"Times New Roman"],
 "\t\t\t\t\t\t\t\t\t\t\t\t\t\t\t(4') \n\nSeparation of variables ",
 StyleBox["N = N1(x,y,z)",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["e", 
    RowBox[{
     SuperscriptBox["\[Nu]", "'"], 
     RowBox[{"t", "/", "\[Tau]"}]}]], TraditionalForm]],
  FontSlant->"Italic"],
 "gives\n\t\t\t\t\t\t\t\t\t\t\t\t\t\n",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["\[Del]", "2"], TraditionalForm]]],
 "N1(r) + ",
 Cell[BoxData[
  FractionBox[
   RowBox[{
    RowBox[{"-", 
     RowBox[{"\[Nu]", "'"}]}], " ", "+", "\[Nu]", "-", "1"}], "D\[Tau]"]],
  FontFamily->"Times New Roman"],
 "N1(r) = 0",
 StyleBox["\t",
  FontWeight->"Bold"],
 "\t\t\t\t\t\t\t\t\t\t\t\t\t\t\t(5)\n\nwith the same solution as (3). The \
boundary condition in the original version of the ",
 StyleBox["Primer",
  FontSlant->"Italic"],
 " is the simplest one possible is used: N1(R) ",
 StyleBox["= 0",
  FontWeight->"Bold"],
 ". This is most certainly not true but gives a first rough estimate of the \
critical radius. This boundary condition is satisfied when  the wavenumber ",
 StyleBox["k=",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Pi]", "R"], TraditionalForm]],
  FontWeight->"Bold",
  FontSlant->"Italic"],
 "or when the \"effective\" neutron number ",
 StyleBox["\[Nu]'=\[Nu]-1 - ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["\[Pi]", "2"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["D\[Tau]/",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["R", "2"], TraditionalForm]],
  FontSlant->"Italic"],
 ". ",
 StyleBox["If \[Nu]' is negative, too many neutrons escape for a chain \
reaction to be started and so the critical radius is\n\n",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["R", "c"], TraditionalForm]]],
 "= \[Pi]",
 Cell[BoxData[
  FormBox[
   SqrtBox[
    FractionBox["D\[Tau]", 
     RowBox[{"\[Nu]", "-", "1"}]]], TraditionalForm]]],
 "\nHere is where kinetic gas theory enters. The time between fission was \
shown above to be ",
 StyleBox[" \[Tau]=",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", 
    RowBox[{"n", " ", 
     SubscriptBox["\[Sigma]", "f"], "v"}]], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["=",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["\[ScriptL]", "v"], TraditionalForm]],
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    SubscriptBox["\[Sigma]", "t"], 
    SubscriptBox["\[Sigma]", "f"]], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[" ",
  FontWeight->"Bold"],
 "where ",
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "t"], TraditionalForm]]],
 " is the transport cross section, and since \n\nD=",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[ScriptL]", " ", "v"}], "3"], TraditionalForm]]],
 "\n\t\n",
 Cell[BoxData[
  FormBox[
   SubscriptBox["R", "c"], TraditionalForm]]],
 " = ",
 Cell[BoxData[
  RowBox[{
   RowBox[{"\[Pi]", 
    SqrtBox[
     FractionBox[
      RowBox[{"D", "\[Tau]"}], 
      RowBox[{"\[Nu]", "-", "1"}]]]}], "=", 
   StyleBox[
    RowBox[{"\[Pi]\[ScriptL]", 
     SqrtBox[
      FractionBox[
       SubscriptBox["\[Sigma]", "t"], 
       SubscriptBox["", 
        RowBox[{" ", 
         RowBox[{"3", 
          SubscriptBox["\[Sigma]", "f"], 
          RowBox[{"(", 
           RowBox[{"\[Nu]", "-", "1"}], ")"}]}]}]]]]}],
    FontFamily->"Times New Roman"]}]],
  FontFamily->"Times New Roman"],
 "=",
 Cell[BoxData[
  FormBox[
   StyleBox[
    FractionBox[
     RowBox[{" ", "\[Pi]", " "}], 
     StyleBox["k",
      FontWeight->"Plain",
      FontSlant->"Plain",
      FontTracking->"Plain",
      FontVariations->{"CompatibilityType"->0,
      "Masked"->False,
      "Outline"->False,
      "RotationAngle"->0,
      "Shadow"->False,
      "StrikeThrough"->False,
      "Underline"->False}]],
    FontWeight->"Plain",
    FontVariations->{"CompatibilityType"->0}], TraditionalForm]],
  FontWeight->"Bold"],
 ", where k = ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[ScriptL]"], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   SqrtBox[
    FractionBox[
     RowBox[{"3", 
      RowBox[{"(", 
       RowBox[{"\[Nu]", "-", "1"}], ")"}], 
      SubscriptBox["\[Sigma]", "f"]}], 
     SubscriptBox["\[Sigma]", "t"]]], TraditionalForm]]],
 "\t\t\t\t\t\t\t\t\t\t(6)\n\nPlugging in numerical values for U 235 gives a \
first estimate of the critical radius (~ 13 cm) and mass. Since the mass \
\[Proportional] ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox[
    SubscriptBox["R", "c"], "3"], TraditionalForm]]],
 "and ",
 Cell[BoxData[
  FormBox[
   SubscriptBox["R", "c"], TraditionalForm]]],
 " \[Proportional] \[ScriptL] \[Proportional]",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Rho]"], TraditionalForm]]],
 ", the critical mass is proportional to ~ ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", 
    SuperscriptBox["\[Rho]", "2"]], TraditionalForm]]],
 ". If a mechanism could be designed to push the active core to a higher \
density, criticality could be reached. This was one of the turning points in \
the Manhattan project. Problems with the high spontaneous fission rates of \
PU239 made the \"shooting\" mechanism (used for the Uranium device in \
Hiroshima) impossible and instead led to the design an implosion mechanism \
used to \"squeeze\" a core of Pu 239 to achieve critical mass. \n"
}], "Text",
 FontFamily->"Times New Roman"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{"Clear", "[", 
   RowBox[{"k", ",", "R2", ",", "M2"}], "]"}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"k", " ", "=", " ", 
    RowBox[{
     RowBox[{"1", "/", "L"}], " ", 
     RowBox[{"Sqrt", "[", 
      RowBox[{"3", " ", 
       RowBox[{"(", 
        RowBox[{"ny", "-", "1"}], ")"}], 
       RowBox[{"sf", "/", "st"}]}], "]"}]}]}], ";"}], " ", 
  RowBox[{"(*", " ", 
   RowBox[{
   "Wave", " ", "number", " ", "from", " ", "Diffusion", " ", "Equation"}], 
   " ", "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"R2", " ", "=", " ", 
    RowBox[{"Pi", "/", "k"}]}], ";"}], " ", 
  RowBox[{"(*", " ", 
   RowBox[{
    RowBox[{"Critical", " ", "radius", " ", "with", " ", 
     RowBox[{"B", ".", "C", ".", " ", "N"}]}], " ", "=", " ", 
    RowBox[{
     RowBox[{"0", " ", "at", " ", "r"}], " ", "=", " ", "R"}]}], " ", 
   "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"M2", " ", "=", " ", 
    RowBox[{"rho", " ", "*", "1000", " ", "*", "4", " ", "Pi", " ", 
     RowBox[{
      RowBox[{
       RowBox[{"(", 
        RowBox[{"R2", "/", "100"}], ")"}], "^", "3"}], "/", "3"}]}]}], ";"}], 
  " ", 
  RowBox[{"(*", " ", 
   RowBox[{"Mass", " ", "in", " ", "kilos"}], " ", 
   "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{"TableForm", "[", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"R2", ",", "M2"}], "}"}], ",", 
   RowBox[{"StyleForm", "->", "\"\<Courier\>\""}], ",", 
   RowBox[{"TableHeadings", "\[Rule]", " ", 
    RowBox[{"{", 
     RowBox[{
      RowBox[{"{", 
       RowBox[{
       "\"\<Critical radius [cm]\>\"", ",", "\"\<Critical Mass [kg]\>\""}], 
       "}"}], ",", "Automatic"}], "}"}]}], ",", 
   RowBox[{"TableSpacing", "\[Rule]", " ", 
    RowBox[{"{", 
     RowBox[{"1", ",", "10"}], "}"}]}]}], "]"}]}], "Input"],

Cell[BoxData[
 TagBox[
  TagBox[GridBox[{
     {
      TagBox["\<\"Critical radius [cm]\"\>",
       HoldForm], "13.75063794679407`"},
     {
      TagBox["\<\"Critical Mass [kg]\"\>",
       HoldForm], "205.8348727770857`"}
    },
    GridBoxAlignment->{
     "Columns" -> {{Left}}, "ColumnsIndexed" -> {}, "Rows" -> {{Baseline}}, 
      "RowsIndexed" -> {}},
    GridBoxDividers->{
     "Columns" -> {False, {True}, False}, "ColumnsIndexed" -> {}, 
      "Rows" -> {{False}}, "RowsIndexed" -> {}},
    GridBoxSpacings->{"Columns" -> {
        Offset[0.27999999999999997`], {
         Offset[0.5599999999999999]}, 
        Offset[0.27999999999999997`]}, "ColumnsIndexed" -> {}, "Rows" -> {
        Offset[0.2], {
         Offset[0.4]}, 
        Offset[0.2]}, "RowsIndexed" -> {}}],
   OutputFormsDump`HeadedColumn],
  Function[BoxForm`e$, 
   TableForm[
   BoxForm`e$, StyleForm -> "Courier", 
    TableHeadings -> {{"Critical radius [cm]", "Critical Mass [kg]"}, 
      Automatic}, TableSpacing -> {1, 10}]]]], "Output",
 CellChangeTimes->{3.38866297646875*^9, 3.38866306775*^9, 3.4074848490625*^9, 
  3.407485487796875*^9, 3.439910639078125*^9}]
}, Open  ]]
}, Open  ]],

Cell[CellGroupData[{

Cell["Critical mass with improved boundary condition", "Section"],

Cell[TextData[{
 StyleBox["In order to yield a better estimate a more accurate boundary \
condition has to be introduced. This was done in P",
  FontFamily->"Times New Roman"],
 "eierls",
 StyleBox["' paper, and he cited a result from astrophysics: \n\n0.71",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"dN", "(", "r", ")"}], "dr"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["+ \[Alpha]N(r)=0 \t\t\t\t\t\t\t\t\t\t\t\t\t\t\t\t(7)",
  FontFamily->"Times New Roman"],
 StyleBox["\n",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["\nThe number 0.71 is an approximation valid for a critical radius \
>> mean free path. The number is an approximation from radiative transfer \
theory, where a scattering integral on the right hand side over the sphere is \
approximated with the average value. On the left hand side, the cosine \
obliquity factor ",
  FontFamily->"Times New Roman"],
 "\[Mu]",
 StyleBox[" ",
  FontFamily->"Times New Roman"],
 "\[Congruent]",
 StyleBox["cos(",
  FontFamily->"Times New Roman"],
 "\[Theta]",
 StyleBox[") is replaced with an \"average\" cosine factor",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   OverscriptBox["\[Mu]", "~"], TraditionalForm]]],
 StyleBox[" . On such choice is a root mean square average\n\n",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{"\[LeftAngleBracket]", 
     RowBox[{"cos", " ", 
      RowBox[{"(", "\[Theta]", ")"}]}], "\[RightAngleBracket]"}], "=", 
    OverscriptBox["\[Mu]", "~"]}], TraditionalForm]],
  FormatType->"TraditionalForm"],
 StyleBox[" = ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   SqrtBox[
    SuperscriptBox[
     RowBox[{"\[LeftAngleBracket]", "\[Mu]", "\[RightAngleBracket]"}], "2"]], 
   TraditionalForm]],
  FormatType->"TraditionalForm"],
 "= ",
 Cell[BoxData[
  FormBox[
   SqrtBox[
    FractionBox[
     RowBox[{
      SubsuperscriptBox["\[Integral]", "0", "1"], 
      RowBox[{
       SuperscriptBox["\[Mu]", "2"], 
       RowBox[{"N", "(", 
        RowBox[{"\[Tau]", ",", "\[Mu]"}], ")"}], 
       RowBox[{"\[DifferentialD]", "\[Mu]"}]}]}], 
     RowBox[{
      SubsuperscriptBox["\[Integral]", "0", "1"], 
      RowBox[{
       RowBox[{"N", "(", 
        RowBox[{"\[Tau]", ",", "\[Mu]"}], ")"}], 
       RowBox[{"\[DifferentialD]", "\[Mu]"}]}]}]]], TraditionalForm]],
  FormatType->"TraditionalForm"],
 ". Assuming that N depends linearly on \[Mu] N \[TildeTilde]C\[Mu], the \
result is ",
 Cell[BoxData[
  FormBox[
   OverscriptBox["\[Mu]", "~"], TraditionalForm]],
  FormatType->"TraditionalForm"],
 "=1/",
 Cell[BoxData[
  FormBox[
   SqrtBox["2"], TraditionalForm]],
  FormatType->"TraditionalForm"],
 "\[TildeTilde]0.71",
 StyleBox["\n\n In the Endnote of the ",
  FontFamily->"Times New Roman"],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[", Serber gives a more intuitive explanation. Imagining a thin \
slice in the material of width 2\[ScriptL] (2 times the mean free path), the \
flow of neutrons going out to the left is ",
  FontFamily->"Times New Roman"],
 StyleBox["1/2 N1 v",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" and to the right",
  FontFamily->"Times New Roman"],
 StyleBox[" 1/2 N2 v",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". The net current is",
  FontFamily->"Times New Roman"],
 StyleBox[" ",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["J = 1/2(N1-N2)v",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". The difference",
  FontFamily->"Times New Roman"],
 StyleBox[" N1-N2 = ",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["dN", "dx"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox["2\[ScriptL]",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[",",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox[" and this gives",
  FontFamily->"Times New Roman"],
 StyleBox[" J = -\[ScriptL]v",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["dN", "dx"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[". Generalizing this to 3 dimensions readily gives ",
  FontFamily->"Times New Roman"],
 StyleBox["J = -",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[ScriptL]", " ", "v"}], "3"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["dN", "dr"], TraditionalForm]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" . At the boundary the neutrons have \"nowhere\" to go except out, \
so this gives them only 1 degree of freedom. This gives the following system \
of equations:\n\t\n",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", "2"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["N v = -",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[ScriptL]", " ", "v"}], "3"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox["dN", "dr"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["\t",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["\t\t(B.C.)\t\t\t\t\t\t\t\t\t\t\t\t\t(8)\n\nN(r) = ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"Sin", "(", "kr", ")"}], "kr"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["\t\t\t\t(Eq 3)\n\nPutting the solution of the Helmholtz equation \
(3) in the boundary condition (8) gives\n\n",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"sin", "(", "kr", ")"}], "kr"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["= -",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"2", "\[ScriptL]"}], 
    RowBox[{"3", "k"}]], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["(",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{
     RowBox[{"kr", " ", 
      RowBox[{"cos", "(", "kr", ")"}]}], " ", "-", " ", 
     RowBox[{"sin", "(", "kr", ")"}]}], 
    SuperscriptBox["r", "2"]], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[")  \n\nx cot(x) = 1-",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox["x", "b"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[" ",
  FontFamily->"Times New Roman"],
 StyleBox["\t\t\t\t\t",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["\t\t\t\t\t\t\t\t\t\t\t\t(9)",
  FontFamily->"Times New Roman"],
 StyleBox["\n\n",
  FontFamily->"Times New Roman",
  FontWeight->"Bold"],
 StyleBox["with",
  FontFamily->"Times New Roman"],
 StyleBox[" b=",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FractionBox[
   RowBox[{"2", "k", " ", "\[ScriptL]"}], "3"]],
  FontFamily->"Times New Roman"],
 StyleBox["and ",
  FontFamily->"Times New Roman"],
 StyleBox[" x=kr.",
  FontFamily->"Times New Roman",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox["\t\t\t\t\t\t\t\t\t\t\t\t\t\n \nThis is identical to equation (10) \
given in P",
  FontFamily->"Times New Roman"],
 "eierls",
 StyleBox[" except that ",
  FontFamily->"Times New Roman"],
 StyleBox["\[Alpha]=k\[ScriptL]",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" and the factor 2/3 is replaced with 0.71. Solving (7) numerically \
for U235  (see below) gives",
  FontFamily->"Times New Roman"],
 StyleBox[" x = 2.26",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" and\n\n",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["R", "c"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox[" = ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FractionBox["x", "k"], TraditionalForm]],
  FontFamily->"Times New Roman"],
 StyleBox["= 2.2698 ",
  FontFamily->"Times New Roman"],
 "\[ScriptL]",
 StyleBox[" ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  SqrtBox[
   FractionBox[
    SubscriptBox["\[Sigma]", "t"], 
    SubscriptBox["", 
     RowBox[{" ", 
      RowBox[{"3", 
       SubscriptBox["\[Sigma]", "f"], 
       RowBox[{"(", 
        RowBox[{"\[Nu]", "-", "1"}], ")"}]}]}]]]]],
  FontFamily->"Times New Roman"],
 StyleBox["\t\t\t\t\t\t\t\t\t\t\t\t\t\t(10)\n\nComparing this to Eq. 6 gives ",
  FontFamily->"Times New Roman"],
 StyleBox["2.26/\[Pi] ~ 0.72",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" smaller radius than the original estimate. Note that the mass is \
~ 80 kg instead of 210 kg and the effect of tamper around the device will \
further reduce the critical radius. Serber states that equations above are \
still empirical and that the full problem with an exact boundary condition \
was solved exactly by two people - Eldred Nelson and Stan Frankel - in the \
Berekely theory group. The ",
  FontFamily->"Times New Roman"],
 StyleBox["Primer",
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 StyleBox[" gives no clue to this, but most probably the solution is based on \
an integral equation like in Peierls' paper. He also states that this \
solution only differed slightly from the one from diffusion theory with the \
improved boundary condition. ",
  FontFamily->"Times New Roman"]
}], "Text",
 CellChangeTimes->{{3.452059757549412*^9, 3.4520598089826317`*^9}, {
  3.4520598460889196`*^9, 3.4520605485565987`*^9}}],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{"Clear", "[", 
   RowBox[{"b", ",", "sol", ",", "M", ",", "y", ",", "R"}], "]"}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{"b", "=", 
  FractionBox[
   RowBox[{"2", " ", "k", " ", "L"}], "3"]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"sol", "=", 
   RowBox[{"FindRoot", "[", 
    RowBox[{
     RowBox[{
      RowBox[{"x", " ", 
       RowBox[{"Cot", "[", "x", "]"}]}], "\[Equal]", 
      RowBox[{"1", "-", 
       FractionBox["x", "b"]}]}], ",", 
     RowBox[{"{", 
      RowBox[{"x", ",", "2"}], "}"}]}], "]"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"y", "=", 
   RowBox[{"x", "/.", "\[InvisibleSpace]", "sol"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"R", "=", 
   FractionBox["y", "k"]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"M", "=", 
   RowBox[{
    RowBox[{
     FractionBox["1", "3"], " ", "rho", " ", "1000", " ", "4", " ", "\[Pi]", 
     " ", 
     SuperscriptBox[
      RowBox[{"(", 
       FractionBox["R", "100"], ")"}], "3"]}], "//", "N"}]}], 
  "\[IndentingNewLine]", 
  RowBox[{"(*", 
   RowBox[{"\[IndentingNewLine]", 
    TagBox[
     FormBox[
      TagBox[
       TagBox[GridBox[{
          {
           TagBox["\"\<Root of Eq. 9\>\"",
            HoldForm], "y"},
          {
           TagBox["\"\<Improved critical Radius [cm]\>\"",
            HoldForm], "R"},
          {
           TagBox["\"\<Critical Mass [kg]\>\"",
            HoldForm], "M"}
         },
         GridBoxAlignment->{
          "Columns" -> {{Left}}, "ColumnsIndexed" -> {}, 
           "Rows" -> {{Baseline}}, "RowsIndexed" -> {}},
         GridBoxDividers->{
          "Columns" -> {False, {True}, False}, "ColumnsIndexed" -> {}, 
           "Rows" -> {{False}}, "RowsIndexed" -> {}},
         GridBoxSpacings->{"Columns" -> {
             Offset[0.27999999999999997`], {
              Offset[0.5599999999999999]}, 
             Offset[0.27999999999999997`]}, "ColumnsIndexed" -> {}, "Rows" -> {
             Offset[0.2], {
              Offset[0.4]}, 
             Offset[0.2]}, "RowsIndexed" -> {}}],
        OutputFormsDump`HeadedColumn],
       Function[BoxForm`e$, 
        TableForm[
        BoxForm`e$, 
         TableHeadings -> {{
           "Root of Eq. 9", "Improved critical Radius [cm]", 
            "Critical Mass [kg]"}, None}, TableSpacing -> {1, 10}]]],
      StandardForm],
     StandardForm,
     Editable->True], ")"}], " ", "\[IndentingNewLine]", 
   "\[IndentingNewLine]", "*)"}]}], "\[IndentingNewLine]", 
 RowBox[{"Plot", "[", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{
     RowBox[{"x", " ", 
      RowBox[{"Cot", "[", "x", "]"}]}], ",", 
     RowBox[{"1", "-", 
      FractionBox["x", "b"]}]}], "}"}], ",", 
   RowBox[{"{", 
    RowBox[{"x", ",", 
     RowBox[{"-", "1"}], ",", "3"}], "}"}], ",", 
   RowBox[{"AspectRatio", "\[Rule]", "Automatic"}], ",", 
   RowBox[{"PlotLabel", "\[Rule]", "\"\<Graphical solution to Eq. 9\>\""}]}], 
  "]"}]}], "Input",
 CellChangeTimes->{
  3.38866301421875*^9, {3.407485030515625*^9, 3.407485152625*^9}, {
   3.4074854631875*^9, 3.407485479546875*^9}}],

Cell[BoxData["0.7862145593496301`"], "Output",
 CellChangeTimes->{{3.407485155671875*^9, 3.407485164984375*^9}, {
   3.4074854665*^9, 3.407485487828125*^9}, 3.439910639171875*^9}],

Cell[BoxData["77.0016920860913`"], "Output",
 CellChangeTimes->{{3.407485155671875*^9, 3.407485164984375*^9}, {
   3.4074854665*^9, 3.407485487828125*^9}, 3.43991063934375*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {Hue[0.67, 0.6, 0.6], LineBox[CompressedData["
1:eJwV13k0Ve0XB3Bd94bMMkWTKJJSSi+hfUIqSUqGREXJ9L5ESDI2ICUZM5Uy
lVlIhniOIZQxc2QeEu49JZUy/Z7fX2d91jrrrHPWs/d37yNl5XTamsbGxqa0
io3t/1dzH4WOlRUKRWvzzzIejMHA7BGBHz8ppHH24UBlxBhctLmkNzFNoX3l
ovGe8WNgbRBX09hNIZU/JYEjGWPguIXzdWwuhYLmbI88eD8G/rWjUXsvUGh8
NGzT0dXjkMYXb2JTwULVIzQBG+9x+PaUq6/Zi4lCn7FvcbKagIx0sY7b1VPI
+QZE2stPgobfoJyjyiSyVuPOuZbyFRy7g7g3lo4jo8CqxVczU/BO/NJQAWMU
DbeX9z6QnoGfW38L+PoOIduFDs6/65nQpnqrM8R1CE32vpke2MyEvBO8cQl2
Q2jafRVnpQwT7N2kpUsNh9B72kWWlwIThmpOKv+UHUIfU856t6sxIarOoKDr
2SAS9ExeWDjLBEKrQDtXcwA9ePtAWiOSCY8PeNhcDOhDByyLLkWyscDYcWn9
2NU+pNx0IaqMnQWiSbfabM71oR0XVXYPrWZBNFeIhtPuPnSmKTNJipcFkT3P
1/r09qKy2aVbd9exIMy9AcXv6kUclh7XPu5hwf38jeLdXT1Inv+FqNdFFuhO
JDeZVfagE7THywaXWLBGYvvtgcweBO+6TKWvsCDYfy9r3K8HBRaYWJMOLAgy
OFozJ9+D9DY2j/e5s+Auy/mqkG83itLehmIesMBX/l39CdkupJx2J25HEQue
Ksht11rbhfZEnt9wp5gFb3fdv6ey0ok8nFwce0tZMK90Sle6uxNZnuD28kMs
cFHrb5gP6ESUs75ufj0Lruj9bE6e6EBznIOfs3pZoOe4tfNPajvy9Du8cWWR
BfZXg5SpsHa03HyhWnCFBfdcpqPGvNvR60Jn9S2rKKh1f2XcYtSOfAxkM9UZ
FICvRk8Kox2FPPfyNeGlQOmRUd9J6za0R1pXbusGCtbl3x1OlfmIdPRCI6XU
KAhS83Se4vuI1KSVPotoUDBf47hK8U8ryl3+Gs8BFPR0m0gVN7Wi3pwzd0Y1
KYhe2m75wb0VPbu0Y9lLl4Kt4QmV/WEtqJ63uHz/WQrEyZ2kU28jcnUSf1rr
SsHvNtGh+ZxGxA65B+3dKegaX1nxv92IWs5XmnN7UBDF3QaRCo2oefiHoM5N
CoRM3FGxTwO6Effh8hN/CnhYFRWrpD+gbk2VpvYQClYkDcojHOrQ2NjelIsp
FCQ0zCmxadUhtiMZFU9SKVC5GfvyX4k69LV9VXZPGgUuvSMRh9/XIoXMUsbR
dArGY9zsf2+rReVDwhk8ORQ0iMSLnRupQS9Kh9p3vqHA9h2E1JXUIMubDJWj
xRTQ3cZo+8JqEONmqPPFEgoOduxk8RA1qGVup3VwGQWvwsmaiifVaM2q1adq
EAWP+Sect5ytQtobRPVb6yjYh4K/PNxdhfrd2oTq6ylodVS0WOCoQoeNdHPK
31OwptnjaEdRJTJd9dYuqYEC7wfcmwKEK5HYo/q5Uy34fIWsYtoTEeKZX3I4
1kUBm0vm7Sx9hNaftUzd2E2B6sc5x7tLFWh5U3rRd+ys0MDD/5yrQCJT2y3D
PlEQzpM9GydajraNVHUXfabgwup5PauQUmRlffps3Ai+P2rpa7tGKVKseBN7
bJSCvzLsgYdZJWjtdbf6X9iRWnyVcidL0HyU4nfdcQrq/WSUKYFiZGXC49z7
hQIRAfm2i5VvUP7z5K03JimwSlR0anN+g3z2/q0X+UrBUvmB9NdtRejzI5HP
R6dwPS4YbPCKfI3OjusrR89Q4BtsXMrUfo2WD+6wlGFS0LjO3OTCz0KkeV3a
Kw/bRsUmTNO4EHWR4b5VLHzebt4MLvECxGEdXt70jQLGt3RmVHweeuT/dqP/
HAWnfXODOfTy0ONuoUMsbNhgf3J2PhflWW0yPPuTAlczv+ri0hx0W9BLZ/sv
Cvo7sjIPH8xCidlSDnm/KfD48WZcLjETcbzg5OGdx/UoVL2Jhy0TeSSOJthg
65z8FNFWmY6SbcMeiv6hIO89w+vC4ReowyYg/9xfCnQnBYo009JQ1u369hfY
46vXf9vKkYae7FQa/44toa10eaY+BaVPhff4LlBQaKXxtGV7CuL/HVlaja3v
f7QnPzgZyd3peMBYpOB2xXm9G3pJqPjUZtpt7I39dgHm2c9Rhfb+1Ars4gVX
Evieo8QNT1TmsXudRe14niWiWRmn3dZLFPDnKvS3WT5F/6KRyBhs7RnNU7HS
T9A7seKv77GzrzipbnsRh9TXPb68bZmCOfGnaYMnYlHtpxMbjLDVGxrXxs49
Rsc55mz9se94L/idjo9Gl9JDjmRiNyrKs7g1o9DG0I0lbdjCI6bn3k1GIP6T
GeW/sc0jA+t9QsMRcWSfkcQKzpPv7w+PloQi4YMnNUyxbY1lNxtzhCD9gfGQ
a9jzpXf+1p8JRn/+3HV+gB20aaRDLSkQpXbsmUzCFr8DuTnUHXT5xvT3IuyX
kwn3pDRuoebSwof12E5xc3VyIr5oRSm0ohv74aVcL768GyhCy8t/DDtHwX7P
nK4rKhDx7mZhx08fayxo+A9xxIWTv7B7vZGy76pLqDqoXH0Re+9m+4NRjgYo
wXL5OJ7zoMLMKD35eB+0Ohix/u+WaCSpUWoG0kVV0svYcqr2MoaytjAoqjPz
B9vV0LO/2sMZOncMHvmBbRSM5Ab63cHCLGT/FPb+SrrrvKYXnG45WTyAzb2Y
coGrww8kfbc2tWLffHzkbnTAbUjcyutGYk8rTWVIqwbAcR1GXja2efOD1rzp
IBBU4PGJwW60U/yl8fQ+2H/dNOiHrc5ok2wweAgGpgfbr2AbG5RuzvJ/BGoS
1ma62HxP7xfM9oaBYkyY6w7s2mlzHdV9EZBOq5Rag+2juuuTb0gkZJ/6ZjmB
z3N/4IpD7UQUFGqvVyaxWR2tyzzEY6iyOhT/GDttS1KYYWwMHP97LvJfbLEK
7TdDx+Nh92tL4Mdu5RbVlU1NgAkz7e+fcf0Fnf3y+b/lJ9Dkw73/JbbXX21/
b9NnwNmz1+Mf7JmoL2lruJNBJcx0QAXXd726wbfsN8lw0K1x2w/cLymjxaqn
LqdAgs78qkxs8z3BjY/LU4G+3+LqWuzGRoVZmasvwe5utEId7teX16LU6ten
w/GnjDWO2Hcklu84vE+HUOb3G0LY6rYtYvlbMuFq4mbBM7jfs9hdNKAzG8wX
BiRLcV6E/zOjaOSeAzX7zUbPYHs4XNniIJYLH157hDBxviTZZd9uUMiDnhCf
e+LY/318bvvb8RUkp5ZvMvpBwcXpActb6BWUyF9eGJjFecVYf45XIB9O5NuU
XMFWUY0+If0qH5rdLnx2+o7nz/NgpZM/CsBs3eJbWwrP09K6Hb1ahSCzwz1j
GOfhVAd9q3VkIcT8s3jDFLuV01fspvJrMAliQ4DzM+HqtcW060Wg9CBncA7n
717CvHZ5oRiCtfY7/8F5vs0sFt3TKwGpkV83DLDFXbuKhZ+UQOcK2+XUMQqW
Xxhkyh8shXqFsr6jeB584D/8yNivDFL/gKn7EM77wZ1mOfQKcBWs8nDopeAo
WvmVE1AJbmL/TZY3UpCocVPYOL8S8izOr7zF8+1X2a89S/2VoGt6MLrkAz6/
Yua/x5WrwGPJVisTz0e2/L6RL6NVIIlexVytoaAk5U3z5kM1cIEjVPw2nrfy
wU5p4Qu1kG4/mFX7ggJ/zulqVdk6yK/9GeKG531PwJXhodN1cJ7ynJbC+0DA
bYv1ihl1UPLolOS1JApGvI6HN5jWQ8LUc93fCRTEOcn60Ivew6XusA8BYbj/
jAaN3K82wsU/86aVeD8xLfOSOBXfCA9d6Yqy13F9S0kM7qhtBEuP25HBbhQQ
rDO2w5JNoKR3mFfXBfd74AdPvbomeBmf+SrTgYKBkteJWza2gIXb5+Q/FhSc
5a91/TLeChYudgKIoMCvSrj7zEIr/Iq5Edh1ED/f7ZJqlcBHQF6dQdPqOI/7
2JYT1D7Ct4ay7XyqOM9eqgWeDvsIPVe1B5X3UFB16FVsuXobfFcJ+cYnRYGC
a0JFREQ79O5NOzuxzII5U1730BftoPM2RqgT75NvNXx23i9rh+ach0mVf1lw
nONiwq3RdnAx+6QT+Qvvm3HSns57O8BGWjFhE4sF6WSG8sn2DpCYsTds/syC
bbwlWWvWdsHvzMqiYLzvbk7rjPMP7wGdGctaVzu8r7O/3rLrZQ9MMG9sSMT7
Na9lZEZveQ8oa+T21eL9+6+kYenerz1wT870Ad95FnSGtX4ah08w7jJU5HMa
7+PeDeK6zE+wtrLh+/AB/P5nKqOFjvbBMU2r+VdcLKhnzw5PXuoHfdlmaZtE
/P8RnuM0IDgAKpxeETIJTJiXytNbt20ANjpk5QzGMEGOKFgdqj8AJZnXJvTD
mRDoVeLpmTgAxwdeGAreZYLOz3dWpw4NgrCH5thWeyZYXXOK9dUcgpOTH5vv
KzFh8HKQ3c6AYZB1u3hX9e0MFFrcE0vWGgPmbKpXQNw0THe+/MdwbgLqa8Xu
qBlMwbFbfItN7l/hrnPlT8Xvk+CvrdVspTgDHqOPDbW0vsDQgQrtqhIWVN93
lW70G4dTPZmiOfHfwFHicq989Chs2jgkGPTzO3z4lH6F58gwsPn8c+GIyg9I
5Xrnea+uH5pt3QvWRs9Bw990KOLpAf2IGsKb+RM8HSLn6C0fYbhbeXPjzt+g
kW6r9KmyFixfPFwtd28eSmO0Rvv0s2Hu0Rqwa/oDna0m/IPHXiHV0qi8I9wL
wNvXFPrz93tk78i5fe/lRXAunjCxde1GuepPf7CnLcHpqcjNXXuG0XCU5Lth
2grUdrwxER38giLo511Py7AR5VMfjWGWhfgGXifRatiIQMvx8tUOPxCbt/8Z
VedVBHp1KK7fbB4tvDyfOyJEI0gWzx228iVEznQ7JorQiMAHj1qOdC2hgN2n
dpmL0whhWxeZ59QSEirWyu7cQCPiNh1SC96yjOTr5DLr5WjE8MA9M5V7y8h0
fDYt5yCNCHaHIlGzFVQoFfjU055GiEn11R+4y0b+F5cbIlRFIxK3njO+ZUQj
8xfjAqtqaMQcXRV5XqaR8+cD/F3qaMT+G+c7fa7RyFvSFm5tjTTC6ZugVkYY
jYzOWmMR3k0jcsPc3oY200iywlpBiEkj9hrmqifpspNCo5INguLsxMzE+CND
PTpZvCOIU9CRndjc9yaw5MZqck+HWFP4VXZi9uzQvtHg1WSG14sw4WvsRIuX
36xEwmrySVOthLgHO8FnIDxZWLGavOXE2LnpFjuRZClZFUnnIE8U3Dq1M5qd
2HgmrUE4goMcOeAdd6yCnaiG+INJZZzkuVGeCw0kO1Gm6dRV0MxJdtxPkD5R
zU4YKGhzdQ1zknWfy7IM6tmJWLHL8XqcXGSWz98KkzZ2ItdPcIFuzEVer7o+
aj3BTuyPdQ54M8dF8ui6KNzioxN8WVesO4CbjK/k5z4gSCfoWalHfxlzk/Kq
2V+/r6UTG5wIT0VHbvKY3GSa5To6ceH+tPlkAjcZsPq81CEZOtGU4UR4LHCT
K1XHRNkO0Ilt5ESkUjkPOasmxeZjTSd8GMF2+87wkf4FFQPKtnRCs2Xu3ENH
PlJgh3k5055OXJqmT9Hv8ZG7JB57WlylE3VGywlHK/hIh3men+o36USCZFmb
ojw/OVY4/3XhEZ1ISx483ckuQHbvbG2/8ZZOtOy6km7PFCDjpgW92RGdsHkp
tRi+IECeTzfcFlJJJ67zRXz6wCVIjsl0ezyrpRMHnJmJdtsEyVnJgQ11rXSC
KyDZ5MpFQZJnzYyN8DidyIiuVnLrFiQPTXAs5vAyCJ39xqvK24RIesqxVBUB
BsGM/ZDNPypE1lve168SYhClxUXv3H4Ikfr9fM87xRiEMJ1V7yW8ljTvENFZ
lGIQFkOlsp4ma0n3KulHx/YzCBSY1m04vJbMfAoyY+cZRGyrVEUKQ4Rkef7g
sLRkEHYjpfqnxUXI3SYvpvsvMQji281nYjtEyCJ+/oJuW3z/EWftKQMRkvQb
ONTgwiBGDVgdhU9EyC4r7wv5AQyCZ2KSxq0uStJky+L8chhEoTs9QiVEjDzM
7uSzlIe/Z9fsOu1kMTJocIuVZwGD0H700cKxRIzkjQne7lrMIAQlMmQ3T4iR
4lxmxTZVDIIruVyD7RD+25+e79TvYhA59JOT3xbFybLjdNnyRQYh75li2+Qv
Qf4PvTb9sQ==
     "]]}, 
   {Hue[0.9060679774997897, 0.6, 0.6], LineBox[CompressedData["
1:eJwB4QQe+yFib1JlAgAAAE0AAAACAAAAXH4s1P//778tKWz94iwCQN7yRxDz
9e+/hkRZ57ApAkBgZ2NM5uvvv99fRtF+JgJAZVCaxMzX77+SliClGiACQG4i
CLWZr++/9gPVTFITAkCBxuOVM1/vv8DePZzB+QFApg6bV2e+7r9SlA87oMYB
QPCeCdvOfO2/eP+yeF1gAUCqqRfUg8PqvzAc+uujggBAPYDfKG446L8fo0i6
N2f/P27XiAseuuW/XgjM6UbR/T/FGGPgswXjv3b6TwnuGPw/9SX3EH9/4L/e
HZok+H36P5c6eGdghtu/Is7kL5rA+D9/KsXITSfWv+N4XohbC/c/GbKF4aUk
0b/0VJ7cf3P1P/8bUL2Ta8e/4L3eIDy58z9fBvhMYo26vxtY5WBbHPI/rWZv
OSanmr/V7BrumYfwP6wDcj0/1a0/0hyi1uCg7T9Iv//2robBP5rCmshTbeo/
C1XfVXZwzD8XgpSa9vTmPyx0nL6Tk9M/kDbsBtiM4z8gptWagZLYP6pN0Gp/
X+A/yAOtkqP93T/w/GpdrdrZP99kiGkthuE/ziNO1Odr0z+0XYkXo0PkPyT8
ZhYE58g/69WoN1P05j9ICVPbZW62P0mCDvzNduk/AIYmzmXYiL+BREPOYi/s
PyDd004Ay76/4Dq+RMK57j9A6r8TYlLMv1HYqxaum/A/AAlPyyR/1L8fHuAR
iPXxP3Tpu8y4Xtu/AH43X0c48z9WAohvYOTgv85odrMTlvQ/PPawmDRe5L9O
k8RAfe31Pyb1eyfKx+e/4dc1IMwt9z9wkbq+mfbqv2GnjgYoifg/BhT4dTlq
7r/0kAo/ac35P/6Z1JqJ0fC/ObqVsEcL+z96L34t12Xyv2puCCkzZPw/HDgn
0Iwc9L+vPJ7zA6b9P24PCnfftfW/4ZUbxeEC/z/mWewtmnH3v5MEXnRSJABA
DnMI6fEP+b8O3rWiAsQAQLaR9VYqpvq/AH2BVDlxAUCGI+LUyl78v/uoXq/i
EgJABIQIVwj6/b9tmq+NEsICQKhXLumtt/+/uCuICBFuA0BomBIXmrYAwAxK
ciyCDgRAU+yqu6uCAcDXLdDTebwEQNL5QmjxXwLAq54/JOReBUCn7veWhS4D
wFivNhEd/gVAPmYVHwr5A8B7haGB3KoGQGeXMq/C1ATAqOgdmw5MB0Dpr2zB
yaEFwISdhWDeTgdAqou3Pl2lBcBgUu0lrlEHQGtnArzwqAXAF7y8sE1XB0Ds
Hpi2F7AFwIaPW8aMYgdA743Dq2W+BcBkNpnxCnkHQPVrGpYB2wXAIIQUSAem
B0ACKMhqORQGwPw4fA3XqAdAwwMT6MwXBsDY7ePSpqsHQITfXWVgGwbAj1ez
XUaxB0AFl/NfhyIGwP4qUnOFvAdACAYfVdUwBsDc0Y+eA9MHQA/kdT9xTQbA
uIb3Y9PVB0DQv8C8BFEGwJQ7Xymj2AdAkZsLOphUBsBLpS60Qt4HQBJToTS/
WwbAunjNyYHpB0AVwswpDWoGwJYtNY9R7AdA1p0Xp6BtBsBx4pxUIe8HQJZ5
YiQ0cQbAKExs38D0B0AXMfgeW3gGwAQB1KSQ9wdA2AxDnO57BsDgtTtqYPoH
QJnojRmCfwbAvGqjLzD9B0BaxNiWFYMGwJcfC/X//wdAGqAjFKmGBsB3tGLi

     "]]}},
  Axes->True,
  AxesOrigin->{0, 0},
  PlotLabel->FormBox["\"Graphical solution to Eq. 9\"", TraditionalForm],
  PlotRange->{{-1, 3}, {-6.626760134734051, 2.2719173239358015`}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{{3.407485155671875*^9, 3.407485164984375*^9}, {
   3.4074854665*^9, 3.407485487828125*^9}, 3.439910639421875*^9}]
}, Open  ]],

Cell[TextData[{
 "The root above was found using a numerical procedure, which was cumbersome \
in the days when numerical calculations required tables and had to be checked \
and double checked for errors.  In the paper by Peierls a simple \
approximative formula for the root is derived, using a Taylor expansion of x \
cot(x) around Pi. Using ",
 StyleBox["Mathematica",
  FontSlant->"Italic"],
 " it's easy to derive an approximation formula for Eq. 9. Expanding around \
Pi to first order gives (remember that ",
 StyleBox["b = 2/3 k \[ScriptL]",
  FontSlant->"Italic"],
 ")\n\n",
 Cell[BoxData[
  RowBox[{"1", "+", 
   FractionBox[
    SuperscriptBox["\[Pi]", "2"], "3"], "-", 
   FractionBox[
    RowBox[{"\[Pi]", " ", "x"}], "3"], "+", 
   FractionBox["\[Pi]", 
    RowBox[{
     RowBox[{"-", "\[Pi]"}], "+", "x"}]]}]],
  FontFamily->"Times New Roman"],
 " = 1-",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "b"], TraditionalForm]]],
 "\nRearranging terms gives\n\n",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "2"], TraditionalForm]]],
 " - x (",
 Cell[BoxData[
  FormBox[
   FractionBox[
    SuperscriptBox["\[Pi]", "2"], 
    RowBox[{"2", "A"}]], TraditionalForm]]],
 "+\[Pi]) + ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "A"], TraditionalForm]]],
 "(",
 Cell[BoxData[
  FormBox[
   FractionBox[
    SuperscriptBox["\[Pi]", "3"], "3"], TraditionalForm]]],
 "-\[Pi]) = 0, A = ",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Pi]", "3"], TraditionalForm]]],
 "-",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "b"], TraditionalForm]]],
 "\n\nThe solution is given below (the first root is physically relevant). \
Compared to the numerical solution 2.26, the error is ~ 3 %. The numerical \
value for b ~0.78  (U235), and the approximation gets better for b < 1 as can \
be seen on the graph below, or when \n\nb = ",
 Cell[BoxData[
  FormBox[
   FractionBox["2", "3"], TraditionalForm]]],
 "k ",
 StyleBox["\[ScriptL]",
  FontSlant->"Italic"],
 " = ",
 Cell[BoxData[
  FormBox[
   FractionBox["2", "3"], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   SqrtBox[
    FractionBox[
     RowBox[{"3", 
      RowBox[{"(", 
       RowBox[{"\[Nu]", "-", "1"}], ")"}], 
      SubscriptBox["\[Sigma]", "f"]}], 
     SubscriptBox["\[Sigma]", "t"]]], TraditionalForm]]],
 "< 1.\n"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{"Clear", "[", 
   RowBox[{"b0", ",", "Eq", ",", "x", ",", "sol", ",", "sol2", ",", "y"}], 
   "]"}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"tayl", " ", "=", " ", 
   RowBox[{
    RowBox[{"Series", "[", 
     RowBox[{
      RowBox[{"x", " ", 
       RowBox[{"Cot", "[", "x", "]"}]}], ",", 
      RowBox[{"{", 
       RowBox[{"x", ",", "Pi", ",", "1"}], "}"}]}], "]"}], " ", "//", " ", 
    "Normal"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{"tayl2", " ", "=", " ", 
  RowBox[{"Collect", "[", 
   RowBox[{
    RowBox[{"Expand", "[", "tayl", "]"}], ",", 
    RowBox[{"{", 
     RowBox[{"Pi", ",", "x"}], "}"}]}], "]"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"Eq", "=", " ", 
   RowBox[{"tayl2", " ", "\[Equal]", " ", 
    RowBox[{"1", "-", 
     RowBox[{"x", "/", "b0"}]}]}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{"sol", " ", "=", " ", 
  RowBox[{"Solve", "[", 
   RowBox[{"Eq", ",", "x"}], "]"}]}], "\[IndentingNewLine]", 
 RowBox[{"numSol", " ", "=", " ", 
  RowBox[{"sol", "/.", " ", 
   RowBox[{"b0", "\[Rule]", " ", "b"}]}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"Plot", "[", 
   RowBox[{
    RowBox[{"{", 
     RowBox[{
      RowBox[{"x", " ", 
       RowBox[{"Cot", "[", "x", "]"}]}], ",", 
      RowBox[{"1", "-", 
       RowBox[{"x", "/", "b"}]}], ",", "tayl2"}], "}"}], ",", 
    RowBox[{"{", 
     RowBox[{"x", ",", "0", ",", "2.5"}], "}"}], ",", 
    RowBox[{"AspectRatio", "\[Rule]", " ", "Automatic"}], ",", 
    RowBox[{"PlotStyle", "\[Rule]", " ", 
     RowBox[{"{", 
      RowBox[{
       RowBox[{"RGBColor", "[", 
        RowBox[{"1", ",", "0", ",", "0"}], "]"}], ",", 
       RowBox[{"RGBColor", "[", 
        RowBox[{"0", ",", "1", ",", "0"}], "]"}], ",", 
       RowBox[{"RGBColor", "[", 
        RowBox[{"0", ",", "0", ",", "1"}], "]"}]}], "}"}]}]}], "]"}], 
  " "}], "\[IndentingNewLine]"}], "Input"],

Cell[BoxData[
 RowBox[{"1", "+", 
  FractionBox[
   SuperscriptBox["\[Pi]", "2"], "3"], "-", 
  FractionBox[
   RowBox[{"\[Pi]", " ", "x"}], "3"], "+", 
  FractionBox["\[Pi]", 
   RowBox[{
    RowBox[{"-", "\[Pi]"}], "+", "x"}]]}]], "Output",
 CellChangeTimes->{3.388662977015625*^9, 3.38866306953125*^9, 
  3.4074848498125*^9, 3.407485487921875*^9, 3.439910640359375*^9}],

Cell[BoxData[
 RowBox[{"{", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", 
     FractionBox[
      RowBox[{
       RowBox[{
        RowBox[{"-", "3"}], " ", "\[Pi]"}], "+", 
       RowBox[{"2", " ", "b0", " ", 
        SuperscriptBox["\[Pi]", "2"]}], "-", 
       RowBox[{
        SqrtBox[
         RowBox[{"3", " ", "\[Pi]"}]], " ", 
        SqrtBox[
         RowBox[{
          RowBox[{
           RowBox[{"-", "12"}], " ", "b0"}], "+", 
          RowBox[{"3", " ", "\[Pi]"}], "+", 
          RowBox[{"4", " ", 
           SuperscriptBox["b0", "2"], " ", "\[Pi]"}]}]]}]}], 
      RowBox[{"2", " ", 
       RowBox[{"(", 
        RowBox[{
         RowBox[{"-", "3"}], "+", 
         RowBox[{"b0", " ", "\[Pi]"}]}], ")"}]}]]}], "}"}], ",", 
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", 
     FractionBox[
      RowBox[{
       RowBox[{
        RowBox[{"-", "3"}], " ", "\[Pi]"}], "+", 
       RowBox[{"2", " ", "b0", " ", 
        SuperscriptBox["\[Pi]", "2"]}], "+", 
       RowBox[{
        SqrtBox[
         RowBox[{"3", " ", "\[Pi]"}]], " ", 
        SqrtBox[
         RowBox[{
          RowBox[{
           RowBox[{"-", "12"}], " ", "b0"}], "+", 
          RowBox[{"3", " ", "\[Pi]"}], "+", 
          RowBox[{"4", " ", 
           SuperscriptBox["b0", "2"], " ", "\[Pi]"}]}]]}]}], 
      RowBox[{"2", " ", 
       RowBox[{"(", 
        RowBox[{
         RowBox[{"-", "3"}], "+", 
         RowBox[{"b0", " ", "\[Pi]"}]}], ")"}]}]]}], "}"}]}], "}"}]], "Output",
 CellChangeTimes->{3.388662977015625*^9, 3.38866306953125*^9, 
  3.4074848498125*^9, 3.407485487921875*^9, 3.439910640578125*^9}],

Cell[BoxData[
 RowBox[{"{", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", "2.317152817965588`"}], "}"}], ",", 
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", 
     RowBox[{"-", "13.81542284106877`"}]}], "}"}]}], "}"}]], "Output",
 CellChangeTimes->{3.388662977015625*^9, 3.38866306953125*^9, 
  3.4074848498125*^9, 3.407485487921875*^9, 3.439910640609375*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {RGBColor[1, 0, 0], LineBox[CompressedData["
1:eJwd1nk41NsfB3Dmq6sZ00aJynpJ4arubaPy+d7kKlERLRSSjCJFkXBL4uJa
QghlX6LFliVD5ivZGYxEGVG2rN8Y+wx+5/7O85znPK/nOX+cz3ne5/McBesb
xlcoQkJCcWj+f71n/csuG3vt+eX/BgkfWuUNzsofg2il6qj/vFxV9lBe/jxc
Z/5qvoS8Orwwh9ZhDaVFGwwFyHFEsNwK+atgSI9zmEdWJ21Clhk3wDGpuHga
uWt1/JKg3RnGz9vvnkBWSqfuXeC7gJJ4HncE2eGQy/VZubswkBCc14/85mNP
6pSOJ4RwFgu+Ii/YG3AnGPdh0Zky+Ak5UZgbMKv7AKJ0nuqykY80hinPtHvD
AaGalg/IP6L1yqeu+oCZhX8gEznIZvECj+8LPXSuczZyK98u4qdcAFjNHKiL
Qnatlt1J5v4L1Y56WgHImx5/rB/TCQK7vL52d2RrVRwbYTyCLjbnqTkyhe/e
NOUZCn/kKr7XR/7YYyozpRsGHCeKpCay2yux4sn2cAg/f+MvCWT98IFfJhMf
w6kzfMXlJRK2uJWbTFyNgDDu+t+GkQmdOz9JfiR4DI3VvkMW7fyuMiYXA+4m
SlvNkTuIdy6jP2LAZBP1Bo78Ij26YiQ3Fl75u/QrIZ9wPmE5rPMMhHjfzEYW
SYiiMqMGGQlQmUAyHJHz7+g735JKhFC3uS2GyM+G90nzPBPB4u6DZVVkx2Zx
xqRuElgVrTHsE5Ag/qyGMtGeDG0bl6iGyPzVhZk3D6YAZjSkqozc9yDl1M/E
FFgVfvyugE9CIeNePHk1FVTlb8ZnIJv9sVtrnJ8Gn7PX02YXSEiqS3AakcsE
BevQKw7zJExFF6h2B2WCi2xh1+/IfzHqeznzmVCw/aTF3BwJo9isKbP1BdAz
OG5eyPsPndAK8HsFzi++ifvOktCcLaCokNmw7daGz5bT6LxAtVilkQMyoxfb
JZEPbI3Q8rPJgeogXlfDFAkhU5k895YcOOfFx/Yi7w5ru2L9MhcU23a+EUyS
cL9OXX+X5RvgfMh4av6ThIzcoq0vI99A9sJ1GQyZE30YU254A3YaRs8zSZRf
xrlSKc18WBVg08obJ6EO89EQkiiAC32XgjzHSJA81CneVFUIjady7lkOk5CV
HdDpoMGEP2MD9O17SRgWHw7Sv8SEd6O/lvK/k6Diqq+9LYIJlo57TQOREw+K
JfXOM8GvqGDx+TcSwmuDbM0qS0BiqEilvRvltffRhN6Fd9B40p9L55KQIFb/
cJ0uARyHz6/bWknQPuRxwcuCgCeb6aLHkbmOanvIOwTwfff5sTgkSHECBxpe
EODQaCtIayEhLNrgmP/acthtmVls10TCQ2X26uXOcjhvV+xRWUcCA1piRp0r
oPetg6CTIEEoeYbZGVgBKT8sHDWRo0W2cOtSK6CMp7gmikVCTa2tbOanCliM
pHcZlJGwzUSQfEXrAxje7sjMKUHv007l1VfhSqD6NlOMCkm4Fv53WXNoFWSl
fgmaeoHe13RyNyuzCl4pf2HvQ449WyOU/b4KPHqWtN0z0X1vkdAJnqqCM28T
MwXPSVB7nlF97Fw1jDolPZ9NRfkpaW1+L1cDl29f/LMqHuW7f3tvflYtdMhH
E5VhJPRYfF5ZXlsLrXon3bihJJzu8Ndo7KsFyjb1/ZOPUN7qB+/2b6qDAOk9
zZtCSMBy09ZK+tWBZu9jK6t/Ub2eCtp3LOvBt2ZwodSbhHIJ6RjNtY2wFL2x
vcSZBJqCj0r3/mZIK3KQYp9E93EwbMc93WZ4Vu9jxDxBwj/n4vfJGDfDZrvU
nDRDVG/oWz1z+2YYWv7Uc/c46ndLo4z2uGbQvrtupZQeCSNc04wmSgsY7aVc
U9MmoSFaZTvR0AJBcQEHstRQvtfWqyZZtYLZzMcFBRHUrzpNWFjDJ9g7sEIy
LHkcnCp9rG+ofAHb5JfOLb+PQ+5N+mmpY13ALK08IlU6Bkke9sM7z/bA/raS
o0fUxsA7QWA0NfoN5mJSxNMjRkGotQcvde8FQ0w/SHNyBLbetnUzUe0HBbr0
G/quEeibfnJG7N0AeEyxB90chkG//9rf4U4/YDoZ913lOwQn+kLVZGWHIYZM
yikI+QHDFdZeNa9HQFo80FL3ySBorPsko2M5BuA9UvvVdQBoMH2JyxsHSXWJ
s8pn+kFvzGxW0fQncMSWfdaL9kHHexthM/UJoN6S+fdZ+nco8P9gIT89AY96
6o9W630D2xUuPWz2JLjLyV01tumGy84eBfRYHrgPNq6ZOciFHQw1tsGlKYhc
vlZFmeqAJnF1MRfNaVDkKRas3d4GJVneT9WxGWD7Hz7lbdMCHsKqpxq6ZkC2
3emrhGQd0LjuRx3zZ8Gm6KTv73ME9KcrRCR4zsGFjvNrYtOzIWKbtEix6Tys
L9V48tojghUaa1t9SX4BqtJXzCstF7PMIi6y6YMLsEVm8wbjwWrW3JG0+YQy
PozMtpnc1uOwevsPu733E4DbXGRRUWwHK/nQrs2NlxYhjyu2a2i2m/Xa12Ud
47clsEnztgvN7mX5mV4+jk0vwaYw13vdUYOsPdykTL/GZQh3LRAtVx1lRW9/
+rrxuBDutufql2Tbn6z7RHD08WIhPMTO4Mi8Io9luNk0rU5JGLf1GTwmtmGG
9adjtfm6AGGcKeokOiI1z0rYYq3UMC2Mu6fVL6XnC1hKFc+jTphTcCvjAtXd
mkJE6G/vVTllFHwj8yHmokMhbHkiLbWqGK5lGPVPHiFCWOdOygc8wvCJi02J
7pdFCQfv0wL9MAxP2tGfv/KmKOFyOr+d/hjDldob1eI9RQn/aZeQ0Ci03yb4
YnekKJGlOc9/Eofh9yMHwktqRIn58uVP6S8x3MtJ59a3HSuJUA49+EM1hifU
fbeWFKESsSnX7f6pxXDKPeW6+HVUIuU2W+doPYZbrAx+sEOOShRKhi7UszF8
WYIufE2LSnDNJOxa2zD82OXhrUecqIRKr7TO914MZzGmjCN6qMTOfHfZ1H4M
f+CQEZs2TiW0fDvnrwxiuFz6329ZAiphoBKXMzSM4atKnsSslqYRzvbyshMT
GL7dzHVoozGN8Dz4YD6Ph+FN8ezPN61ohO+q7x9vT2M4OzdqkONII6KzUwLn
5jD8Y6SIV14gjUj2EmEwFzC8OmutzM5YGvHS6MphTwGG25w6LyjMoBH5ilUy
2ksYvlJml4peEY0o421FXzcM73KdTe+ppBH/A0NkJq8=
     "]]}, 
   {RGBColor[0, 1, 0], LineBox[CompressedData["
1:eJwVjns0lGkcgL/Pp1Q4B51ts8ToQsYlzq4jNn5vaXdzXVbKysZMuRRNpZld
oVISKq0yRy5ZK9R2Wpes9qRa7xtKxCDbsGHDGCSMMQbjNtv+8Zznr+ecx4x7
/LswDYqivD/yv/PPclfaH45ydTew7FWrZVDfwfLaz3KHu2JVicacDNQvapJY
rO9hZV4Tb8WkDPJJuukK1hGo3xRNaffJINpFcGzW9DQgh/JN64gMOhYihZOm
aeAY6HDUOlEGWt0DFuOmOWD+y87ZQEoGhU0FJz+Y3gOGHaRfSU1AWXlad7Tt
Y2iVX7Hnij9ABLTnjMXUAfwk0BG1vQee1FJSVdYIzQ//0eK8GYY1Zhct3m1v
g6hiV5m8TArX9F6xC0M7oMz3wHXXUglwu/diplkMGqWxpmmJ/XDy+UXucYu3
MHTR7Jbyyb/w4ISO/3r3XtCzesbzJW+hMD5q1G5/HwTbcw3v8sRwoWDRb3qs
HyYrY5wE4a+B6uhDT+MkkOgt/2KTtBnM+eGxe9lSOKxU3/H3q4dB5c192n8N
wSOLbq2hXx+Bh/TomRsnR4AIZZd8jfPBZzDDysRkFPq2Gp1fnVSER+u4iS9L
P8Aky7AFiquxrb54g1vIOJzIan2afPw5XgNKTo9iArpKPv9t8XQL/mY8aHZj
wCQkxjXJXdZ34K7aw3SQtRy6ypN9toV04oep9QdZSjkMHml3uBfcjcNXCPpE
oiloqzP6yqLiHT4UE/9QJ1cBzn7xsnW9/XhbhJXIizMNPK0JI0WnBLcaWGsL
nJRwQzc0uH1Jip+UXcizZmZAb99nA4rlYRxPs32be2fg7ERH/HqrUbymJ24P
r2oW2loMlJ94j2HpHTNhQcIcJMhy/tjvMYGFWw01qwNUcNVOr87baRJn5IY3
cFjzkHp5VqT/pRwHCX8Q6QzPw3RKoUet4xSe212iKqhZgEAbzzaFuwJLpLti
a1MWYUyRHdvgPo1vu9gbtXCW4EfuENs4RIlLkwX6ETbLYDvVqTbgzuCUgEOe
jHIZArdc3T0TMosdegrvpbSoIY5CmUGn5nC2ZV5piyeFmhZvanL4KnyOpGd7
VlNIr1B2buTyPPY2Cihp2kyjlyOro4JvLeCdvIYD+mk0+t3vZWxM3iIuMOZu
blbSKJV7KWN1xRLeXHc3y+eABhIKR5p2VS7jDJta9usaDaR3VFcQWa7G4QrN
9kY2g2oV2UItDkW4D6ZYaT8zaK3bktqqnyLRF/wXPa4zKFrHyzZSQhGBf1Wn
TiaD6t/kBhdLKZKqFFzLyGKQIHJ7tfEoRcqcVAs38xnUlX4qRldBEdUztfjO
fQbd6no/OKFJk4zXOun1DQwy54kbH2ylSW7RschLjQxKdNwyN86mSRFf5Lbn
1cee4puzbWjy57qM+VciBqVl6ifdtqdJT9DayI43DBp75LUj05kmFhJDtwEJ
gyqYujK+F03squJMiqUMWtWi31vhQxPn5G5V2DCDOFmh2uO+NPGyyK94P8og
A8vliLAAmsREsUzkcgbxfZzM9h2kScKO86pKBYNEn6Z+eyOUJsm6A3/zlR9/
+8VnRFyaZJcXXZmbY1DnKf7bryNocjtRM+LxPIPsXOpWJR2hyX2/sF0JiwxK
XWngiKNoUrXxxQbXZQb1t4aGLRyjSY3CXKVWM8g5pzzT8QRN/gM5nmva
     "]]}, 
   {RGBColor[0, 0, 1], LineBox[CompressedData["
1:eJwVknk0FHgcwGf8GDNjei1t14oZrZdIOl7bYzu+36XX4SjSbKXTOHdHck3Z
CBWKqKmsJIwrKktLap+y5rfYFJKS44XN0uT2DAbj3PaPz/v8/3kfI9Hp/R4a
DAbD4Qv/OzVMxNrgLt7O3aI+lOPExcoGgf1BwR6QCV+4XdjHxfkXZZcEgsOQ
o3fr+FEHLqbSeL6W4Ceo8Ocb6+3hos82yakJ/i/wTeyjoVDkYsO0d8IwPwaO
tbDD963jonZrp8kg/w5Y5grP9vC4mFEt8+/nPwBPh/PLe6s4WPAoptXH4hls
H5K4em7goBe8vTMQUAFO3IDwwiw2+ipMu4oLXoHHoIEyZSkbuUaRJh8t62Ev
d3iBQ5g2XvuqxizjZANom5r7dQyyUNR6QE5qm0BwMiWyfh8L/f+OFJ02+QBi
hXrzdJEWFvrxnJftaQeD6ksiAV8LM0LEfesPdkBbsX3eXKgmXpTNOI0N/AuC
M0mFzV0EGQ0dWHquC8LHLcb77QiuCvIMPmCmAI/rNf6T2Rr4SXX7R50/P0Nb
hOk/X2troK3i5/M3/XvA5s0x29DDTNz7SbrG0LAPQnStTol/Y2BfhSjiZX4/
pO/kDIexGWih22Rgc2IQXNtTtkbnzgEXVK5to0PAzI9I2+I4C7sGXSZWCoeh
3Ktl1cHGaWgpd2e6mCshk3V0pOHwFDy5UnlcoFJCl3jZ3NXhSfDUknTU1Y3A
gGtvmvDkBLgFhDzhJY9CZtibYZ1OFazzWlNn7zoGcZmPS6NFY/BGz1xHYqWC
wEWx+u1/jMDzgot3zck47Np6NrLbUAkhTDPH2vZx2G/L1MnNHwJu27ndvsUT
ULQ5x/K5bj8ocowSZKGTYKhVy1ud1A0Jq5drlgjV4Kv/luoXdYE02bPKVTAF
5tWOntPNH8El4Vgdr3sK/ExSemPMWmByxz21rGwaOj51x6mm6qBLYR1cfnkG
/Gx1FW7GpZC5bYP+a9dZiEi/FnUrLlueHyXR9Vo7B/X3CzcFzlfJLwvd7Ihq
DoLTqzTFzxvl37VlPLj8eh4cdz8NNArqkCeZ3s1/bcfAYLsw60vzCnk4jU+y
K2Gg+1mrDEvxgNxBX3iv2piJG914LzclK+U/+FYd0Y1hYv2NANmvmiq5bIXI
uFbFxLSSlIw8rlpuXJGbuPeIBhaV2omMTWbl0rXlZu/KNHCcMyPTq2ZQz1HN
t6/MCAqbU2KbgggVFY4IYq4TbPFfHPgMWNTnovOM7Q2CihKxy7g9i0qci5t5
twieseGdsXJh0SsqyTVpIsGPeNOjMYhFC6zU07dTCS5M9+h2f8ii6r/mm3Ly
CHIez6b3LNGm0ne8+MoqgvHSyNVPVdo0OeuUd/QrgmPMA7GfCZtmBdXZ7K4h
aKXGBXw9Nn26RDpVU0cwozY78b4Fm7a5LPJuaCS40dHmrqE3m5p0Lbfp/PK9
t8WK/NQ2Nl1ffM4wW0EwT0dj1rmfTb+PalV7dBN8956O6k2xqb1J6u+9fQQX
V1dOPlzKoQFigaFSSbB584lDO5w5NHTrBXXRKMEP+puaUMShUQs63wepCJpm
LrTe6c+hSY+yrk5Ofulz/BDb5zqHZkZoej2bIjiT+K0kLo1D85w8rENnCPqk
FY49yefQ4pUvDLbPEXysdJT2lHJo2egq9fw8QV7tUSfjWg79Dxc5TEc=
     "]]}},
  Axes->True,
  AxesOrigin->{0, 0},
  PlotRange->{{0, 2.5}, {-3.3466198963422884`, 3.289868064027705}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.388662977015625*^9, 3.38866306953125*^9, 
  3.4074848498125*^9, 3.407485487921875*^9, 3.4399106406875*^9}]
}, Open  ]],

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

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{
   RowBox[{"A", "=", 
    RowBox[{
     FractionBox["\[Pi]", "3"], "-", 
     FractionBox["1", "b0"]}]}], ";"}], " ", "//", 
  "N"}], "\[IndentingNewLine]", 
 RowBox[{"Approx", "=", 
  RowBox[{
   RowBox[{"N", "[", 
    RowBox[{
     SuperscriptBox["x", "2"], "-", 
     RowBox[{"x", " ", 
      RowBox[{"(", 
       RowBox[{
        FractionBox[
         SuperscriptBox["\[Pi]", "2"], 
         RowBox[{"3", " ", "A"}]], "+", "\[Pi]"}], ")"}]}], "+", 
     FractionBox[
      RowBox[{
       FractionBox[
        SuperscriptBox["\[Pi]", "3"], "3"], "-", "\[Pi]"}], "A"]}], "]"}], 
   "\[Equal]", "0"}]}], "\[IndentingNewLine]", 
 RowBox[{"sol2", "=", 
  RowBox[{
   RowBox[{"Solve", "[", 
    RowBox[{"Approx", ",", "x"}], "]"}], " ", "//", " ", 
   "FullSimplify"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"z", "=", 
   RowBox[{"x", "/.", "\[InvisibleSpace]", "sol2"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{" ", 
  RowBox[{"Plot", "[", 
   RowBox[{
    RowBox[{"z", "\[LeftDoubleBracket]", "1", "\[RightDoubleBracket]"}], ",", 
    RowBox[{"{", 
     RowBox[{"b0", ",", "0", ",", "1"}], "}"}], ",", 
    RowBox[{"AxesLabel", "\[Rule]", 
     RowBox[{"{", 
      RowBox[{"\"\<b=2/3kL\>\"", ",", "\"\<Numerical value of root\>\""}], 
      "}"}]}], ",", 
    RowBox[{"BaseStyle", "\[Rule]", 
     RowBox[{"{", 
      RowBox[{"FontSize", "\[Rule]", "12"}], "}"}]}]}], "]"}]}]}], "Input",
 CellChangeTimes->{{3.407485212140625*^9, 3.4074852791875*^9}, {
  3.40748530965625*^9, 3.40748536784375*^9}}],

Cell[BoxData[
 RowBox[{
  RowBox[{
   FractionBox["7.193832906510146`", 
    RowBox[{"1.0471975511965976`", "\[InvisibleSpace]", "-", 
     FractionBox["1.`", "b0"]}]], "-", 
   RowBox[{"1.`", " ", 
    RowBox[{"(", 
     RowBox[{"3.141592653589793`", "\[InvisibleSpace]", "+", 
      FractionBox["3.2898681336964524`", 
       RowBox[{"1.0471975511965976`", "\[InvisibleSpace]", "-", 
        FractionBox["1.`", "b0"]}]]}], ")"}], " ", "x"}], "+", 
   SuperscriptBox["x", "2"]}], "\[Equal]", "0"}]], "Output",
 CellChangeTimes->{
  3.407485270515625*^9, {3.407485319859375*^9, 3.40748533209375*^9}, 
   3.407485368671875*^9, 3.407485488109375*^9, 3.439910640828125*^9}],

Cell[BoxData[
 RowBox[{"{", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", 
     RowBox[{"3.141592653589793`", "\[InvisibleSpace]", "-", 
      FractionBox["1.5000000000000002`", 
       RowBox[{"0.9549296585513721`", "\[InvisibleSpace]", "-", 
        RowBox[{"1.`", " ", "b0"}]}]], "-", 
      RowBox[{"0.5`", " ", 
       SqrtBox[
        FractionBox[
         RowBox[{"9.869604401089358`", "\[InvisibleSpace]", "+", 
          RowBox[{"b0", " ", 
           RowBox[{"(", 
            RowBox[{
             RowBox[{"-", "12.566370614359172`"}], "+", 
             RowBox[{"13.159472534785806`", " ", "b0"}]}], ")"}]}]}], 
         SuperscriptBox[
          RowBox[{"(", 
           RowBox[{"1.`", "\[InvisibleSpace]", "-", 
            RowBox[{"1.0471975511965976`", " ", "b0"}]}], ")"}], "2"]]]}]}]}],
     "}"}], ",", 
   RowBox[{"{", 
    RowBox[{"x", "\[Rule]", 
     RowBox[{"3.141592653589793`", "\[InvisibleSpace]", "-", 
      FractionBox["1.5000000000000002`", 
       RowBox[{"0.9549296585513721`", "\[InvisibleSpace]", "-", 
        RowBox[{"1.`", " ", "b0"}]}]], "+", 
      RowBox[{"0.5`", " ", 
       SqrtBox[
        FractionBox[
         RowBox[{"9.869604401089358`", "\[InvisibleSpace]", "+", 
          RowBox[{"b0", " ", 
           RowBox[{"(", 
            RowBox[{
             RowBox[{"-", "12.566370614359172`"}], "+", 
             RowBox[{"13.159472534785806`", " ", "b0"}]}], ")"}]}]}], 
         SuperscriptBox[
          RowBox[{"(", 
           RowBox[{"1.`", "\[InvisibleSpace]", "-", 
            RowBox[{"1.0471975511965976`", " ", "b0"}]}], ")"}], "2"]]]}]}]}],
     "}"}]}], "}"}]], "Output",
 CellChangeTimes->{
  3.407485270515625*^9, {3.407485319859375*^9, 3.40748533209375*^9}, 
   3.407485368671875*^9, 3.407485488109375*^9, 3.43991064184375*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {Hue[0.67, 0.6, 0.6], LineBox[CompressedData["
1:eJwVjmk8lIsegIejEiIpa44z1FhHFBHlX+mgcGwVRZJQsktIRZJUSnGRNOJI
13YY2RLl/Tdki0QRZaeZd+ZVUeEo6XY/PL/n+fjQPYMcvYVpNNr5X/zfe7zJ
LuTvM/uVxu2yvvUi5orBwwoGQKsJtFdfaklsys5gVSjsAprVlpMBos6Edceb
8iyFvUBb3XTEYsUxwtlr7kmKghfQAtjhanIRxOCGBPPrCmFAa6gbW6qdQOSs
sbt0VeEi0CJNCX2nW0SsaIX6FYVUoI2FRfzIyydcDvzNnprOg7J1A8uuHX1E
iP+MdbP8nQ3nj/WxnY83EOsdPpINiZUQ/NeClvvhNiLixRbI+PoIgn3CFVZt
7CLkKqWHhIQIwGaXIhWdHqJ6zLy9cysHJuujNK5y+gnOiSylA5ONgEatC+XJ
Q8QMY5tUSUAzFMXazgqVjBLC88Urlda3gUCpLX40ZILY/NOj2qm0Hca+G7x4
Ysgj2i8avGjS6oSv06FhLYoCQuUmNV/48SWUdV9hB0d+IKi54aEevW5IoZ6v
7PSbIvIEdW2SSa9g39w9SeWjnwlJ5fxPExOvgRXXXdQc+ZW45lvlk7OuF5be
5NyPSZ4lZrcXTxyPeQO1Un72jUvnibCNjZWyDX0w29rZox33nRgubUzsor8F
N5/dAf7Bi8T3E6/vaEa9g4twgN3ZSkPtHcsWfeoGIG/bidojnsJomUs/9kBu
CLgnzaXlLEUwQvmw6qjbMLQnF52bM1qKOpP7xN/YjsAPJc4Zq79EUUZqSrTx
2wisruaoWr8QQxHiDluNNQoq40y3bK8V6GdQ17/cegxiBhNrgnauxOLdnjn1
s2PQ1pTJkFqQxhlKcSU3fRxK+3+74/lcBqucLH9L3TUB7qY3e++LyGKMVXx2
JXcC1hccibE/Ko/0wT199+LfQzXT5QFLXAlNjC6mKmzkQl3OflbJ0d/R7ffO
eaFXXBiK3k+YLqWjpEJoVV0PFxRWBBT2StLxqcyakFN9XMgdsGiJkKMjY7kr
SQ5yYbX1nfp2DTp++sLt7SS50LZ5sxXHmo4XWhcrsha5UDvJ4Iym0LEgTDdw
ixYPXp9qmVVjqOJM2/Wx4Fge2NwcvOLjqYYD4j3HF+J4cL/kTJ28vxpybNZ+
SrjEA7kqS6nuU2p480XR97uJPPg3LnDrwStqqNXdvLojjQfi9NTLT8vU0L1f
yFKziAfR8hrUZ9o6bOKdKh7t5sHL27f9lfLXYYaIe5jjOhIOiB+6EifBwMWj
coZdDBKszssbicoz0Kvh5YydJgmyZ0t+pqkxUO+CeYStLgm+H2UlOkwY2PJT
I8rKmISD9qe3F/kycO7blxgzGxJCEm80arUycP/0lUStMBLubu92Truhjo/t
zW0Kw0k4fePYISeWOqqWLUhonCahF4dnlArV8WNgUNL6aBKipe7KdHDUMX5y
b/Ifl0losloZVjarjpU8lVuyLBIs+gKyEzw0UGaoKle4kYTAmiXK77ZqYpLE
uXMPm0gItkszNbHRRDHTXS7+rSS0dBVT+a6aKJzRvaL3BQlJOXmO7ChN/Ozw
KbKwn4S1REdc8CNNfPlMw87+EwkTHBGjzi1aeK3kzvcsRT58u898zLLQRtEB
z15HZT54FBl8+NdZGy+KaT1Y9gcfaFO3U3x8tfHMsRqfkPV8cF2rcy/8mjae
oPd07dLjQ8Hjigb/V9polSZZQP3Jhz0n6XMtXjq45Fzs3i0hfNAoHGBLpDOx
cnpynjjJh9GRSYezBUz08nHJtgjnQ/r14cOLtUxstNsgcDrDh6QYV6bJCBMv
qg7GBMbz4dCIJes/2roo0mJcfO82H0wZQdPcZ7ooLDMlJMnhw9bquaTPInpY
nuCan9rIh3evkR8gpYeeC002Ss18GGOO+X1V0EPO+6xbGu18uMrRZKpt0MML
NdZM814+5PbTspYf1EMh93yXSAEfnGvtzoqV6yGt0J09tkoA5y1/qB/z08em
yTeB2WsEMIWlHwIj9fGanoOum7wAFguSmi/F66NcjXlJj7IAFG1qdHnZ+qjb
rFHcoiGAEeM5xcAefXR7//m/pWYC0LK7x2ebb8RH9IS7UScE0CSR5zbF3ITR
3jR3owAB5PYeUTfZtgl3FZ5W/hokgJ1d+a4ZNpvwpZ4/K+CUAOJq7fNS/DYh
CQ6ZHucFsNmvzyKuaBPKuiulW6QLgJQ++zJKxwBDM9nXV3EEsHealyWx0xAf
LmQmcBoF4J+69y7XyRAX3C/FhjYLwOrdBeW33oaYoHboVHe7AD46jodKXjVE
1j9ih1LeCCC09PKt8leG+KzeW2fVBwF0NOr2KvlvRvlxpefS8hQINMlHGyqM
8NCfy549VaQgx27QrazFCHPzP9eHKFPQ1bSm1HHICHX8W8u7VCmIHijN+LHc
GHfMhGcmMynYNW11csbTGE8s6/aVNqegbyi17aHSFnyifVlUOpCCxvb+ELli
E9z8Wq4jJZgC4xUewlUNJlh2Nj959UkKvrrVekYNmGBuR5OifCQFZfXZstGS
png5aAlT5QIFq6BXKyzMFB0rLjgw03/9yevZhVttRZ7Juczd9RTs0JddMy9j
hh7jEoefIwVBrqy/fTaa4dtElpptAwVpKtdGxRzMsH2g7h/7Fgrs115un08y
w/Lob/XO3RQ8sh3/HrsC8CwnYtybS0F15i3zCbHtOL02Mcd2hgL26eEvzuHb
8X8nCBMa
     "]], LineBox[CompressedData["
1:eJwV0XtI01EUB/AdLZyaZiFp0AMjcUpPKrGs7mk2E1/g0BaZaKY4rRFl5iTT
hCYm0U9TU1eKMl3D8hEqmUw0ylBIzVca6SwDmaC7RYthpev0x+Xy4Xw55x6u
T/JVeaqDSCTypfP/1pZJRbOChbXmzP1U3MQ+1wmPKhPZW3V82DsPMCV/rKL0
IXnIsabuNqAxoLzsdAV5UqmWkD0/xpXaKsm6nqS2XEBVgZfQpCWv/o0LIu+Y
1hZ71JEdJWlhtwDvaOoLZposrKrPb3d6DqBbtlP4xWbynqAbK2pArVK1eaGV
8sfyRu6R2yOO6ngH1QsNT5uyARc2jb516CU3JJqWsgAjakVi/3HKl7WnZF4H
3PIy8UHWHws719A1UnQZsNm8/bF7JGfT/c+tEwmARTWaJFM0Z0WpquR4cop8
2bclhurui7NfLwBuMxrbohSciVc2rFniAYuF+Hf3k6kumOedzwMqA6t/uORw
pvSznZCeBdx11zPMSU/9rPqDPdGAa0G5blMGzgb09bEy8qflb2P6Z5x16Q9p
hqIASxTtCaEvOKsSb91oigS0B8gzNUbKawb77eGAM6NCreM4Z+ap1z2yM4Bd
hbZL45OcSXoNr4ZDAcuDE/110zS/xfONghzRuK9DauIszCV7LV0G2K0eGixY
pHyazEsIAXy094gQs0T9Qw7ne5OvzT+J9eGcGcy11nopoCTyylyfld73fd6p
8xTgOphoKLHRPtDdeJL8pTM4I+k3Z+pAccwA0v9n6PYfWKV5ae+d5eTKna6/
7HbaJ2P9h88M8B9jaxge
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesLabel->{
    FormBox["\"b=2/3kL\"", TraditionalForm], 
    FormBox["\"Numerical value of root\"", TraditionalForm]},
  AxesOrigin->{0, 0},
  BaseStyle->{FontSize -> 12},
  PlotRange->{{0, 1}, {-32.924559101218044`, 2.186570750445423}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{
  3.407485270515625*^9, {3.407485319859375*^9, 3.40748533209375*^9}, 
   3.407485368671875*^9, 3.407485488109375*^9, 3.43991064196875*^9}]
}, Open  ]]
}, Open  ]],

Cell[CellGroupData[{

Cell["Peierls' formula for critical mass", "Section"],

Cell[TextData[{
 "In Peirls paper a simple formula is derived for both the low and high beta \
regimes and he produces a graph of the critical radius valid for both regimes \
of approximation: \n\n",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Beta]r"], TraditionalForm]]],
 "= ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    FormBox[
     SqrtBox["3"],
     TraditionalForm], "\[Pi]"], TraditionalForm]]],
 "x + 0.71 ",
 Cell[BoxData[
  FormBox[
   FractionBox["3", 
    SuperscriptBox["\[Pi]", "2"]], TraditionalForm]]],
 " = 0.552 x + 0.216",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "2"], TraditionalForm]]],
 ", for x = ",
 Cell[BoxData[
  FormBox[
   SqrtBox[
    RowBox[{"1", "-", 
     FractionBox["\[Alpha]", "\[Beta]"]}]], TraditionalForm]]],
 " \[LessLess]1 (from Taylor expansion of x cot(x))\n",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Beta]r"], TraditionalForm]]],
 "= 0.78 - 1.02 (1-x) , for x ~ 1 (from Eigenvalue calculation)\n\n(remember \
that",
 StyleBox[" \[Beta]=n(",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "s"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["+\[Nu]",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "d"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[")",
  FontSlant->"Italic"],
 " ). The first part of the formula can be derived as follows: From eq (7) \
and the solution of the diffusion equation (3)\n\nrk cot(rk) = 1-",
 Cell[BoxData[
  FormBox[
   RowBox[{"(", 
    FractionBox[
     RowBox[{"\[Alpha]", " ", "r"}], "b"], ")"}], TraditionalForm]]],
 " where b = 0.71 \n\nTaylor expansion around Pi yields:"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{"Clear", "[", 
  RowBox[{"x", ",", "r", ",", "k", ",", "b"}], "]"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"Ptayl1", " ", "=", " ", 
   RowBox[{
    RowBox[{"Series", "[", 
     RowBox[{
      RowBox[{"x", " ", 
       RowBox[{"Cot", "[", "x", "]"}]}], ",", 
      RowBox[{"{", 
       RowBox[{"x", ",", "Pi", ",", "2"}], "}"}]}], "]"}], "  ", "//", 
    "Normal"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{"Ptayl2", " ", "=", " ", 
  RowBox[{"Collect", "[", 
   RowBox[{
    RowBox[{"Expand", "[", "Ptayl1", "]"}], ",", 
    RowBox[{"{", 
     RowBox[{"Pi", ",", "x"}], "}"}]}], "]"}]}]}], "Input"],

Cell[BoxData[
 RowBox[{"1", "+", 
  FractionBox[
   RowBox[{"\[Pi]", " ", "x"}], "3"], "-", 
  FractionBox[
   SuperscriptBox["x", "2"], "3"], "+", 
  FractionBox["\[Pi]", 
   RowBox[{
    RowBox[{"-", "\[Pi]"}], "+", "x"}]]}]], "Output",
 CellChangeTimes->{3.388662980078125*^9, 3.388663218390625*^9, 
  3.407484855328125*^9, 3.407485489671875*^9, 3.439910642734375*^9}]
}, Open  ]],

Cell["\<\
Now discard the quadratic term to get a managable polynomial equation of \
second order and set x \[Rule] r k. Together with the RHS of the equation \
above this gives: \
\>", "Text"],

Cell[CellGroupData[{

Cell[BoxData[
 RowBox[{"PExpr1", " ", "=", "  ", 
  RowBox[{
   RowBox[{
    RowBox[{"Part", "[", 
     RowBox[{"Ptayl2", ",", "2"}], "]"}], " ", "+", " ", 
    RowBox[{"Part", "[", 
     RowBox[{"Ptayl2", ",", "4"}], "]"}], " ", "+", 
    RowBox[{"\[Alpha]", " ", 
     RowBox[{"r", "/", "b"}]}]}], " ", "/.", 
   RowBox[{"x", " ", "\[Rule]", " ", 
    RowBox[{"r", " ", "k"}]}]}]}]], "Input"],

Cell[BoxData[
 RowBox[{
  FractionBox[
   RowBox[{"k", " ", "\[Pi]", " ", "r"}], "3"], "+", 
  FractionBox["\[Pi]", 
   RowBox[{
    RowBox[{"-", "\[Pi]"}], "+", 
    RowBox[{"k", " ", "r"}]}]], "+", 
  FractionBox[
   RowBox[{"r", " ", "\[Alpha]"}], "b"]}]], "Output",
 CellChangeTimes->{3.388662980375*^9, 3.388663218484375*^9, 
  3.407484856484375*^9, 3.407485489703125*^9, 3.43991064325*^9}]
}, Open  ]],

Cell[TextData[{
 "Solving for the reciprocal of r, and using the wavenumber from Eq. (2) ",
 Cell[BoxData[
  FormBox["k", TraditionalForm]],
  FontSlant->"Italic"],
 " ",
 StyleBox["= ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SqrtBox["3"], TraditionalForm]]],
 StyleBox[" ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox["\[Alpha]", TraditionalForm]],
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SqrtBox[
    RowBox[{
     StyleBox["1",
      FontSlant->"Italic"], 
     StyleBox["-",
      FontSlant->"Italic"], 
     FormBox[
      FractionBox["\[Alpha]", "\[Beta]"],
      TraditionalForm]}]], TraditionalForm]]],
 StyleBox[" = ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   RowBox[{
    SqrtBox[
     RowBox[{"3", " "}]], "\[Alpha]", " ", "x"}], TraditionalForm]]],
 " ",
 StyleBox[" ",
  FontWeight->"Bold"],
 "and Taylor expanding around x=0 yields:"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[
 RowBox[{"\[IndentingNewLine]", 
  RowBox[{
   RowBox[{
    RowBox[{"Psol", " ", "=", " ", 
     RowBox[{
      RowBox[{"Solve", "[", 
       RowBox[{
        RowBox[{
         RowBox[{"PExpr1", " ", "\[Equal]", " ", "0"}], " ", "/.", " ", 
         RowBox[{"k", " ", "\[Rule]", " ", 
          RowBox[{"\[Alpha]", " ", 
           RowBox[{"Sqrt", "[", "3", "]"}], " ", "x"}]}]}], ",", "r"}], "]"}],
       " ", "//", " ", "FullSimplify"}]}], ";"}], "\n", 
   RowBox[{
    RowBox[{"Psol12", " ", "=", " ", 
     RowBox[{"r", " ", "/.", 
      RowBox[{"Psol", "[", 
       RowBox[{"[", "1", "]"}], "]"}]}]}], ";"}], "\n", 
   RowBox[{
    RowBox[{"InvR", " ", "=", " ", 
     RowBox[{"1", "/", "Psol12"}]}], ";"}], "\n", 
   RowBox[{"TrueForm", " ", "=", " ", 
    RowBox[{
     RowBox[{"InvR", " ", "/.", " ", 
      RowBox[{"\[Alpha]", "\[Rule]", " ", "1"}]}], "//", " ", "Simplify"}]}], 
   "\n", 
   RowBox[{
    RowBox[{"ser", " ", "=", " ", 
     RowBox[{
      RowBox[{"Series", "[", 
       RowBox[{
        RowBox[{
         RowBox[{"1", "/", "Psol12"}], " ", "/.", 
         RowBox[{"\[Alpha]", "\[Rule]", " ", "1"}]}], ",", 
        RowBox[{"{", 
         RowBox[{"x", ",", "0", ",", "2"}], "}"}]}], "]"}], "//", " ", 
      "Normal"}]}], ";"}], "\n", 
   RowBox[{"PeierlForm", " ", "=", " ", 
    RowBox[{
     RowBox[{
      RowBox[{"ser", " ", "/.", 
       RowBox[{"\[Alpha]", "\[Rule]", " ", "1"}]}], "//", " ", "Simplify"}], 
     " ", "//", " ", "Expand"}]}], "\[IndentingNewLine]", 
   RowBox[{"Plot", "[", 
    RowBox[{
     RowBox[{"{", 
      RowBox[{
       RowBox[{"PeierlForm", "/.", 
        RowBox[{"b", "\[Rule]", " ", "0.71"}]}], ",", 
       RowBox[{"TrueForm", "/.", 
        RowBox[{"b", "\[Rule]", " ", "0.71"}]}]}], "}"}], ",", 
     RowBox[{"{", 
      RowBox[{"x", ",", "0", ",", "1"}], "}"}], ",", 
     RowBox[{"PlotStyle", "\[Rule]", " ", 
      RowBox[{"{", 
       RowBox[{
        RowBox[{"RGBColor", "[", 
         RowBox[{"1", ",", "0", ",", "0"}], "]"}], ",", 
        RowBox[{"RGBColor", "[", 
         RowBox[{"0", ",", "0", ",", "1"}], "]"}]}], "}"}]}]}], 
    "]"}]}]}]], "Input"],

Cell[BoxData[
 FractionBox[
  RowBox[{"6", " ", "x", " ", 
   RowBox[{"(", 
    RowBox[{
     SqrtBox["3"], "+", 
     RowBox[{"b", " ", "\[Pi]", " ", "x"}]}], ")"}]}], 
  RowBox[{
   RowBox[{"\[Pi]", " ", 
    RowBox[{"(", 
     RowBox[{"3", "+", 
      RowBox[{
       SqrtBox["3"], " ", "b", " ", "\[Pi]", " ", "x"}]}], ")"}]}], "+", 
   RowBox[{
    SqrtBox[
     RowBox[{"3", " ", "\[Pi]"}]], " ", 
    SqrtBox[
     RowBox[{
      RowBox[{
       RowBox[{"-", "12"}], " ", 
       SqrtBox["3"], " ", "b", " ", "x"}], "+", 
      RowBox[{"\[Pi]", " ", 
       RowBox[{"(", 
        RowBox[{"3", "+", 
         RowBox[{"2", " ", 
          SqrtBox["3"], " ", "b", " ", "\[Pi]", " ", "x"}], "+", 
         RowBox[{
          SuperscriptBox["b", "2"], " ", 
          RowBox[{"(", 
           RowBox[{
            RowBox[{"-", "12"}], "+", 
            SuperscriptBox["\[Pi]", "2"]}], ")"}], " ", 
          SuperscriptBox["x", "2"]}]}], ")"}]}]}]]}]}]]], "Output",
 CellChangeTimes->{3.38866298634375*^9, 3.388663069359375*^9, 
  3.388663141796875*^9, 3.388663228609375*^9, 3.407484863375*^9, 
  3.407485495203125*^9, 3.439910650734375*^9}],

Cell[BoxData[
 RowBox[{
  FractionBox[
   RowBox[{
    SqrtBox["3"], " ", "x"}], "\[Pi]"], "+", 
  FractionBox[
   RowBox[{"3", " ", "b", " ", 
    SuperscriptBox["x", "2"]}], 
   SuperscriptBox["\[Pi]", "2"]]}]], "Output",
 CellChangeTimes->{3.38866298634375*^9, 3.388663069359375*^9, 
  3.388663141796875*^9, 3.388663228609375*^9, 3.407484863375*^9, 
  3.407485495203125*^9, 3.439910650875*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {RGBColor[1, 0, 0], LineBox[CompressedData["
1:eJwVkms4lAkARpEtXdzb1VSmVUxLJBpjlbyk3GkjyZN7ulCtWBVtinZIIrTz
RC7JLVkxiWkrirCKiDGaoYZhTObD8FVrUtpl2x/nOb/Ov6MXEuF5UElBQcH9
K//b5SDBbRzztnlupupXbOhko2y/8oSIQofBPcFhrhENmwuy82ooO8CUbDy5
awMDrp2Ce/mUPTjm/2wZy9gBPqEzj69SQuHhe/yMneleDJhetE+jRGMV/+/o
Kvoh3Px2V1IKhQlrrRuxmTiFBJWa9ZcoLDze8HSJn3cS9vkWst+9L0FmZfGA
+uksLJ1P8HOksmHtrDx7MaMQBruniObLtbAOjvKd9C/D6ZdWyJ5+CDfNfum5
p5XQqdUcVFRswLI4rS++0mrcF9t3dFk3gV8W77xlKwdN4fmrfGUtKO9WX7Iu
9QHktG3qlcef4V5h7l+W7Doofa7QWGXQDuZJs2SXtQ1gzAfd96rqQGBkgUxn
81N0MOkvW426sCm7OtiipwlrMiY+l091Y1s5NaA0pwUTM6LBV5t6YEatuE13
aUXJeF272hUePHcuvGyi/xxqumWkRNILb3kDt+5LG1LDOIdu6vPxOrMjSufh
C3y0rZAcOS+AfhyNb1nUiWjzltrvmvvQFde7L/ZEF0RVLZe5eq+x/+HZTy84
3fgS3ptreOYN/EwWxEuqudhgt2juUJ0Q9Bi6cKSgB45FeoerdQYRkPuWqlDK
w2ndwLXDfiKYai/m3MnthbHMe6nAfQj0eFOmRcYraKu/U2mZHYJWKcsxMY8P
5YZc9rq8YQRvXr2i8boAR+l1/YtdxSizeFLkzu5DhXPIzScfxfjeS9Ama+qH
fGKlxui1ESRJwx91tr4Gx8txAWuHBK6Gxg7skTc475RYUDsqwbQl1+fsOyH0
Blz6ihPfok618VzV7AC2WDJZFPNR5LuO81KoIvhRuz4r8kZhy58n+SuGIG9P
E59IkOJcmpXNb+whZCsHRHvqExiT27LKtg1De5BTpNRCYJjWPu/fPoxv4hL2
WEWOYabZNYgZKoZCeQBbrDUOht1VhW65GFE57DStpnHUTKdncVkjXz9MVtH8
eQKs7iuK8eYSaLpEGV9Qk+HDn/mSuVYJhCbdvNh6GbiHBb0pYW9x9wb0JQGT
kG10WVyvNAp/NkNTd24ShIVGmX3OKOITCxOEf0zBrS/SzsRWikpCN1fNjURW
MPv3MKEUyfmJQYMeJMJjwnbRRFKEek4aVO0m8SkugjYyLMXq+vq77j4kbLwD
af5SKVLS97emhpCIpr887vlBiiOM6++XxJK4xiFWblchsJa53GnRLRIlTRpT
FnQCcz+eVRXcJpFz10n7I4NA/+RIz60KEk5pHV73rQhk+NT4O1ST4Go9oFmC
wLyR5y+J9SQyPajLt7gQEHLTbyzgkWhMXq+6M4jAg6SZA7xXX/sQqnThAQKs
rYGGxX0kko2Chc8PEnAt3Vi7fZCE+cX2Fe5HCTyK6WxLGCPBiDTu33uKwDUT
i/TdMhKJjRfSKLEEIsV5e/RIEmlz+v7CXwn84HZM1DhNos3B3SMkgYCyYm9J
xgyJf37iRRgwCQxxtoYHzZIYcy69QyQRqA8vNt30L4kMs2blO5cIZK1ZKp+f
JzG1yCQmIpXAf/Gfbc8=
     "]]}, 
   {RGBColor[0, 0, 1], LineBox[CompressedData["
1:eJwVzHk41HkAx3HjWGRT0jUVGqFBqNYqOj5byhVZR6RcIVq7bZJV2o6dMjpE
tkSHcS+V1Sj0CEUaQhj3OGcY0y9lzFe7Ie2T2frj/bz+ezOCjrgfVFRQUHD5
0ledDo62Vb/ds7Vh3VzfHGOHrcq2yyJEdEv002WBbSZG+C7jZloxfQcsZi77
uJpaYVez4BGH7glPXfv062vs4B0y/fQaPQRsQ1E1LLwwaHHBNoEeBbeILdfy
LUORucg17jI9FsOCC3OiEQ2WWvHqS/RknJGHGxPPOOz1yeJOvM/F9s9K1qa/
pUJDzvK11+XC8ZvmwJiELBi6yUZfxJeApp1W9tkrH8dbrHHzwxOES95/LC8u
xJISLSGNVoVsq+KJhv6HeCy2beJvrgFr/9Qol1mKmnDOch8pD7NzFnknRJRh
0mjLvMLDL7F6ZrZc6VIFFGcK5i83bIR/zm7+1PgzWMkDH3s8aEI8o36ZktZz
NMVattSZ8MHxNeCr5ddAL2ls5p6sFfyqPzNUQngYmxYJu9a242KGvNtDtQ65
7yoaNRM78G+NquqG9pfQ1MknEkkngjnlt+beasCVn0pDMw26IaX1rdDY+QpT
PxRIDp0V4HR2H8+N3oyo9bySxS96oKYa5sdsb4HoAS++jdEHs/SYKZlPK/4L
77xjfLIfjSWnbA45tMF0m+psaMUA2NvHPXOZ7bDPZoQ9XCLE3y1Lk/J0OnBc
J0B/2FeEKjd1BlO5E2ukezQELkN4abGNVL/phPa8CTXepyHU2trxI4e6oFx1
h7sqbRhZjMNy88Zu/GxZ0au+SwzrsLBmxwYBChyDMp9NiaE0G7D3fG0PJseW
zadSRjB4xI7m9agXpR72Ssk7JCg+eWBjb0kfzjqwM0ooCbz4jvMq8/vBGHTq
yWG/RtLzA2ecUgdgsyE2mb6eQna735aetEH46vJnaB0UaKt5UwHXhJhsTBBH
sN4gourcPWa6CDeV/aPcDUZx0milQr/zELSFpdmKvFGkWSuuHKwfgspplqf1
0bforA12NfpxGAr3/LniBe/gzPlgoc8fRuRtbsKCmncQ3B3ayPYR46npRTWt
X8cg3JdSaCIWQ8spcs05TSmOED+pR9AIBsxaO2IqpWi1UFfRkI2gKB0GEv9x
aK4qWmh+VAI/rpWWzuw4+g+cEQb/I8Ef7CzWwH0Z9LibZp3PvkbhqM4dTWcC
PUdBf7QShYscdqBwN0HbsaYTdSoUQtzHDR+4ETDLivQXq1FYUVlZ5OJN0Lhz
J6f0WwqXr+6vuxJEkMK6YTa5iMIhq1vv58QQ1OvHl0czKejHLnRQzSOoCkg+
f96FwuzGU3MFdwki/PlPO10p9I6PtOcVEFj5EWUjdwpJ3sV+dg8JNHQGn9R7
UZCbuB9jVxJs/Rj6TDOQwkDb1XSlji+/Cdk6TiSFsrjp4I4ugqXdnjYkikLy
pgDjnB6CuKbbrtuOU9j1l3nJdiFBUnFt5uvfKZSfaG5gvSWwESUOmrMppJh9
f9VNSpD4XOrMukDhqDjNk0EIunYwGzouUWA6/yKq/kCQN201eiKRgjKtMzdp
mkD7inr8qyQKQ6WbwgM/EWyOL7PRvU6hMjzHYu1ngup82+mIGxRS9TQm5XKC
fcn3q1+kUvgfVbBwew==
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesOrigin->{0, 0},
  PlotRange->{{0, 1}, {0., 0.7671429965196219}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.38866298634375*^9, 3.388663069359375*^9, 
  3.388663141796875*^9, 3.388663228609375*^9, 3.407484863375*^9, 
  3.407485495203125*^9, 3.43991065103125*^9}]
}, Open  ]],

Cell[TextData[{
 "\nUsing \[Alpha] = ",
 Cell[BoxData[
  FormBox[
   StyleBox[
    RowBox[{"\[Beta]", " ", 
     RowBox[{"(", 
      RowBox[{"1", "-", 
       SuperscriptBox["x", "2"]}], ")"}]}],
    FontSlant->"Italic"], TraditionalForm]]],
 "~ \[Beta] for x << 1, yields the desired forumla. Note that the \
approximation is valid ",
 StyleBox["only",
  FontSlant->"Italic"],
 " for x \[LessLess] 1 typically useful in the range x = [0, 0.2]. The blue \
curve shows the deviation from square root solution above, and it follows \
closely up to values ~ 0.4.\n\nThe second part of Peirls formula  (x ~1) is \
derived as follows. Starting from Eq. 1', in this approximation, the sphere \
is small compared to 1/\[Alpha] and the exponent is thus negliable, leading \
to\n\n",
 Cell[BoxData[
  RowBox[{
   StyleBox[
    RowBox[{"N", 
     RowBox[{"(", 
      RowBox[{"x", ",", "y", ",", "z"}], ")"}]}],
    FontSlant->"Plain"], 
   StyleBox["=",
    FontWeight->"Bold",
    FontSlant->"Italic"], 
   StyleBox[
    RowBox[{
     StyleBox[
      FractionBox["\[Beta]", 
       RowBox[{"4", "\[Pi]"}]],
      FontWeight->"Bold"], 
     RowBox[{"\[Integral]", 
      RowBox[{"N", 
       RowBox[{
        StyleBox["(",
         FontWeight->"Plain"], 
        StyleBox[
         RowBox[{
          RowBox[{"x", "'"}], ",", 
          RowBox[{"y", "'"}], ",", 
          RowBox[{"z", "'"}]}],
         FontWeight->"Plain"], 
        StyleBox[")",
         FontWeight->"Bold"]}], 
       StyleBox[
        FractionBox["1", 
         SuperscriptBox["r", "2"]],
        FontWeight->"Plain"], 
       StyleBox[
        RowBox[{"\[DifferentialD]", 
         RowBox[{"V", "'"}]}],
        FontWeight->"Plain"]}]}]}],
    FontSlant->"Plain"]}]],
  FontFamily->"Times New Roman",
  FontSlant->"Italic"],
 "\t\t\t\t\t\t\t\t\t\t\t\t\t(11)\n\nAssuming that there is no anglar \
dependence in N, using the cosine theorem for ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["r", "2"], TraditionalForm]]],
 "\[Rule] ",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["r", "2"], TraditionalForm]]],
 "+",
 Cell[BoxData[
  FormBox[
   SuperscriptBox[
    RowBox[{"r", "'"}], "2"], TraditionalForm]]],
 "- 2 r r' cos(\[Theta]) yields \n"
}], "Text"],

Cell[TextData[{
 "N(r) = ",
 Cell[BoxData[
  FormBox[
   RowBox[{
    FractionBox["\[Beta]", 
     RowBox[{"4", "\[Pi]"}]], 
    SubsuperscriptBox["\[Integral]", "0", "r"]}], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", 
     RowBox[{"2", "\[Pi]"}]], 
    RowBox[{
     SubsuperscriptBox["\[Integral]", "0", "Pi"], " ", 
     RowBox[{
      RowBox[{"n", "(", 
       RowBox[{"r", "'"}], ")"}], " ", 
      FractionBox[
       RowBox[{
        SuperscriptBox[
         RowBox[{"r", "'"}], "2"], 
        RowBox[{"sin", "(", "\[Theta]", ")"}], 
        RowBox[{"dr", "'"}], " ", "d\[Theta]", " ", "d\[CurlyPhi]", " ", 
        "dr"}], 
       RowBox[{
        SuperscriptBox["r", "2"], "+", " ", 
        SuperscriptBox[
         RowBox[{"r", "'"}], "2"], "-", " ", 
        RowBox[{"2", " ", "r", " ", 
         RowBox[{"r", "'"}], " ", 
         RowBox[{"cos", "(", "\[Theta]", ")"}]}]}]]}]}]}], TraditionalForm]]],
 Cell[BoxData[
  FormBox["=", TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Beta]", "2"], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "r"], 
    RowBox[{
     RowBox[{"n", "(", 
      RowBox[{"r", "'"}], ")"}], " ", 
     RowBox[{"log", "(", 
      FractionBox[
       RowBox[{"r", "+", 
        RowBox[{"r", "'"}]}], 
       RowBox[{"|", 
        RowBox[{"r", "-", 
         RowBox[{"r", "'"}]}], "|"}]], ")"}], " ", 
     FractionBox[
      RowBox[{"r", "'"}], "r"]}]}], TraditionalForm]]],
 ". Putting ",
 StyleBox["x = r/R",
  FontSlant->"Italic"],
 " and",
 StyleBox[" f = nr",
  FontSlant->"Italic"]
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{"Clear", "[", 
  RowBox[{"x", ",", "r"}], "]"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"1", "/", "2"}], 
   RowBox[{"Integrate", "[", 
    RowBox[{
     RowBox[{
      RowBox[{"x", "^", "2"}], " ", 
      RowBox[{
       RowBox[{"Sin", "[", "t", "]"}], "/", 
       RowBox[{"(", 
        RowBox[{
         RowBox[{"x", "^", "2"}], " ", "+", " ", 
         RowBox[{"r", "^", "2"}], "-", 
         RowBox[{"2", " ", "r", " ", "x", " ", 
          RowBox[{"Cos", "[", "t", "]"}]}]}], ")"}]}]}], ",", 
     RowBox[{"{", 
      RowBox[{"t", ",", "0", ",", "Pi"}], "}"}], ",", 
     RowBox[{"Assumptions", " ", "\[Rule]", "  ", 
      RowBox[{
       RowBox[{
        RowBox[{"Re", "[", 
         RowBox[{
          RowBox[{"r", "/", "x"}], "+", " ", 
          RowBox[{"x", "/", "r"}]}], "]"}], " ", "\[GreaterEqual]", " ", 
        "2"}], " ", "&&", " ", 
       RowBox[{"r", " ", "\[NotEqual]", "x"}]}]}]}], "]"}]}], " ", 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{"%", "/.", " ", 
  RowBox[{
   RowBox[{
    RowBox[{"-", 
     RowBox[{"Log", "[", 
      RowBox[{
       RowBox[{"(", 
        RowBox[{"r", "-", "x"}], ")"}], "^", "2"}], "]"}]}], " ", "+", " ", 
    RowBox[{"Log", "[", 
     RowBox[{
      RowBox[{"(", 
       RowBox[{"r", "+", "x"}], ")"}], "^", "2"}], "]"}]}], " ", 
   "\[RuleDelayed]", " ", 
   RowBox[{"2", " ", 
    RowBox[{"Log", "[", 
     RowBox[{
      RowBox[{"(", 
       RowBox[{"r", "+", "x"}], ")"}], "/", 
      RowBox[{"(", 
       RowBox[{"r", "-", "x"}], ")"}]}], "]"}]}]}]}]}], "Input",
 CellChangeTimes->{3.3886631645625*^9}],

Cell[BoxData[
 FractionBox[
  RowBox[{"x", " ", 
   RowBox[{"Log", "[", 
    FractionBox[
     RowBox[{"r", "+", "x"}], 
     RowBox[{"r", "-", "x"}]], "]"}]}], 
  RowBox[{"2", " ", "r"}]]], "Output",
 CellChangeTimes->{3.38866303628125*^9, 3.3886632131875*^9, 
  3.388663284765625*^9, 3.407484912515625*^9, 3.407485543453125*^9, 
  3.439910686875*^9}]
}, Open  ]],

Cell[TextData[{
 StyleBox["\n\n",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["2", "\[Beta]R"], TraditionalForm]]],
 "f(x) = ",
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{
     SubsuperscriptBox["\[Integral]", "0", "1"], 
     RowBox[{
      RowBox[{"f", "(", 
       RowBox[{"x", "'"}], ")"}], " ", 
      RowBox[{"log", "(", 
       FractionBox[
        RowBox[{"x", "+", 
         RowBox[{"x", "'"}]}], 
        RowBox[{"|", 
         RowBox[{"x", "-", 
          RowBox[{"x", "'"}]}], "|"}]], ")"}], " ", 
      RowBox[{"dx", "'"}]}]}], ",", " ", 
    RowBox[{
     RowBox[{"f", "(", "0", ")"}], " ", "=", " ", "0.", " "}]}], 
   TraditionalForm]]],
 "This is an Eigenvalue problem, ",
 StyleBox["\[Lambda] f =",
  FontSlant->"Italic"],
 StyleBox[" L ",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox["f,",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" with",
  FontVariations->{"CompatibilityType"->0}],
 StyleBox[" \[Lambda] = ",
  FontSlant->"Italic",
  FontVariations->{"CompatibilityType"->0}],
 Cell[BoxData[
  FormBox[
   FractionBox["2", 
    RowBox[{"\[Beta]", " ", "r"}]], TraditionalForm]],
  FontSlant->"Italic"],
 ". The kernel is not bounded, but it's square is which means that one can \
get the highest Eigenvalue from the Rayleigh-Ritz quotient\n\n",
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[CapitalLambda]", "1"], TraditionalForm]]],
 "< ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    StyleBox[
     RowBox[{
      StyleBox["L",
       FontWeight->"Bold"], "f"}]], "f"], TraditionalForm]]],
 "< ",
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[CapitalLambda]", "2"], TraditionalForm]]],
 ". \n\nThe idea now is to iterate the integral using eigenfunctions to get \
an estimate of \[Lambda]. The simplest choice here is f =  x, which leads to\n\
\n",
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{
     SubsuperscriptBox["\[Integral]", "0", "1"], 
     RowBox[{
      RowBox[{"x", "'"}], " ", 
      RowBox[{"log", "(", 
       FractionBox[
        RowBox[{"x", "+", 
         RowBox[{"x", "'"}]}], 
        RowBox[{"|", 
         RowBox[{"x", "-", 
          RowBox[{"x", "'"}]}], "|"}]], ")"}], " ", 
      RowBox[{"dx", "'"}]}]}], " ", "=", " ", 
    FractionBox["1", "2"]}], TraditionalForm]]],
 "(1-",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "2"], TraditionalForm]]],
 ") ",
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{"log", "(", 
     FractionBox[
      RowBox[{"1", "+", "x"}], 
      RowBox[{"1", "-", "x"}]], ")"}], " ", "+", 
    RowBox[{"x", ".", " "}]}], TraditionalForm]]],
 "The value of \[Lambda] = ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    StyleBox[
     RowBox[{
      StyleBox["L",
       FontWeight->"Bold"], "f"}]], "f"], TraditionalForm]]],
 " for ",
 StyleBox["x=0",
  FontSlant->"Italic"],
 ", and",
 StyleBox[" x=1",
  FontSlant->"Italic"],
 " is in the range of [1,2]. An obvious improvement is \n\n",
 StyleBox["f = x - b",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "3"], TraditionalForm]],
  FontSlant->"Italic"],
 ", which gives  (see below)\n\nLf = ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "2"], TraditionalForm]]],
 "(1-",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "2"], TraditionalForm]]],
 ") l",
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{"og", "(", 
     FractionBox[
      RowBox[{"1", "+", "x"}], 
      RowBox[{"1", "-", "x"}]], ")"}], " "}], TraditionalForm]]],
 " - b { ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "4"], TraditionalForm]]],
 "(1-",
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "4"], TraditionalForm]]],
 ") l",
 Cell[BoxData[
  FormBox[
   RowBox[{
    RowBox[{"og", "(", 
     FractionBox[
      RowBox[{"1", "+", "x"}], 
      RowBox[{"1", "-", "x"}]], ")"}], " "}], TraditionalForm]]],
 "+ ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "2"], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "3"], TraditionalForm]]],
 "+ ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "6"], TraditionalForm]]],
 "x}\t\t\t\t\t\t\t\t\t\t(12)\n"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[
 RowBox[{"\[IndentingNewLine]", 
  RowBox[{
   RowBox[{
    RowBox[{"Clear", "[", 
     RowBox[{"b", ",", "x", ",", "z"}], "]"}], ";"}], "\[IndentingNewLine]", 
   RowBox[{
    RowBox[{"func", " ", "=", "  ", 
     RowBox[{"Log", "[", 
      RowBox[{
       RowBox[{"(", 
        RowBox[{"z", "+", "x"}], ")"}], "/", 
       RowBox[{"Abs", "[", 
        RowBox[{"(", 
         RowBox[{"z", "-", "x"}], ")"}], "]"}]}], "]"}]}], ";"}], 
   "\[IndentingNewLine]", 
   RowBox[{
    RowBox[{"intRaw2", " ", "=", " ", 
     RowBox[{"Integrate", "[", 
      RowBox[{
       RowBox[{
        RowBox[{"(", 
         RowBox[{"x", " ", "-", " ", 
          RowBox[{"b", " ", 
           RowBox[{"x", "^", "3"}]}]}], ")"}], " ", "func"}], ",", 
       RowBox[{"{", 
        RowBox[{"x", ",", "0", ",", "1"}], "}"}], ",", 
       RowBox[{"Assumptions", "\[Rule]", " ", 
        RowBox[{
         RowBox[{
          RowBox[{"Im", "[", "z", "]"}], " ", "\[Equal]", " ", "0"}], " ", "&&",
          " ", 
         RowBox[{"0", "<", " ", 
          RowBox[{"Re", "[", "z", "]"}], " ", "<", " ", "1"}]}]}]}], "]"}]}], 
    ";"}], "\n", 
   RowBox[{
    RowBox[{"res2", " ", "=", " ", 
     RowBox[{"Collect", "[", 
      RowBox[{"intRaw2", ",", 
       RowBox[{"{", "b", "}"}]}], "]"}]}], " ", ";"}], "\[IndentingNewLine]", 
   RowBox[{"res3", " ", "=", " ", 
    RowBox[{
     RowBox[{"Collect", "[", 
      RowBox[{"res2", ",", 
       RowBox[{"{", 
        RowBox[{"b", ",", 
         RowBox[{"Log", "[", 
          RowBox[{
           RowBox[{"(", 
            RowBox[{"1", "+", "z"}], ")"}], "/", 
           RowBox[{"(", 
            RowBox[{"1", "-", "z"}], ")"}]}], "]"}]}], "}"}]}], "]"}], " ", "//",
      " ", "FullSimplify"}]}], "\[IndentingNewLine]"}]}]], "Input",
 CellChangeTimes->{{3.388663416828125*^9, 3.388663422375*^9}}],

Cell[BoxData[
 RowBox[{
  RowBox[{
   RowBox[{"-", 
    FractionBox["1", "6"]}], " ", "z", " ", 
   RowBox[{"(", 
    RowBox[{
     RowBox[{"-", "6"}], "+", "b", "+", 
     RowBox[{"3", " ", "b", " ", 
      SuperscriptBox["z", "2"]}]}], ")"}]}], "+", 
  RowBox[{
   FractionBox["1", "4"], " ", 
   RowBox[{"(", 
    RowBox[{
     RowBox[{"-", "2"}], "+", "b", "+", 
     RowBox[{"2", " ", 
      SuperscriptBox["z", "2"]}], "-", 
     RowBox[{"b", " ", 
      SuperscriptBox["z", "4"]}]}], ")"}], " ", 
   RowBox[{"Log", "[", 
    RowBox[{
     RowBox[{"-", "1"}], "+", 
     FractionBox["2", 
      RowBox[{"1", "+", "z"}]]}], "]"}]}]}]], "Output",
 CellChangeTimes->{3.388663065375*^9, 3.38866331421875*^9, 
  3.38866345178125*^9, 3.407484946421875*^9, 3.407485573234375*^9, 
  3.439910740546875*^9}]
}, Open  ]],

Cell[TextData[{
 StyleBox[" ",
  FontSlant->"Plain"],
 "[",
 StyleBox["Incidentally, there is an error in Peierls' paper: The polynomial ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", "4"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[" (1- ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "2"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox[") in the second term should be ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["1", "4"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["(1-",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   SuperscriptBox["x", "4"], TraditionalForm]],
  FontSlant->"Italic"],
 StyleBox["). I verified the calculations in the paper with mathematica using \
the correct formula, and the numerical results are identical to those in the \
paper so it must have been a typesetting problem",
  FontSlant->"Italic"],
 "].  ",
 StyleBox["The constant b is determined such that ",
  FontSlant->"Plain"],
 "Lf/f",
 StyleBox[" assumes the same value for x=0 and x=1. Using the fact that ",
  FontSlant->"Plain"],
 Cell[BoxData[
  FormBox[
   UnderscriptBox["lim", 
    RowBox[{"x", "\[Rule]", "1"}]], TraditionalForm]]],
 StyleBox[" ",
  FontSlant->"Plain"],
 "(1-x) log(1-x) = 0 ",
 StyleBox["gives the value for b = 0.633975 (see below). The first plot shows \
the Eigenfunction and the second plot shows the ratio ",
  FontSlant->"Plain"],
 Cell[BoxData[
  FormBox[
   FractionBox[
    StyleBox[
     RowBox[{
      StyleBox["L",
       FontWeight->"Bold"], "f"}]], "f"], TraditionalForm]]],
 StyleBox["for this function in the same interval. Peierls probably tabulated \
this function around x ~ 0 to get a robust estimate of \[Lambda]. It is \
interesting to note that results by ",
  FontSlant->"Plain"],
 StyleBox["Mathematica",
  FontSlant->"Italic"],
 StyleBox[" tally almost exactly with the numerical values in the paper, \
which gives some credit to enormous work that scientists in Peierls days had \
to put into evaluating integrals and functions. ",
  FontSlant->"Plain"],
 "\n"
}], "Text",
 FontSlant->"Italic"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{
   RowBox[{"Clear", "[", "b", "]"}], ";"}], " "}], "\[IndentingNewLine]", 
 RowBox[{"ser1", " ", "=", 
  RowBox[{
   RowBox[{"Series", "[", 
    RowBox[{
     RowBox[{"res3", "/", 
      RowBox[{"(", 
       RowBox[{"z", "-", 
        RowBox[{"b", " ", 
         RowBox[{"z", "^", "3"}]}]}], ")"}]}], ",", 
     RowBox[{"{", 
      RowBox[{"z", ",", "0", ",", "0"}], "}"}]}], "]"}], "//", " ", 
   "Normal"}]}], "\[IndentingNewLine]", 
 RowBox[{"ser2", " ", "=", 
  RowBox[{
   RowBox[{"Series", "[", 
    RowBox[{
     RowBox[{"res3", "/", 
      RowBox[{"(", 
       RowBox[{"z", "-", 
        RowBox[{"b", " ", 
         RowBox[{"z", "^", "3"}]}]}], ")"}]}], " ", ",", 
     RowBox[{"{", 
      RowBox[{"z", ",", "1", ",", "0"}], "}"}]}], "]"}], "//", " ", 
   "Normal"}]}], "\[IndentingNewLine]", 
 RowBox[{"sol", " ", "=", 
  RowBox[{
   RowBox[{"Solve", "[", 
    RowBox[{
     RowBox[{"ser1", "\[Equal]", "ser2"}], ",", "b"}], "]"}], " ", "//", 
   "N"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"numsol", "  ", "=", " ", 
   RowBox[{"b", " ", "/.", "sol"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{"bval", "=", " ", 
  RowBox[{"numsol", "[", 
   RowBox[{"[", "1", "]"}], "]"}]}], "\[IndentingNewLine]", 
 RowBox[{"\[Lambda]", " ", "=", " ", 
  RowBox[{
   RowBox[{"res3", "/", 
    RowBox[{"(", 
     RowBox[{"z", "-", 
      RowBox[{"b", " ", 
       RowBox[{"z", "^", "3"}]}]}], ")"}]}], " ", "/.", " ", 
   RowBox[{"{", 
    RowBox[{
     RowBox[{"b", "\[Rule]", " ", "bval"}], ",", 
     RowBox[{"z", "\[Rule]", " ", "0.0001"}]}], 
    "}"}]}]}], "\[IndentingNewLine]", 
 RowBox[{"EigenFunction", " ", "=", " ", 
  RowBox[{"Plot", "[", 
   RowBox[{
    RowBox[{"res3", " ", "/.", 
     RowBox[{"b", "\[Rule]", " ", "bval"}]}], ",", 
    RowBox[{"{", 
     RowBox[{"z", ",", "0", ",", "1"}], "}"}]}], 
   "]"}]}], "\[IndentingNewLine]", 
 RowBox[{"RayleighRatio", " ", "=", " ", 
  RowBox[{"Plot", "[", " ", 
   RowBox[{
    RowBox[{
     RowBox[{"res3", "/", 
      RowBox[{"(", 
       RowBox[{"z", " ", "-", " ", 
        RowBox[{"b", " ", 
         RowBox[{"z", "^", "3"}]}]}], ")"}]}], " ", "/.", " ", 
     RowBox[{"b", "\[Rule]", " ", "bval"}]}], ",", 
    RowBox[{"{", 
     RowBox[{"z", ",", "0", ",", "1"}], "}"}]}], "]"}]}]}], "Input"],

Cell[BoxData[
 RowBox[{"2", "-", 
  FractionBox[
   RowBox[{"2", " ", "b"}], "3"]}]], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.439910741640625*^9}],

Cell[BoxData[
 FractionBox[
  RowBox[{
   RowBox[{"-", "3"}], "+", 
   RowBox[{"2", " ", "b"}]}], 
  RowBox[{"3", " ", 
   RowBox[{"(", 
    RowBox[{
     RowBox[{"-", "1"}], "+", "b"}], ")"}]}]]], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.439910741703125*^9}],

Cell[BoxData[
 RowBox[{"{", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"b", "\[Rule]", "0.6339745962155614`"}], "}"}], ",", 
   RowBox[{"{", 
    RowBox[{"b", "\[Rule]", "2.3660254037844384`"}], "}"}]}], "}"}]], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.43991074178125*^9}],

Cell[BoxData["0.6339745962155614`"], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.439910741859375*^9}],

Cell[BoxData["1.5773502682965752`"], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.4399107419375*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {Hue[0.67, 0.6, 0.6], LineBox[CompressedData["
1:eJwV13k8VVsbB3CkSIYMlalBRYSovArp140Gc0qlMiVXkiQRaXB1KUXmSpkS
RanILFJJUZGZBqEzrXOO4aCiCe+6f+3P9+xz9tlrr2f9nrXV3I9u/1tESEho
qbCQ0H9Hy7+5Lc94O9crPDr9YdPW/etFzZT9epUMMEtsFaurGFidkZxapGQO
zXV++xKK7WDV2FWYpuSAWSHrdXyLXbHbY/xJgpIHrGoOHtpVfBSf9S6YXVYK
gEKF/KfVxaG4Ocfu/CWlcCi/SeurK4pFmHjRsotKSRj6K0h5WWEGHPdk5g+P
ZKNpWfrmkdACzJoKc9qyIB/lQnb8+pByqNsPcV9EFaM9QsT3ndozBL0zQvK3
CgRwJkzNd9RiXrFsj7DwU7TPlrrubVCPUoZZQ9O6Gji91tiwy7kBNd5pKnsG
anF2Z5F+4apmfNcwlXlwpA58o7FZelKtEPmZN1tF/Q3CN7qdtxJrh+GUW+mO
hw3wkGkYlejqQEO4wbtXy5uwND5aqyipCwvj+n/eHWqG35r48p1HPqB/vLen
Q78Vhc9WidupfEI2v/KNdEwbbugNhgX1dUN6fo6AxWqHokbDfvfYHkQfKvG8
ubQTjEbTomkSfRjbkMfyCu3C1RLjGeodfQhYVVs898V7RK5Sc1kd9QW9D2uj
WtQ+4ojNF6muLQz89m5P0Qr5hN3undpknAHtv8QmPSu74cdXvZCZyMSWW2oH
H83rgaGspkLjKhaC5rsu/uLUi+bSoOyJJyzoDOyc1WXTBzf9D0G2pmzIywyL
1/7qg0TrWY1PVWyIPk3JX5L6BQcUtmTlanFw2KDyw0wrBnRapgqVkjjIs3C/
WT3GgFqb89+5wxx871eezbnKRIHZ9eBZ6wiS0khd6nUmds0zf2dsSrDarjh0
RyoTwwU7HA+tJzhWaCN4lslE+HP7mbUbCAaD/3mXcp+er0jrPrqJgD2dE21f
w0Rt1KI1mbYEHQsfzXwyyIRld9FKbzeCkh1bpiWZs/Ct0eFbXBhBpFV5hOoW
FlJSG6VPnyNwMtMSv23Bwvnm5OOe/xKIrJ4lVWrLglItUTGKILCVa5r73pGF
9OpH7I5IAm7zruWqR1i49tJy22AsgbKtp332FRbud8QqJKbR+9vU1aaTzMKr
7/P3HUwneGa6dVfJDRaOJUjyjTMIPHWX73uVwcJ8/8a+3psERZJDHty7LBQG
6eqrZRNYvw0M1qlmIVqNtzXsHkHo1oiMYg4LcZuKxkpLCVpdvJz0eSwMy7se
9ysjWBpopXS/n4WVMVWmWuUEbzJlE7OGWSid2l6ZXEGg8DvtfMIvFh777XM8
WkWQ+6D0iJ80G8rhVlK/nxM0y3FNdAzZKB30m4xqIFis9fZHzlo2tmZGz9Fr
JAjEw5IlJmwUqccFt1Ar+wToqWxgg100Z7FCE4HHS6ElEpZsuN2JFo5vIfgR
pCTBdWaj92Al26uTQO2z5fusCDakl/b0xvUSXKoxyd8ZSf/fbqe5Vh/BaI7O
ebEoNsqkc1ufUdf4SxscjmMjcGWyz+AXen3x1thVKWw4xd5LXMei9796z5bn
BWwMl3elF3IJVlw6WNrziY2xqCjtiWGCq0cdL8f1sLHINC03ZIRgysHCY+MX
NhRLrGzGqZsXasvd4bAx4rRGamiU1lfJ0JEjI2wk+Dr83fKNoLgvUP3PDA6k
v056Hv1BYLwmPElpFQceZVOrxKcIOt2VDNMMODAoWTR6gPp4zMOuRWs48A8e
Ha2mvs9+r6y5joOde30r/IW4WJCkm2m4mZ5nrnRsEeZCZKTzwY69HDgYmDR7
i3KRrupj1+XEgZpaetZjapOtwiN7Xek6kXpTNnM6FwEZ2gbuHhyo/th46zY1
x+afx36+HPSdXXK2fQYXb+9p1V0+x8G6+FFz6ZlceHVWe8lGcGCzUWXBdmpR
EQeJKxc4YGmpr7lCvX7PWZvUaA7WZN/SpDOFArG2tntXOVh9+KGY4iwukjxO
99Xd42Ca+567I5Jc6MfJnrN6wMH6x69Paklx0Vh5Z0lTPgf1G1ZfdKOeId/i
2VnMgW59mUcDdfBz9UFWNQfnVwxLpkhz4bSg6adwGwe1rINzFGZzIa3kX1LZ
wcGZs9n+ZtTP5eccC3zPwTKl8f5j1Boz93G5nzmIX5K0qoFa8JXT2cSl8yF9
R+qELBe3hi4mXOzngLkmyiKd2oGnY2s2xIErb+TmS+ryHv+XpV85CI79mSkr
x8W515NFaZMcpG7n8dOpDWozj+4WJujSaRV6Rk2qzbVlRQn+fZGzso/auvjS
rfCZtI6z3zHny3Mx+VDXZb0kzYm8Xx4m1AV3m5V+SBOYs+uFHKnnZMyNP6xA
ILlHODKGuv56hbX6PILgvzYfz6UOSXIS71UiOPqYEfycuufirdDtC+n3c7ta
h6njwjeZSC4mWPtEYpmYAhdmodyxl0sJxrIsr6hS5was8DVaTvBK7s5Tc+q9
R1u0vuoQuOukntpNLeUdwL6vRxBj7GB/iPqZx7xMz1UEqmGVG0Ko/V0fOy36
H0H7X08tLlEv3eus+HENgXeVtc916k4HofZEY4J/Vtjn5lCbWG62EqM528s6
7vOcesicN+P5RoK46z7fGqhvIromhObuOY34K13U2431zhpspTlZ/c7uC7Xo
/1qNhizpuomXW8KnLtUL/J5jQ2DWaCEzSn1oueKj/dsIvqftV/hJrape6aOy
gz5/nS2rpqjfLXTR7NhJkJVBvETncPGPsjArxpHguop+uTj16jnZGVv3ESyo
kVokSc2R2bJPxIXgYOGRm9LUyRL8uVU098VnmhrNpracfrk18ADBZK8X/z9P
TOnF6HkSpNr2FP3ngl+tFjwvOp5Dl67IULt/D5yedZjgysa9sVLUCsOKz518
Cb5xjTMlqOv4lafnHiMY8VJ6O4M6hO2ytvk4wfQGhqQItU6f8LeLJ2gOK4Z7
/qHj6fmYnW92kmBgO7/rO3Vcx5bDE6cIFp0S2T9EvbGZr1F2luZubIUwh/r7
m8sMP9qnbl78U9FN3T2rw+sP7UsPv9w2baOusVYVXDhPICr9z1g99d3LBwLl
L9I+c+kU68l/139373d6FM3BlkixQuoTMqNhy2MIXu5Jdr1N7bzNSLw0jkAi
+cbANerlrXUKjVcI/p5Ye/u/epktL5PimEywV7i02Zt6fMcuNdYNAsVAor+X
uraDpfub9rEZFlqxhtQuH4S3aNE+JeJ4/Dmb1rO5skVj8X2CqsyeuCZq7X1x
OzbkU9uGJZVT/+he4La7mI533oq1kdQJfSYnI6oJgtKOCP5bb8Fq/07Npn3J
0L7bRZja1f1NROoLWl/ZoiIsuj51WI4JRfUEyZJmzDvUr0hg3pdW2gf9qmSW
UN/XfLLySAdBfrWpihB14iHR8h9ddLxqv/Z003xw60+olflMcDGV5ZdA/Wso
/7MpIWBJvhD5TvOmT2/8wGsenf981ztvqev81vMdBgikRMpPZ1InjTaMHaZ9
RO/ZvueW1HpjPJkbvwkS5+TnJspwMXfNyqvqk/R51Xr1HaD+ExSs+oj2AYV2
Z8vV1G9+ztCqpzlemVNR3ETz0WNiycYx+nn4soAvkzQ/k0VdArYv5WLEuC7a
kObv5IF5/2vR4EJxzjz53zSfPV40f7fT4mJsxp3Wamr9c2ZBNivoOj9hNmFG
XT+lGbJ1LT0efztpTvN9/NfX0PXW9L57g2+ridM83vtgwxNbWrcageUfxGg9
VXgKr7Pn4nJi2bR46pjgD+eMdnHBEmreOEH7i8ZYdYSBG/39jyv7G+h97xq5
GLU8gIs1SJXfNI2Lqm1m1ndP0FzdlpA1LMLF4oI/kpon6bGy2TuFesj3aIz6
WS68j3HSh2j/ixhwiF8UyYU6KUyJos+lmCy8NjeVC82/7n7ImCCQ7ym5JVJL
59Hj2jf9cZprkmfOlL2idepv+LtijM6jibmjz2v6XK78Jb+RWiS5VarzHc29
ZZuD7L7T/YS9IPjuB9oHjt+XdP1K+/1LTbttAi52d/fvsRcQRD9I+Z2mzEPu
60FmGIfmRLd75/b5PJweuzU4wCYIl1j+SGwRDz/qp6bvpj51sNzzmDoP/5Y7
uGvS/Ye3WkeLuT4Pe53Iq6d0f7L1inRu/yYe9MrctJ930xw4E+ZgdIyH0i7L
Kilal8UjAz+fHufhcPBFb0+6P/LwdMzYfIKHibifZtXNBLV2evwdp3hwL5sf
fJjup8IXfw71jeBBfxTMSrofE61fm5d1nQeXzXJFBnV0vPLDwtI1PCgyduZE
VBIUXtiXk1TLg8Hd5sQXj2lf+PPKWqWOhxq/3lwh6hp22jXNBh4mQ/ZbB9P9
4blyK12zTh5ubIs94VxCIOyS4xjMp+PlK5tPPSQQuuuSz5Djo+x3ydKbmbTv
DHT5Zszhoznjj08B3a9G69uvcFLkY50Lv/gpzYF55WYPOubzwXLbYPaJ7n9X
1Gnm1WvyYUkCds2gueHEHr3zcD0ftplXTHQSCCrULqSHePPBjTZdWEpzzv9G
/mW5Gj4qBt+0PXYmKPtz40JNLR/tVXLxsU4Ef1zOh/nX8fE0hGt7gOb+hSXO
ga0NfCQV5LwV20Nz/L6Ec0IXH8GL7V+aO9Ccq/5bR26Qj7zam9diaB9SZKq8
lVXsR4m/aaabIcET7UhxWd9+ZF77yP0mQXOmfV5jgl8/eu4WTs+kfb/gdE68
wvF+MIJfGFuLE9xqfKWsGNyPsTPcgYzp9P3h6HTdhef6IX4+22edEMH2onP2
ulf7YVZlV2zxnQNifOaGRXU/FPLLtIt7OJC19Nc5Jz0AY48PC6ofctCt29x2
smoAiVsdcgONOChIx1KWyyAOFaR0+zxmwznfUHb+5CCUg6+6hWiz8U9EZlj3
vSH4ib5SLg1iQSpIzHL/gyGIaF52/hzAwg2vI3Kc/CHMH2h8LuLPQpGVUZag
eAjGDSLMLT4scGRbakWeDuHVl0OrK91YsEoXEtdqG8KKom1GZvS9Z26Za0zg
7yFU/6ia9lSRhazcl7t+TQzBRF61IGsOC/o3tBeGCgmgFVosfF6OBYsz4/mR
0wWI0gwNNZNk4bRZbEuKjAB7A5KtMqeYYDRVK9QsEYAU5Qww2Ew84M5PkbYW
YNf8wR3WBUxEpkW49dgK0DnKrzlG3+M8tg+qP7QXQDxLTDoxlwnVqqoCm90C
vFl8crz+JhOXYve9inYXYN2h9hVC8Ux4GV4fkTgpALcsOe+SHxPm/UJlH08J
MHOtw0cPHyYW3vQ6fe+sAEfUedFGXkx0zlwrZhlOf3/LenqrKxObPneqXowV
IFpRN6fSlonF4Qpbxe4IEL7Os8ZWm4nJtaelunLp9fQLvJ9qMPFhkNl6J08A
k9a8xmWLmYjbXeS8+ZEATfUqt3sUmZhavv14RJUAhqwXMZ+mM/Gxt2LtzqcC
ZBjNVRYWZqI0SW1yaY0AU9EqZ1X+MOA7ORxZWydA4A6bCp1RBrpbYtOntQkQ
Fpu1JP8zA+Xnxw+0dQjwtans7aH3DCSZuGplvRfAZlTGQr6NAavbK4o39giQ
fi8mRbeeAY29V0/KfRHA80OfQ0oNAyIyk+sZTAFSI71jRqoYeBzc+DqMJ0Bo
TJyfzSMGrur+L9Z+QIDEF68jrPMYOMZIdVATCGje+Ptq3GbA+pqo8siIAAsr
D35sT2dA09qn99k3AZYtzWPbJjMgKtyeHTcuANxuKF6IZ6CvxMTb7ZcAr7V2
xnleYqDKO0tPf0KAkxx92SfnGLi2cNb3qSkBSv/8jM05ycD/AZ3UP2g=
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesOrigin->{0, 0},
  PlotRange->{{0, 1}, {0., 0.752571542579062}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.43991074203125*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {Hue[0.67, 0.6, 0.6], LineBox[CompressedData["
1:eJwV2Xk0VV8UB3BDIoVQmVOmEPITTWJXSqLIPCRTiiZCQkURkjmhFEWSkMgc
yRwNCKESyhtN75IGN0O//f6yPutZ795zzvfss8966928LY5zcXBwDHBycLD/
Gh9ndtWPWuuzVqZk9ewkYYmB5LlhCW0QK4jMZHvzgzvppRJ7QT21MY1tqrRC
aICSKQh1XYpn26S9vyRDwgrsHUtC2FaIlcn5J2sPfPonTrE9b3zstquSE0Sf
v23O9lj4x2bXVFdYzWGozbat+5/aJAl3mKL4iLDdbCBZ2ZR+AvjlpCa60Zry
esU/ZU9CY83uBrbTuVzyFLNPg0sJ7SbbfCNhD22UvGAtP68T29XfutubrM5B
b3OOItvnDrsbO6b6gAytgdmFHtx03SBO4jzYnfB1ZTv5gXi94hF/qGtfKsa2
sVD+zlfpF0ApT7ztA7qC9V6HkA2C8o4yabbPODmVXHe9COProL4TLddBaKzL
vgRr/ho6sx1fKKJsrhQCJ6w6EzrQ+2QePRr1uAIS38zWsx3y+u7KBqur4MRn
+qwdrb1FLc4+NRREDZuq36MzV5tFRkuEw17NpIy36Evtt5ttqsOh0JpjOds2
Ed+45I9EwEPTCb836BW/fK7UpEfC8Be9bW3owI/JgROyN8DTpC26BW0ZO1hR
VX8DeEjFj81ojb1Kv8Jdo2FfPacE29TSSh+Z7Bhwa96c0og2v/XllKlSPEi1
hNnWodVM5PMkW+OB71j35Vfs+eY+w6B7JEDNucyMWrTzi8XndVaJ8HhZVVcN
erVKooZd6k2Ij+njrkKH8pVuuCGRDM+Ofk0oQqvuvxlVci4ZHK29rZ6huyO8
Rwdak8HrWqVoITs/3GoFGhdSQK/HIywf/XrhkcbH7lRwDzNel8Me70yqjmxs
GoQabnRKQ5f/559qREmDc6VSHbfRTucs//jsuAsDLxu3paKLJoVeNDPvQia/
zmIS2oIZtfPUvnQw9ig8EIu+M3jRoHzxPtSPL6EEo3dL2z8asn4A3nmnFy+x
8+qwlYev8AHYvopedRG9diLoWKF6Jih7/tC8gPY+uUWR42smuCqdVvVCi7gX
PcnZ9hBCK9euO4q2s88qmpp+BHJTYWJb0cSAFf9DyRzY+D6/XBsd6cR33NIg
B2pX7TTTQpce85asSM4BB59759XRgl56EZe2PgZpYfMoeXRL2GeHpcG5UHK9
okcQrfVUeKk0Xz48vHVm7rsuCW80WlzaNfNBrp6Xdxjt8jywJsQ+HzYsFRb8
io6rGD73PT8fTCK8BPvQzPrCgceHCuDXpP2nN+j7vQee/3frKfxRE79ZjF7+
L9Rx/9oiUF/mdfsi+s5JWyc1wyKwN/omF4BW+KjmIuxVBD5xM3l+aP28/mMD
tUVQQVvIP4N+kZdaXqlRDNKcuYFO6EnLiZcmHcUgZeC/Zzfa+sntd34CJfAt
yfUmNzr9aNjBE+tLQGGP1c1/O0gYETnbbqdTArtMWbFzaO/gPZ16R0tA0U4i
YAZ9w2Kym7ewBKL0zMUo6Nr5PV/uHSyFMN+g6w1oRXMWsymmDE777la+hJa8
sn+HR1YZHAnu3HEBvbIwM4a/sgy8x42MfNB/+Sw1zEfKYFjzn6UH+kN9pd/Q
tnKY+S9MwxJ9UTN0cZZWDqIVKlMq6A6hVas0dlUCj7Dpk57tJAR0bIc7P18A
417UPzG0oEY57Q1fNXDLBG4XQT+K04ydk66G/TcEvVeguw5u+Oy0rxqCvP9r
50BvfLfKTzG1Gr6k2FmPbiNhuJV4XLK1BhyMLqS8QO9vyBFsv/gSZD0O9dug
xcqEhzg564Cvcpd/+FasVyMG7zt3NsKNrl8QpkNC46kMKfuJZoA6i2K9zST8
UtITKjzbCsLpqwuHNEngIgtWSim+BSm9hTZeDRK2/HOpsHz2Hv563Lwrv5GE
9+HaHa9VOyH51QtLMWUSZBPHyTzWBxjbv949QYGE8T/DQ72a3bB34O5mkfU4
vrGat4LxPbCR7yNTXAbHL5NLUKkfYWJU32KpBAmxJ8tPZCr0QUVndqz4ahJ+
7yqgel7phwua7hurVpJwXqu5bE3TJ1gjkf81bgWO/1lzTNf6LyBgOMua5yVh
7tTHeyoXB8CfFMoO5sb52s27eKLmK+xeVbDT4t8s7H+43uO52BBw5Xtu/vZ3
FgJknOW+Ow5Dr88Ork+/Z0Ftwnp5/6FvQKW+MxCYmQVRoSm+5r/fQEeiQkCD
mIUldfeK5NO/w5Fr/ckD47NwWrvm8zKTETA7vS8tlzkL9R96n1YdGgFhyuLv
e+jVZ6eueBwegdy5GMdEdF2OolKL9Qj8ubMdgtCiYom+V11GgDUUf8wEXUO6
L//jPwIdNnHq04xZWF4nqEfNHIHB6dYfe9AuR1RW3soeAbvG5bbb0OV/DCi7
H4/ArpqJJnW0s2bQjQcFI3BeJLNSHF2aRel1qBiBlFNhvAR9Fo6EV3l1vRuB
1e8J0wx0wQG3zFe/R2CgvjefA62ttl4fyBHoXhrZ/Yc2C68Evw3UzeF4bpzl
ItAfeo6KNXBQICdo9uIg+udR+/gmfgoMCtztrEbr+x6+3LaWAvfTZpYHsP//
rr5dtyEFvA0qZn9RZ8E+eOGXxQEK3Di0U4uFHnF+eavHhAKe5UXn6eifCjs6
Px6mwOI//1X9aIkiHcN+BwpwPFeor0K7Nalpf/WigMzo5vZg9K9xyZX0VAo8
N6YyBNHJGYzW9DQKNOjxvuRFbzYru2KZToEeSsw9DrRPySGiPosC7QWfL/+g
zMJk4NWOe08pUJl17E0fmsZDjzVvpIDBoXU3s9DhlSV7+VooMKAnKnsPLX/y
yvyrVgpocTGqktGu78XPqLVTIOOA+coo9GCSiQnvJwo89RM180b3yj5fVjtJ
gfGfXqK70Oe7ghv8pihw2Ff/2A606DXjINUZCuz08anTRpvTKczbsxTIkpHP
UEG3P13T5stNhQy5+yxRdMv2y5HKElQ4G2g4PTYyC+7jRvrDUlTYlLBtOR3N
nbH6d8paKlgHrdL6jt7NUXScW4EKPwYLHveja19/3zukQYVcKyO1FnS55X7u
5L1U+FZcV/wQHWVSFSG9nwoCTuYj99GOBip8OQeokBR+T/YemmvzcoEKUypI
1/jX3kKbinSu+WRHhYDhZdcj0Ov5d6W5HKGCk9UGsTD0T87nUqNHqWBxY0tp
CPredNK6v25owf1CgWjmBxtV6bNUeK84Z3saXdPW+vSRNxXaNvgYeaLj67dt
Uvelgu6EyL7jaO1iSW39ACr8+fXbzRkdmjCk5xJKhda9glpWaEnTE+aPUqgw
e2RSfhd6cl9/j9odKridTvqih67XM7Ipv0sF2+DPGbroE+qqR14/oILCw9PG
W9GlK1juzDwqRJzy3bkJHbnEmXHuKRWUf5J26miH+c6T5DMqLDQ4h21Ec0w8
9+IvowIRTp3bgD74zj9Q7RUV1Dtzo9ejZZvof8vq8X1ifguuQ/+otg3Wa6KC
aZJL5lp0Wv72MNM2KnzdljMmhT7zMH9J/1sqnIaEbEn0rrtS153bqbDTsfSk
BJp+YyHmXDd+3tomI4Z+EeolSH7E+Uh+s3wNOjZoODG0nwpyAfK8q9HOPodF
+b9Q4ZDOkNAqtNbJhpSkr/g+aZxKomgeVy1xqWEq1PmmHBBBf7LLvpv9nQpe
iVkXhdFXjCIelNExX5O7eVaiu508HTVHqbCSd7mbEFrB30Ti6TgV9k5S3gmi
A2I0+jawqGC//OMett9mCd/KnqJCV+fIawG0dNVPM9kZKgg+FbNj27ujf8W9
X1Q4turi7xXoRmr1mzWzVBC+KvyQ7VVzGZFJf6mw5yTNge0TwqEGggtU4D+0
sI7tFxvcOaL/UcG/wOXncvQK/f21PFw0uMgr1su2k5XqxdAlNJijb2hiu/iU
wNaFpTTo6L77im3u0KmZwGU0aLUJaGXb5nZP8c/lNLjF2TLI9pPCirPnBGkw
7h7DyX7e36Y01YmVNPBd3a3F9qEvlxkeojRozM/0YTtzyvkRZTUNPOoWa9n+
sdTA1VmcBqNTM2Ls8e6VUVo7IEmDpuHrV9i+vXnZgI0MDbw3Nf5ke+zAxO1u
WRo4ROYHsOdvp0unlakcDVhxRnzs+Y6/UCL8VoEGN39kP2b7e2xKx74NNMg6
1mzBXh/t7MCYBhUaXLn1fAV7/SJfHDHSU0Nv9fvI9qdOfZ4XGjTg4hPOZ6/3
Rvr6Ru3/aGBfnxTHzsMHEaau2hZ8X6pLCDs/cirvZnO30aCn72UMO1/+8Kxc
XpcGszVSuez8SZ45v0lqFw28Yjn52Hk9E2Y7kbqHBnGr7xxi57nuzo48kX00
vKcczpJm16sWDnl+Y5xfJSM/WXTFAGU44iANqoNTp9n7Y9mP1+mcZjR45K0U
Iod+tjZuDWlJA+Xnqo2K6NkACX7mURqEaD38rYY2jp9/fcyFBoe6vzA00BmP
hq8Nu9HgxYVkhiZ6T1fOQp8HDY69zBfRYe8PVa3pFh8abPQ2ENdHD+1a/Wz3
eRq8tXAPYdeL/2xnT9VeoEGpp+HMHnTftVfUsks0yLvstcyIXd8GjT9lR9Bg
xcb/mBbo6EbdIusoGqT0my+xYa93rlokbwwNEvVltezZefYV1D6dSAMF3rFG
dr1y5+tO0LpHg0mu49Vn0e8nG09QM2igTlt53Aet01Oql5pJg7tRtev90Uvv
p4yTOTRI8DTouszO42b7/Q3FNBDuV+WPQ6+UMF7rV0qD+ywew5vooMUdvxQq
aHCNsiU5hT0/b6Szo2pooGfY4sGu5+NO3xbNXtPAnL9hohitEe1RMTRAg7Ov
CwI+o1O97eISh2jgKvtzcAj9z+qA+57vNPjcC7ZUdj5kN4o8ptPAz1rgGoH2
KWedPTtNg4Wz5fVL8Xwq++avOL+UDpE8xOA2tPTrE/OFy+iwWJIRBOiIAtse
5xV0mHzCqWKItrmwPbRJmA4F4WuaLNF/li98jZamAzOpmMY+P3dsDU+W0KKD
sW8aqwDd5yaxJUObDmVDVwxK0X7xz/rXbaWD1nG1/Gr0U9onSeWddHA2XHj6
Br02WT1riyEd0hsjbjDRXNN9hZYOdCiUCnysgv3AfekzZv2OdJCfVJPWQusa
cU47ONPh/GV4tAN9/sFGbTd3OqzVeTRugqYfulp9zosOQT5Bdt7od/kqrXFh
dKg68abkBdqz75WncAQdpIO+WTShl3BZ8adcp8MauymedrS+fcih9Fg6GFSu
z/mGLubt6clPpYOiw1AbH/ZLye6Xv7Xm4/M9jz9zRDuu7SQ5e+jgUZKwVRT7
M0EJ3/KaXjpwx5aarkM3iK728f9EB7+w9GB1tNKyI0zmIB2siWk5IzQxQ+/r
ZNLBbFo8LgQd9maxNGORDjKKn8bY/eCT8xpe21UZcM1WVIWG/aKDd5fKjBoD
nGsvWPxCC5w6T3u6iQEhxw6k8GB/6utc7bhOhwG6TF8PJbSusaEJ7y4G5K6x
hpPoDlkn5V5rBsQM2LyZQf96GzdyLpQBedWQKz02C1+X93rOX2OAa1Bthwa6
8aA0cT2SAdZi4wK70Ykd+XP3YxiwTzi85zhatbt1VXsKAx4nTV0uRjt95tyv
ks8Asa9fLIyx/37N8C/43s2Afu0PjokTuP7Ktf+d7WWAR/3XxBz0rZNLqmb7
GbBZcuVQNdplPKlZaJABxUOSzTT0X1bRoB6DAT29gvU7J2dh0+9RobtzDJBc
GAsl0HeWOJ23UGCCHXdfqCf2/4vHxHS6lJhQXyAbcBXt3vThl5kKE84cWBZ5
B60ZZhBwSIMJHyLGptrQbf+ULxptY8JLN+V51SnM99+ZK/oHmRDpX/pjBm0z
fSNG9TwT5uwuPU/7MQsvDxsczLvABFbq7ukytFzx/ArlICaE35A++AHN8vKO
VwxhgvttsOTF+0nEhNXNdVFMKJu+9uQCuowhe3tNOhPkaDZBDj+xvxwqf8jV
zIRXupNLt+L9Jn5FcHDlayZ8NX1yxA7Nr7vX7swbJiR3LHYFobnudAv0deDz
66vla9E/zInAvM9M+GsqdnTfH9z/LcpmhwkmGLsWaDvOYj0tvDeXITkKMBv+
PRvvUzzBoVbbfUZhpVlB6H28b5VNT5B1fqPArKaMtKDdT9g9MLwwCvln8t0n
0c1mm8YsL43CrZcib3dykBAuN3jFK2IUdgg9kBlEL2nbVpCdNgqNw5I/5Lnw
Pio6xSnYOAp0kbPZb5eQwJHnVDQiMgbGrRNPQ/hJeD3R7/Vg9RhQix4klqBj
Nc01HMXHIGUrM4OOFqsyKOyVGYNQf28ds+UkaLQqF7Qpj4HNYJu9HN4nHWk/
Hj/Tx+/TWXbtgwAJL9Zfv3/x1BhcMr5RdFCYBN+7RXEijWPgVbFm6rU4CbUb
o/iEvcZhv/bTyEsqJAgb+6qFCU6AsVhlboguCV/VP/QEvZyAtoCuV+fNSSi+
DwpUp0loWPWgOsGDhJmLM7yurpMglV3WzEBvsc0dHzw2CQqRhzvBk4QaIaHS
fs9JuL5njEKgX18d2v3OdxLC1jH+mp0iYcAt2LkkchLGIhuPCJ0lYemGmrtX
n02Ce/24XaQvCQe4vUMWiiehIOW33iA6dljO7WLpJOSWcQhq++H73olWOV81
CQsmdae/o6WXOVR5NE7C5Pd5Hl1/ErTGZ3tN+ybh6BO7CWYACUeLtgjLLE6C
+t470crBJFQ4zAp+4GCB8Lzrr1C04NLqFde4WZD74s66AXS9ox7fKB8LDgTX
Po8NIUGOf++/MlEWzH+MTpi4QgLd3XzyoAoLwl4yczPCSNBfKTr+byPa75Hh
D/Ttmo/MEg0WGG7yjDK8RoKRiB1VXJsFBQUPJybR+XVOA1R9FiwO685tjSDB
S/LMm8tWLDjInaheeR3nq0W9dZMtCya3bmQtiSJB1odoHrFngZowh5oF+kOr
b72RMwusqgJXjaO1/IMqV51igVHtwoE10ST86YzMeXqFBfa1aU/MYkk4fMko
2zmMBRLqXIO30E+U+LNEIligLHLkaD/aITguPTCaBSmXkzc7xpHwUjX51t4U
FrRzhHQ4xpNwNSIr9Gs+C4ChZrglkQSBAF5j10IWePXp1Z9E3/U8K0IvYsF7
7c/u6ehSk+3ZRBkLVmbKr/2Hpgt3NXPVscAvR/vbi5uYR+6tcdcbWMDnq6dP
R//7mW69opkFMbIfcoSTSBD/5Elf/Qa/73bKleNok/scfCo9LNg59ekLxy0S
PiWc+FDYy4J4gVkuRbR76Ps7Wp9Y4H/VStAIHeJ+R2XnIAtOMjpeRKP5bRan
64dZwPNWwywffXu/e/W+ERY4+pmWtKGLVTVNzBgs6LCUeseZTIKedKrox1EW
VDgp+kuj3wrMD9hNsCBbJ7dXB23zz/XRIIsFf/VD6YfQlKnWM27TLPgq7pfl
jj43oq7DmGHBsqeSHBfR8z23Fk7/ZkGbzsxcHPpGC9kyNcuCLCOr+Ez0mkrn
eP85Fug3nSt+js5+0mLzd4EFG2qtjzSgNe9ulL3CQcCE5+j1TvTLmJsMbm4C
+sZGNb6iDwT/KYriIYBfYcKIge7zOhoowEfg/dyhYwrt5tK0K4mfAFdHv8pZ
NGGuskxMgICaXUbzi+jLBgld94QISP1Wk8ydQsIynV9p60QI+DR+NXApOlXp
iFvOKgIqL4g/4EXLiTeoqooR8ML5wz/250XLNsw8kyDgV3bDzSXonXOxNZul
CWjN7DDhQLdN/LhWtZYAYneWwl98nvWQ3UG99QR0SuaJ/ECPdL5a1ShPwM+T
EiuYaO8GhUFDJQJseRhL2eObL4nOeadMwG6N17Pt6KhHU2cPbyRA3ES5vxa9
OtVmS686Aasjt98uQD+8/nLRXpOAaNMh9dtojSC51iEtHF9NccJVdM2pqIRj
OgQIdywUe6CNHFm2zK0EVCU0RhxEfzxkte7sDgKuCuRxbkK7QjVzeicBGao8
4kJolua65xeAgOYDXTUTmJ+LcpFBc7sJ+PEmf7gVzbtqYvfVvQTUWnb4Z6KT
eSz4efYToPhT9NwFdCFT5p7gQQIeN1vbSKKjMiJchkwJON58bO0I5t3dYlLx
mTkBPKN/bB+jpV++LD5kS8Bo2/4mefTsOYULMg4E/LfjbfNX3D89irG6k44E
bDtq0n8THZ1w5HWsGwG/O4RaCNx/J/Y2xTgeJ8DjdJ3nHfQeUtVczZOA9DnV
ZzvRpNvfgXdnCQiIXJIdkIDjF3PLuncO8+TzTlgQXfz+zYnTfvh51kDmA6wH
nlvSpvmDCIiTuZxQivVi7zhH5ZdLBDxjcGfroGUzPS/nhxCQ/+emfgnWl75l
23iNwwnQ04l0uxtDQsmrB+8lrhNwZY9DIg863o83afQGAU+l3yaewvq1b7BP
+kYCAdLrvXTX3yBhXZL+iF0SAUUrVgf6Yr2bN3ycq5xCQJ5Pf8FLrJdlxf5a
bXcJ7N/s5ndEYj7DVxnxPsb1/v02uh/r8eK2ywL9Twjg3F3zlRlKwudJSvfj
AgLkvHobpq6SkGhbetTwOQE2rW+LRrDen1khJb+mDOfveOXAOzwPjBrCmLQK
Al7T7YfzLmO9UrXwi3hJwGapT+Z7LpLwZfjFNus6AqTCe1PnAvF8SV6/qNBI
wNtjY465eP54LU5FNbcScPnN40s950kwLrUzTX5LwAmPC+oWeH4petaLurcT
MEMZWlHvg+dvV8J97h4CLMeW/HT3IqEq8s+xnl4CJC+d9E87g3nSdVbJ/oT7
Q1eMvwrPU5McjbI9QwTsK36v+/QECUoOqUEi3wk4WWERHuyO/YfQov4IBfMs
Wiio6UZCdWD7m9BR3B9VRYTqUdzf6joJ5hMElKebqJ5yIMFnJN1qPUFArHlB
RYAtCQdvL5GcniZgPo9nyNiKBOWDZ4brfxKQ+zDfuf8w9jucHx8l/sG8fvpQ
yn+IhG/luqdc/hKQkGZ5qsUI69Wp7E2aCwQkebg4txhg/ZVd/uvfP6xnefIz
7N9H/weACIq9
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesOrigin->{0, 1.56},
  PlotRange->{{0, 1}, {1.555186981725615, 1.5988842215223664`}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.3886630655*^9, 3.388663154859375*^9, 3.3886633143125*^9, 
  3.4074849465*^9, 3.40748557340625*^9, 3.43991074265625*^9}]
}, Open  ]],

Cell[TextData[{
 "The eigenvalue at x = 0 is then \[Lambda]=",
 Cell[BoxData[
  FormBox[
   FractionBox["2", 
    RowBox[{"\[Beta]", " ", "r"}]], TraditionalForm]]],
 "= 1.57735. This approximation was done assuming  \[Alpha] = 0, and this is \
now improved using a first order Taylor expansion of Eq. 1': \n",
 StyleBox["lf = ",
  FontSlant->"Italic"],
 StyleBox["L",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox["f + ",
  FontSlant->"Italic"],
 StyleBox["M",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox[" f = ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   RowBox[{
    FractionBox["\[Beta]", 
     RowBox[{"4", "\[Pi]"}]], 
    RowBox[{"\[Integral]", 
     RowBox[{
      RowBox[{"n", "(", 
       RowBox[{"r", "'"}], ")"}], 
      FractionBox[
       RowBox[{" ", 
        RowBox[{"(", 
         RowBox[{"1", "-", 
          RowBox[{"\[Alpha]", " ", "r"}]}], ")"}]}], 
       SuperscriptBox["r", "2"]], " ", 
      SuperscriptBox[
       RowBox[{"r", "'"}], "2"], " "}]}]}], TraditionalForm]]],
 "sin(\[Theta]) d\[Theta] d\[CurlyPhi] dr = ",
 StyleBox["L",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox["f  - ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", "2"], TraditionalForm]]],
 StyleBox[" ",
  FontSlant->"Italic"],
 Cell[BoxData[
  FormBox[
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "r"], 
    RowBox[{
     SubsuperscriptBox["\[Integral]", "0", "pi"], 
     RowBox[{
      RowBox[{"sin", "(", "\[Theta]", ")"}], " ", 
      FractionBox[
       SuperscriptBox[
        RowBox[{"r", "'"}], "2"], 
       SqrtBox[
        RowBox[{
         SuperscriptBox["r", "2"], "+", 
         SuperscriptBox["r", "2"], "-", 
         RowBox[{"2", " ", "r", " ", 
          RowBox[{"r", "'"}], " ", 
          RowBox[{"cos", "(", "\[Theta]", ")"}]}]}]]]}]}]}], 
   TraditionalForm]]],
 "d\[Theta] dr . The \[Theta] integral turns out to be\n-2\[Alpha] (|r + r'| \
- | r - r'|) ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"r", "'"}], "r"], TraditionalForm]]],
 "and so \n\n",
 StyleBox[" ",
  FontSlant->"Italic"],
 StyleBox["M",
  FontWeight->"Bold",
  FontSlant->"Italic"],
 StyleBox[" f ",
  FontSlant->"Italic"],
 " = -\[Alpha]r",
 Cell[BoxData[
  FormBox[
   RowBox[{
    SubsuperscriptBox["\[Integral]", "0", "1"], 
    RowBox[{
     RowBox[{"f", "(", "x", ")"}], " ", 
     RowBox[{"{", 
      RowBox[{"|", 
       RowBox[{"x", " ", "+", " ", 
        RowBox[{"x", "'"}]}], "|", " ", 
       RowBox[{"-", " ", 
        RowBox[{"|", 
         RowBox[{"x", " ", "-", " ", 
          RowBox[{"x", "'"}]}], "|"}]}]}], "}"}], " ", 
     RowBox[{
      RowBox[{"dx", "'"}], ".", " "}]}]}], TraditionalForm]]],
 "Now we turn to first order perturbation theory which states that a small \
perturbation of the value is given by\n \n",
 StyleBox[" \[Delta] l =",
  FontSlant->"Italic"],
 " - ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"\[Integral]", 
     RowBox[{
      StyleBox[
       RowBox[{"f", 
        StyleBox["M",
         FontWeight->"Bold"], "f"}]], 
      RowBox[{"\[DifferentialD]", "x"}]}]}], 
    RowBox[{"\[Integral]", 
     RowBox[{
      SuperscriptBox["f", "2"], 
      RowBox[{"\[DifferentialD]", "x"}]}]}]], TraditionalForm]]],
 "\n\nUsing the eigenfunction in Eq 12 gives \n"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{"Mf", "=", " ", 
  RowBox[{"Integrate", "[", 
   RowBox[{
    RowBox[{
     RowBox[{"(", 
      RowBox[{"z", "-", " ", 
       RowBox[{"bval", " ", 
        RowBox[{"z", "^", "3"}]}]}], ")"}], 
     RowBox[{"(", 
      RowBox[{
       RowBox[{"Abs", "[", 
        RowBox[{"x", "+", "z"}], "]"}], "-", " ", 
       RowBox[{"Abs", "[", 
        RowBox[{"x", "-", "z"}], "]"}]}], ")"}]}], ",", 
    RowBox[{"{", 
     RowBox[{"z", ",", "0", ",", "1"}], "}"}], ",", 
    RowBox[{"Assumptions", "\[Rule]", " ", 
     RowBox[{
      RowBox[{
       RowBox[{"Im", "[", "z", "]"}], " ", "\[Equal]", " ", "0"}], " ", "&&", 
      " ", 
      RowBox[{"0", " ", "<", " ", 
       RowBox[{"Re", "[", "z", "]"}], " ", "<", " ", "1"}], " ", "&&", 
      RowBox[{
       RowBox[{"Im", "[", "x", "]"}], " ", "\[Equal]", " ", "0"}], " ", "&&", 
      " ", 
      RowBox[{"0", " ", "<", " ", 
       RowBox[{"Re", "[", "x", "]"}], " ", "<", " ", "1"}]}]}]}], 
   "]"}]}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"denom", " ", "=", " ", 
   RowBox[{"Integrate", "[", 
    RowBox[{
     RowBox[{
      RowBox[{"(", 
       RowBox[{"x", " ", "-", " ", 
        RowBox[{"bval", " ", 
         RowBox[{"x", "^", "3"}]}]}], ")"}], "Mf"}], ",", 
     RowBox[{"{", 
      RowBox[{"x", ",", "0", ",", "1"}], "}"}]}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"nom", " ", "=", " ", 
   RowBox[{"Integrate", "[", 
    RowBox[{
     RowBox[{
      RowBox[{"(", 
       RowBox[{"x", "-", " ", 
        RowBox[{"bval", " ", 
         RowBox[{"x", "^", "3"}]}]}], ")"}], "^", "2"}], ",", 
     RowBox[{"{", 
      RowBox[{"x", ",", "0", ",", "1"}], "}"}]}], "]"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{"\[Delta]l", " ", "=", " ", 
  RowBox[{"denom", " ", "/", " ", "nom"}], " "}]}], "Input"],

Cell[BoxData[
 RowBox[{"4.156376983941085`*^-10", " ", 
  RowBox[{"(", 
   RowBox[{
    RowBox[{"1.643288625`*^9", " ", "x"}], "-", 
    RowBox[{"8.0198051`*^8", " ", 
     SuperscriptBox["x", "3"]}], "+", 
    RowBox[{"1.52530581`*^8", " ", 
     SuperscriptBox["x", "5"]}]}], ")"}]}]], "Output",
 CellChangeTimes->{3.388663066625*^9, 3.388663329296875*^9, 
  3.40748494771875*^9, 3.40748557471875*^9, 3.439910744921875*^9}],

Cell[BoxData["0.7960118233770359`"], "Output",
 CellChangeTimes->{3.388663066625*^9, 3.388663329296875*^9, 
  3.40748494771875*^9, 3.40748557471875*^9, 3.439910745546875*^9}]
}, Open  ]],

Cell[TextData[{
 "This means that \n\n",
 Cell[BoxData[
  FormBox[
   FractionBox["1", 
    RowBox[{"\[Beta]", " ", "r"}]], TraditionalForm]]],
 "= ",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Lambda]", "2"], TraditionalForm]]],
 "-",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Delta]l", "2"], TraditionalForm]]],
 "= 0.788 - 0.398 \[Alpha] r . Now Peierls uses a trick: For small \[Alpha] \
use the fact that x = ",
 Cell[BoxData[
  FormBox[
   SqrtBox[
    RowBox[{"1", "-", 
     FractionBox["\[Alpha]", "\[Beta]"]}]], TraditionalForm]]],
 "\[Tilde]1 - ",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", 
    RowBox[{"2", " ", "\[Beta]"}]], TraditionalForm]]],
 " and thus \[Beta] = ",
 Cell[BoxData[
  FormBox[
   FractionBox["\[Alpha]", 
    RowBox[{"2", 
     RowBox[{"(", 
      RowBox[{"1", "-", "x"}], ")"}]}]], TraditionalForm]]],
 ". Now use the value for r when \[Alpha]=0 on the rhs: r \[Tilde] ",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "0.788"], TraditionalForm]]],
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Beta]"], TraditionalForm]]],
 "= ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"2", 
     RowBox[{"(", 
      RowBox[{"1", "-", "x"}], ")"}]}], 
    RowBox[{"0.788", " ", "\[Alpha]"}]], TraditionalForm]]],
 " and so\n\n",
 Cell[BoxData[
  FormBox[
   FractionBox["1", "\[Beta]r"], TraditionalForm]]],
 "= 0.788 - ",
 Cell[BoxData[
  FormBox[
   FractionBox[
    RowBox[{"0.398", " ", "*", "2", 
     RowBox[{"(", 
      RowBox[{"1", "-", "x"}], ")"}]}], "0.78"], TraditionalForm]]],
 "= 0.788 - 1.01(1-x)  "
}], "Text"],

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

Cell["\<\
Graphing this gives the following plot (fig. 1 in the original paper). Since \
none of the formulae are valid in the intermediate range, some interpolation \
must done.   Thus Peierls formula must be seen as a back-of-the-envelope \
estimate and more accurate values must be calculated numerically as in eq. 9\
\>", "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{"p1", " ", "=", " ", 
  RowBox[{"Plot", "[", 
   RowBox[{
    RowBox[{
     RowBox[{"0.552", " ", "x"}], " ", "+", " ", 
     RowBox[{"0.216", " ", 
      RowBox[{"x", "^", "2"}]}]}], ",", 
    RowBox[{"{", 
     RowBox[{"x", ",", "0", ",", "0.7"}], "}"}], ",", 
    RowBox[{"PlotStyle", "\[Rule]", " ", 
     RowBox[{"RGBColor", "[", 
      RowBox[{"1", ",", "0", ",", "0"}], "]"}]}]}], 
   "]"}]}], "\[IndentingNewLine]", 
 RowBox[{"p2", " ", "=", " ", 
  RowBox[{"Plot", "[", 
   RowBox[{
    RowBox[{"0.78", " ", "-", " ", 
     RowBox[{"1.01", " ", 
      RowBox[{"(", 
       RowBox[{"1", "-", "x"}], ")"}]}]}], ",", 
    RowBox[{"{", 
     RowBox[{"x", ",", "0.6", ",", "1"}], "}"}], ",", 
    RowBox[{"PlotStyle", "\[Rule]", " ", 
     RowBox[{"RGBColor", "[", 
      RowBox[{"0", ",", "0", ",", "1"}], "]"}]}]}], 
   "]"}]}], "\[IndentingNewLine]", 
 RowBox[{"Show", "[", 
  RowBox[{
   RowBox[{"{", 
    RowBox[{"p1", ",", "p2"}], "}"}], ",", 
   RowBox[{"AxesLabel", "\[Rule]", " ", 
    RowBox[{"{", 
     RowBox[{
     "\"\<x = \!\(\*SqrtBox[\(1 - \*FractionBox[\(\[Alpha]\), \(\[Beta]\)]\)]\
\)\>\"", ",", " ", 
      "\"\<\!\(\*FractionBox[\(1\), SubscriptBox[\(\[Beta]r\), \
\(c\)]]\)\>\""}], "}"}]}], ",", 
   RowBox[{"AspectRatio", "\[Rule]", " ", "Automatic"}]}], "]"}]}], "Input"],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {RGBColor[1, 0, 0], LineBox[CompressedData["
1:eJwVjXk4lAkAh5nW1hYbg7aabMtMLGpdFcL8yIMSKiLJMSnK8eSehGi3xJa2
VFuI1FaUoR2hcZZrlsY9yD2fY+brkKMiQtr6433e5/3rVfUOcvShSElJ2X/l
m2+6WSVl5zkw/5KfvG04Yc40qmW0GjE2wvW4UFGVuh55nM7webopwnqDWIlU
HRxb/1tdP90O9pGz18oVAasDbsFNdHe0XOc27F+1G5vZk7V8eiBMraYys2le
4Jw1maiix+C5m2TeST0IL1cc+q6SngQz7xi5ApM4DO3c4c5zzECHmuga3fMy
1Bo1NHzq7yPpJPN0R2gKDFYUEpyaXGxIYZav3XYbJZYRXobaj0FZaueflHEP
7COyRS8LivCwW3uJftUDCO7OqrgklmCQpsKv2JmHlToP1UdsKzAbk3yrRMxF
ZHPOeS/fSlDSvEYPLSlA4TXOhY321XiY5V8ldC4C7HK+3FasxdU/piuPnuOh
WvNCDaOND+/xD43M/BIkX1+VnFBYB+WP8j7aCuV4NJXZlnziObYWlsbN2j/F
T8PttGhWA27u2RbEZlUicvp1UZh5E1I+eoxa3aiCxtznakWpFlyeOzdKlFVD
wUYU+Dy0FTHFpS2GCzV48md5vrtzGwafVVXKaPEhrRQcTdgKkfJhpTwt/D/s
yFcae23SDpMlIWyrsjpstTFaXb+uA8Yui8cxUQ8BQ6/+b/lOpEpfFGibCyAz
/0DPcKoTgiYr4sqpBpiKFSxZ7S/QrCt3wIfXiDXKe3oi+V34KMySbvuxGRlM
OU2f7G6oyKTaKDu2oO7fzM/m6T3wHjvypuJAK8YM6FH3knoBjgu1alkbSrkb
uH2BfbDR/fWcR14bNkV89qew+vEPS3CYdUgIG/19xML2AexVkoFweTs8KOnz
saoieFxd7oXH7ViztiR9kkogNzyD6uvXgRgt55ktUwR+Vjw6maPcie0WhwO2
hA9igSlbMdbWiTS/Pp4xOYgNA8WbGWdeIJhm1y3tMoQvvxvRtup3YU2IRplF
yxBmdHKjU991QdmyOirAfBip8VTnzPvdEBfXLciWDWNmpN/RxL8HesV8x1CD
EZwWd31q3NgLocGZ+Kq7I+D80PQ+b6QX1k6cqEQVMcJTE6KjOH2Qm7hiZ3zp
a7daTwcc7YejUtj7vEUx6h/tcln8ZQD+AxVuaSckmHtCud4yOoAXavWZwa8k
WMzQbQnPEcFsvIZ97I0Es+82qTbkiXCv5Jk9660E/q7vLqrmixCyhze/e1IC
rruA18wTQTY2a7/OrAS92/l2WnwRLHrOrhxfSuJivCiaJETIvWQRF6BBIk58
azxUmYDiQTOXw5okfBM8UxtXE4hSN950UJtEranYU30dAZsK3T5bHRJikcW2
XjUCw2/WG2oZkngqZa1vqUtgtfXi+CtrEp2lIRprdxGIU5jjD+4goR9/JTbC
gQDZP53ebUvie++Tk617CRSEjdnWO5CgOqUyEl0JONzpz8p2+fqjqjXP+BAo
DOw+lelKoi84bN7JjwDNqGPfDTcS/kPnLbmBBF41N1ASPEm8/bSC4RdGYHda
XU8si8TwnYhyPpvAE58aLtubhIP7/SC1KALr9J4lHD9CQlfzplncKQJnFko9
fX1JXFi2X73/NIH/ATXoZ0g=
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesOrigin->{0, 0},
  PlotRange->{{0, 0.7}, {0., 0.49223998779428574`}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.38866306684375*^9, 3.38866333346875*^9, 
  3.407484947953125*^9, 3.407485574921875*^9, 3.43991074603125*^9}],

Cell[BoxData[
 GraphicsBox[{{}, {}, 
   {RGBColor[0, 0, 1], LineBox[CompressedData["
1:eJwVkms0lHkAh1EjaSpFxUm7FIqkkkslfui0ym4MySXUZG0kuqglCRGtRaU6
bKdtTS617cklVmlsZehCl01DMS6D8c478zfMvLEiWWbbD895Pj/nPKZhR3x/
0NLQ0Nj5hf+da5HnuGEDhcbtZf5heiKctduXbGdHoa03bXbUQhHi3C2f2dtT
oOMCHx4zECFk7yPfjY4UtG7OWJ5qKMKqfDrG2YmCs+YeJc9EBAHLocRjK4Uq
vna6eJ0IjLRNL2QXhetW3MpAHxG8ig0V6ccpHGXrs70vi+BUNPUb+BQKjeLi
Lb/uwHM2Jyz4MwWXgIXB5UUdOKGYY7tusxQ+J5NfFazohK+n5YzDKVKo+uXD
pys68feiailbIMVYUlYuy7YL8d/faDTQppERZkUn1nVhm9nYvfeeNHKsrj0P
9ezGXQ6/ejyLxvTyrSaKV934Vi2uTntLYyJ8YG14oBhnjY4+O7xABm57U01m
lxi9tw5J9H1kcFlz8rVhSA+i7jdNm+bLUPza7Jd8eQ+urL9ifalNhrjrup2P
D/Zi0KYu+KCxHMJ7ypZ3Y704wPPL/8yVgxPEsb7j2oeUxADJ0A05VEe+We2W
2ofGjr/WbJfJ4bPaieXV2Iey4vNZM1cSKJK/W3pXR4IGRkAsYwg84kPYP3tL
MFke6qtTRbBv/7snS/IkEHZEvPYYISijImsm30vQdq7V44PjAOp0/83wN+5H
Aq+0deLkAGa9vJOrvbcfWM1EPXwwAN3EJ/ujC/phv+nKDDI1gJa3HXN3SPvR
VFVSesJFgebH7LFT5hSyy5cERKUrcNnHhb8imsKxZaNsoxcKzGtrd0+rpPBx
3Ll11exBNIyc13EYoXDbb+BaodcgNlW7/8HdKMUnY63jiRcHIZtrcK4vQYro
wHT/ec2D8Araabm4Vor5iyOddA2GkLgj4/cGtRQc8NdF+A2hR9ix/6krDYHi
0Mo114bAUluPb/mJxhvNHBu/riGEPb3epPuChjp7kZvdUiVOJZ7gXpgtw76r
GruTuEr4t3u2ZnrLEKoXlLCyUInIV85polwZsnT1eZvlSiB8LZvXLkNXvMNb
bWsVhKrMkKXL5HjBaVxgF63C1QVqMhQqh/wO3+9lpQpmHcc59iVyFMUa3a4f
UUG2kTXaQsuRWtqvHWvLIDxrOtTCguBi5NyI7HgGN8+MFX04RLA592Y76wED
4/lzzPTKCKxj7QXmfAZtBbeUtuUEy/ye3d5Wy+CCjXvN7goCDUM6IeMhAw2v
BM9fKwkaeWbGrHoGdI78qMV9Ar+KYu7MlwwqdJ8+2lJHcPhN4YBWNwM37aSA
g0IC7t31LcvFDCbzDE1zWgh8LtfXuvcw+NO8WlHeSmDnL8lO62Ngvm0wafQ9
waTYxEaLZjAnPagkpZMgU8mL1VQxeKfpMJxHEZxqXhtsyjA4f0lY+0BKEF1Z
t9XtAwMP05j0LprA+8deg9QRBnzXkiWmhGDx1Fc1GuMMClIWupQOEej0lvNM
PjEInF8+q1lJMCFwyXSdYKDH2yEcVhF0p+8NOjPJIO3xmXDHYYI3B1RuN/5l
sMXL2GbPl68F25OtBFMMPoprxk//Q1BlNU+/b/pLf8yuet4oQTG7YFKtZhAx
pcpq+EjwHwpff+s=
     "]]}},
  AspectRatio->NCache[GoldenRatio^(-1), 0.6180339887498948],
  Axes->True,
  AxesOrigin->{0.6, 0.4},
  PlotRange->{{0.6, 1}, {0.37600000824489793`, 0.779999991755102}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.38866306684375*^9, 3.38866333346875*^9, 
  3.407484947953125*^9, 3.407485574921875*^9, 3.4399107461875*^9}],

Cell[BoxData[
 GraphicsBox[{{{}, {}, 
    {RGBColor[1, 0, 0], LineBox[CompressedData["
1:eJwVjXk4lAkAh5nW1hYbg7aabMtMLGpdFcL8yIMSKiLJMSnK8eSehGi3xJa2
VFuI1FaUoR2hcZZrlsY9yD2fY+brkKMiQtr6433e5/3rVfUOcvShSElJ2X/l
m2+6WSVl5zkw/5KfvG04Yc40qmW0GjE2wvW4UFGVuh55nM7webopwnqDWIlU
HRxb/1tdP90O9pGz18oVAasDbsFNdHe0XOc27F+1G5vZk7V8eiBMraYys2le
4Jw1maiix+C5m2TeST0IL1cc+q6SngQz7xi5ApM4DO3c4c5zzECHmuga3fMy
1Bo1NHzq7yPpJPN0R2gKDFYUEpyaXGxIYZav3XYbJZYRXobaj0FZaueflHEP
7COyRS8LivCwW3uJftUDCO7OqrgklmCQpsKv2JmHlToP1UdsKzAbk3yrRMxF
ZHPOeS/fSlDSvEYPLSlA4TXOhY321XiY5V8ldC4C7HK+3FasxdU/piuPnuOh
WvNCDaOND+/xD43M/BIkX1+VnFBYB+WP8j7aCuV4NJXZlnziObYWlsbN2j/F
T8PttGhWA27u2RbEZlUicvp1UZh5E1I+eoxa3aiCxtznakWpFlyeOzdKlFVD
wUYU+Dy0FTHFpS2GCzV48md5vrtzGwafVVXKaPEhrRQcTdgKkfJhpTwt/D/s
yFcae23SDpMlIWyrsjpstTFaXb+uA8Yui8cxUQ8BQ6/+b/lOpEpfFGibCyAz
/0DPcKoTgiYr4sqpBpiKFSxZ7S/QrCt3wIfXiDXKe3oi+V34KMySbvuxGRlM
OU2f7G6oyKTaKDu2oO7fzM/m6T3wHjvypuJAK8YM6FH3knoBjgu1alkbSrkb
uH2BfbDR/fWcR14bNkV89qew+vEPS3CYdUgIG/19xML2AexVkoFweTs8KOnz
saoieFxd7oXH7ViztiR9kkogNzyD6uvXgRgt55ktUwR+Vjw6maPcie0WhwO2
hA9igSlbMdbWiTS/Pp4xOYgNA8WbGWdeIJhm1y3tMoQvvxvRtup3YU2IRplF
yxBmdHKjU991QdmyOirAfBip8VTnzPvdEBfXLciWDWNmpN/RxL8HesV8x1CD
EZwWd31q3NgLocGZ+Kq7I+D80PQ+b6QX1k6cqEQVMcJTE6KjOH2Qm7hiZ3zp
a7daTwcc7YejUtj7vEUx6h/tcln8ZQD+AxVuaSckmHtCud4yOoAXavWZwa8k
WMzQbQnPEcFsvIZ97I0Es+82qTbkiXCv5Jk9660E/q7vLqrmixCyhze/e1IC
rruA18wTQTY2a7/OrAS92/l2WnwRLHrOrhxfSuJivCiaJETIvWQRF6BBIk58
azxUmYDiQTOXw5okfBM8UxtXE4hSN950UJtEranYU30dAZsK3T5bHRJikcW2
XjUCw2/WG2oZkngqZa1vqUtgtfXi+CtrEp2lIRprdxGIU5jjD+4goR9/JTbC
gQDZP53ebUvie++Tk617CRSEjdnWO5CgOqUyEl0JONzpz8p2+fqjqjXP+BAo
DOw+lelKoi84bN7JjwDNqGPfDTcS/kPnLbmBBF41N1ASPEm8/bSC4RdGYHda
XU8si8TwnYhyPpvAE58aLtubhIP7/SC1KALr9J4lHD9CQlfzplncKQJnFko9
fX1JXFi2X73/NIH/ATXoZ0g=
      "]]}}, {{}, {}, 
    {RGBColor[0, 0, 1], LineBox[CompressedData["
1:eJwVkms0lHkAh1EjaSpFxUm7FIqkkkslfui0ym4MySXUZG0kuqglCRGtRaU6
bKdtTS617cklVmlsZehCl01DMS6D8c478zfMvLEiWWbbD895Pj/nPKZhR3x/
0NLQ0Nj5hf+da5HnuGEDhcbtZf5heiKctduXbGdHoa03bXbUQhHi3C2f2dtT
oOMCHx4zECFk7yPfjY4UtG7OWJ5qKMKqfDrG2YmCs+YeJc9EBAHLocRjK4Uq
vna6eJ0IjLRNL2QXhetW3MpAHxG8ig0V6ccpHGXrs70vi+BUNPUb+BQKjeLi
Lb/uwHM2Jyz4MwWXgIXB5UUdOKGYY7tusxQ+J5NfFazohK+n5YzDKVKo+uXD
pys68feiailbIMVYUlYuy7YL8d/faDTQppERZkUn1nVhm9nYvfeeNHKsrj0P
9ezGXQ6/ejyLxvTyrSaKV934Vi2uTntLYyJ8YG14oBhnjY4+O7xABm57U01m
lxi9tw5J9H1kcFlz8rVhSA+i7jdNm+bLUPza7Jd8eQ+urL9ifalNhrjrup2P
D/Zi0KYu+KCxHMJ7ypZ3Y704wPPL/8yVgxPEsb7j2oeUxADJ0A05VEe+We2W
2ofGjr/WbJfJ4bPaieXV2Iey4vNZM1cSKJK/W3pXR4IGRkAsYwg84kPYP3tL
MFke6qtTRbBv/7snS/IkEHZEvPYYISijImsm30vQdq7V44PjAOp0/83wN+5H
Aq+0deLkAGa9vJOrvbcfWM1EPXwwAN3EJ/ujC/phv+nKDDI1gJa3HXN3SPvR
VFVSesJFgebH7LFT5hSyy5cERKUrcNnHhb8imsKxZaNsoxcKzGtrd0+rpPBx
3Ll11exBNIyc13EYoXDbb+BaodcgNlW7/8HdKMUnY63jiRcHIZtrcK4vQYro
wHT/ec2D8Araabm4Vor5iyOddA2GkLgj4/cGtRQc8NdF+A2hR9ix/6krDYHi
0Mo114bAUluPb/mJxhvNHBu/riGEPb3epPuChjp7kZvdUiVOJZ7gXpgtw76r
GruTuEr4t3u2ZnrLEKoXlLCyUInIV85polwZsnT1eZvlSiB8LZvXLkNXvMNb
bWsVhKrMkKXL5HjBaVxgF63C1QVqMhQqh/wO3+9lpQpmHcc59iVyFMUa3a4f
UUG2kTXaQsuRWtqvHWvLIDxrOtTCguBi5NyI7HgGN8+MFX04RLA592Y76wED
4/lzzPTKCKxj7QXmfAZtBbeUtuUEy/ye3d5Wy+CCjXvN7goCDUM6IeMhAw2v
BM9fKwkaeWbGrHoGdI78qMV9Ar+KYu7MlwwqdJ8+2lJHcPhN4YBWNwM37aSA
g0IC7t31LcvFDCbzDE1zWgh8LtfXuvcw+NO8WlHeSmDnL8lO62Ngvm0wafQ9
waTYxEaLZjAnPagkpZMgU8mL1VQxeKfpMJxHEZxqXhtsyjA4f0lY+0BKEF1Z
t9XtAwMP05j0LprA+8deg9QRBnzXkiWmhGDx1Fc1GuMMClIWupQOEej0lvNM
PjEInF8+q1lJMCFwyXSdYKDH2yEcVhF0p+8NOjPJIO3xmXDHYYI3B1RuN/5l
sMXL2GbPl68F25OtBFMMPoprxk//Q1BlNU+/b/pLf8yuet4oQTG7YFKtZhAx
pcpq+EjwHwpff+s=
      "]]}}},
  AspectRatio->Automatic,
  Axes->True,
  AxesLabel->{
    FormBox[
    "\"x = \\!\\(\\*SqrtBox[\\(1 - \\*FractionBox[\\(\[Alpha]\\), \
\\(\[Beta]\\)]\\)]\\)\"", TraditionalForm], 
    FormBox[
    "\"\\!\\(\\*FractionBox[\\(1\\), SubscriptBox[\\(\[Beta]r\\), \
\\(c\\)]]\\)\"", TraditionalForm]},
  AxesOrigin->{0, 0},
  PlotRange->{{0, 0.7}, {0., 0.49223998779428574`}},
  PlotRangeClipping->True,
  PlotRangePadding->{
    Scaled[0.02], 
    Scaled[0.02]}]], "Output",
 CellChangeTimes->{3.38866306684375*^9, 3.38866333346875*^9, 
  3.407484947953125*^9, 3.407485574921875*^9, 3.4399107463125*^9}]
}, Open  ]],

Cell[TextData[{
 "Note that the definition of the cross sections in Peierls and the ",
 StyleBox["Primer",
  FontSlant->"Italic"],
 " differ, and \[Alpha] contains the scattering, absorption and the \
disintegration cross sections. In the ",
 StyleBox["Primer",
  FontSlant->"Italic"],
 " the transport cross section is the sum of the fission c.s and an angular \
integral over the scattering cross sections. Thus the values given above \
cannot be directly used in Peierls' formula. Nevertheless, if one tries to \
identify alfa with the reciprocal of the average mean free path, this gives:\n\
\n",
 Cell[BoxData[
  RowBox[{
   StyleBox[
    FormBox[
     SuperscriptBox["k", "2"],
     TraditionalForm], "Text"], 
   StyleBox["=",
    FontSlant->"Italic"], 
   StyleBox[" ",
    FontSlant->"Italic"], 
   RowBox[{
    StyleBox["3",
     FontSlant->"Italic"], 
    StyleBox[" ",
     FontSlant->"Italic"], 
    StyleBox[
     FormBox[
      SuperscriptBox["\[Alpha]", "2"],
      TraditionalForm], "Text"], 
    RowBox[{
     StyleBox["(",
      FontSlant->"Italic"], 
     StyleBox[
      RowBox[{
       StyleBox["1",
        FontSlant->"Italic"], 
       StyleBox["-",
        FontSlant->"Italic"], 
       FormBox[
        FractionBox["\[Alpha]", "\[Beta]"],
        TraditionalForm]}], "Text"], 
     StyleBox[")",
      FontSlant->"Italic"]}], 
    StyleBox[
     StyleBox[
      RowBox[{
       StyleBox[" ",
        FontWeight->"Bold"], 
       StyleBox[" ",
        FontVariations->{"CompatibilityType"->0}]}]], "Text"], 
    StyleBox[
     RowBox[{"(", "Peierls", ")"}],
     FontFamily->"Times New Roman",
     FontVariations->{"CompatibilityType"->0}]}]}]],
  FontFamily->"Times New Roman"],
 StyleBox["\n",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  StyleBox[
   RowBox[{
    FormBox[
     SuperscriptBox["k", "2"],
     TraditionalForm], "=", 
    RowBox[{
     RowBox[{
      FormBox[
       FractionBox["1", 
        SuperscriptBox["\[ScriptL]", "2"]],
       TraditionalForm], 
      FormBox[
       FractionBox[
        RowBox[{"3", 
         RowBox[{"(", 
          RowBox[{"\[Nu]", "-", "1"}], ")"}], 
         SubscriptBox["\[Sigma]", "f"]}], 
        SubscriptBox["\[Sigma]", "t"]],
       TraditionalForm]}], "=", 
     RowBox[{
      FormBox[
       SuperscriptBox[
        RowBox[{"(", 
         RowBox[{"n", " ", 
          SubscriptBox["\[Sigma]", "t"]}], ")"}], "2"],
       TraditionalForm], " ", 
      FormBox[
       FractionBox[
        RowBox[{"3", 
         RowBox[{"(", 
          RowBox[{"\[Nu]", "-", "1"}], ")"}], 
         SubscriptBox["\[Sigma]", "f"]}], 
        SubscriptBox["\[Sigma]", "t"]],
       TraditionalForm], 
      RowBox[{"(", "Primer", ")"}]}]}]}], "Text"]],
  FontFamily->"Times New Roman"],
 "\nIdentifying ",
 StyleBox["\[Alpha]",
  FontSlant->"Italic"],
 " and ",
 StyleBox["\[Beta] ",
  FontSlant->"Italic"],
 "gives \n\[Alpha] \[Rule] n ",
 Cell[BoxData[
  FormBox[
   SubscriptBox["\[Sigma]", "t"], TraditionalForm]]],
 "\nand\n",
 Cell[BoxData[
  RowBox[{
   RowBox[{"1", "-", 
    FormBox[
     FractionBox["\[Alpha]", "\[Beta]"],
     TraditionalForm]}], "\[Rule]", 
   FormBox[
    FractionBox[
     RowBox[{
      RowBox[{"(", 
       RowBox[{"\[Nu]", "-", "1"}], ")"}], 
      SubscriptBox["\[Sigma]", "f"]}], 
     SubscriptBox["\[Sigma]", "t"]],
    TraditionalForm]}]],
  FontFamily->"Times New Roman"],
 " or\n",
 StyleBox["\[Beta] \[Rule] ",
  FontFamily->"Times New Roman"],
 Cell[BoxData[
  FormBox[
   FormBox[
    FractionBox[
     RowBox[{"n", " ", 
      SubscriptBox["\[Sigma]", "t"]}], 
     RowBox[{"(", 
      RowBox[{"1", "-", 
       FractionBox[
        RowBox[{
         SubscriptBox["\[Sigma]", "f"], "(", 
         RowBox[{"\[Nu]", "-", "1"}], ")"}], 
        SubscriptBox["\[Sigma]", "t"]]}], ")"}]],
    TraditionalForm], TraditionalForm]],
  FontFamily->"Times New Roman"],
 "\n\nUsing the values above for the cross sections for U 235 one sees that \
Peierls original formula underestimates the critical radius as compared to \
Eq. 10, which is not surprising given the number of approximations involved. \
Also more accurate values for the different cross sections are probably \
needed to make an honest comparison beteeen the formulae\n"
}], "Text"],

Cell[CellGroupData[{

Cell[BoxData[{
 RowBox[{
  RowBox[{"alfa", " ", "=", " ", 
   RowBox[{"nd", " ", "st"}]}], " ", ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"beta", " ", "=", " ", 
   RowBox[{"nd", " ", 
    RowBox[{"st", "/", 
     RowBox[{"(", 
      RowBox[{"1", " ", "-", " ", 
       RowBox[{"sf", " ", 
        RowBox[{
         RowBox[{"(", 
          RowBox[{"ny", "-", "1"}], ")"}], "/", "st"}]}]}], ")"}]}]}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"x", " ", "=", " ", 
   RowBox[{"Sqrt", "[", 
    RowBox[{"1", "-", 
     RowBox[{"alfa", "/", "beta"}]}], "]"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"r", " ", "=", " ", 
   RowBox[{"1", "/", 
    RowBox[{"(", 
     RowBox[{"beta", " ", 
      RowBox[{"(", 
       RowBox[{
        RowBox[{"0.552", " ", "x"}], " ", "+", " ", 
        RowBox[{"0.216", " ", 
         RowBox[{"x", "^", "2"}]}]}], ")"}]}], ")"}]}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"cm", " ", "=", " ", 
   RowBox[{"rho", " ", "4", " ", "Pi", " ", 
    RowBox[{
     RowBox[{
      RowBox[{"r", "^", "3"}], "/", "3"}], " ", "/", "1000"}]}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"rho235", " ", "=", " ", 
   RowBox[{"rhoList", "[", 
    RowBox[{"[", "1", "]"}], "]"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"rho239", " ", "=", " ", 
   RowBox[{"rhoList", "[", 
    RowBox[{"[", "2", "]"}], "]"}]}], ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"MmodernU235", " ", "=", " ", "56.0"}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"MmodernPU239", " ", "=", " ", "11.0"}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{"RmodernU235", " ", "=", " ", 
   RowBox[{
    RowBox[{
     RowBox[{"(", 
      RowBox[{"3", " ", 
       RowBox[{"MmodernU235", "/", 
        RowBox[{"(", 
         RowBox[{"rho235", " ", "*", "1000", " ", "4", " ", "Pi"}], ")"}]}]}],
       ")"}], "^", 
     RowBox[{"(", 
      RowBox[{"1", "/", "3"}], ")"}]}], "*", "100"}]}], 
  ";"}], "\[IndentingNewLine]", 
 RowBox[{
  RowBox[{
   RowBox[{"RmodernPU239", " ", "=", " ", 
    RowBox[{
     RowBox[{
      RowBox[{"(", 
       RowBox[{"3", " ", 
        RowBox[{"MmodernPU239", "/", 
         RowBox[{"(", 
          RowBox[{"rho239", " ", "*", "1000", " ", "4", " ", "Pi"}], 
          ")"}]}]}], ")"}], "^", 
      RowBox[{"(", 
       RowBox[{"1", "/", "3"}], ")"}]}], "*", "100"}]}], ";"}], 
  "\[IndentingNewLine]"}], "\[IndentingNewLine]", 
 RowBox[{"StandardForm", "[", 
  RowBox[{"TableForm", "[", 
   RowBox[{
    RowBox[{"{", 
     RowBox[{
     "x", ",", "r", ",", "cm", ",", "RmodernU235", ",", "MmodernU235", ",", 
      "RmodernPU239", ",", "MmodernPU239"}], "}"}], ",", 
    RowBox[{"TableHeadings", "\[Rule]", " ", 
     RowBox[{"{", 
      RowBox[{
       RowBox[{"{", 
        RowBox[{
        "\"\<x = \!\(\*SqrtBox[\(1 - \*FractionBox[\(\[Alpha]\), \(\[Beta]\)]\
\)]\)\>\"", ",", " ", "\"\<Peierls critical radius[cm]\>\"", ",", 
         "\"\<Peierls critical Mass [kg]\>\"", ",", 
         "\"\<Modern value of radius (U235)\>\"", ",", 
         "\"\<Modern value of mass (U235)\>\"", ",", 
         "\"\<Modern value of radius (PU239)\>\"", ",", 
         "\"\<Modern value of mass (PU239)\>\""}], "}"}], ",", "None"}], 
      "}"}]}], ",", 
    RowBox[{"TableSpacing", " ", "\[Rule]", " ", 
     RowBox[{"{", 
      RowBox[{"1", ",", "10"}], "}"}]}]}], "]"}], 
  "]"}], "\[IndentingNewLine]"}], "Input",
 CellChangeTimes->{3.388663346171875*^9}],

Cell[BoxData[
 TagBox[
  TagBox[GridBox[{
     {
      TagBox["\<\"x = \\!\\(\\*SqrtBox[\\(1 - \\*FractionBox[\\(\[Alpha]\\), \
\\(\[Beta]\\)]\\)]\\)\"\>",
       HoldForm], "0.680881781221968`"},
     {
      TagBox["\<\"Peierls critical radius[cm]\"\>",
       HoldForm], "5.817031522676345`"},
     {
      TagBox["\<\"Peierls critical Mass [kg]\"\>",
       HoldForm], "15.583128978612795`"},
     {
      TagBox["\<\"Modern value of radius (U235)\"\>",
       HoldForm], "8.910030783813173`"},
     {
      TagBox["\<\"Modern value of mass (U235)\"\>",
       HoldForm], "56.`"},
     {
      TagBox["\<\"Modern value of radius (PU239)\"\>",
       HoldForm], "5.096290033186184`"},
     {
      TagBox["\<\"Modern value of mass (PU239)\"\>",
       HoldForm], "11.`"}
    },
    GridBoxAlignment->{
     "Columns" -> {{Left}}, "ColumnsIndexed" -> {}, "Rows" -> {{Baseline}}, 
      "RowsIndexed" -> {}},
    GridBoxDividers->{
     "Columns" -> {False, {True}, False}, "ColumnsIndexed" -> {}, 
      "Rows" -> {{False}}, "RowsIndexed" -> {}},
    GridBoxSpacings->{"Columns" -> {
        Offset[0.27999999999999997`], {
         Offset[0.5599999999999999]}, 
        Offset[0.27999999999999997`]}, "ColumnsIndexed" -> {}, "Rows" -> {
        Offset[0.2], {
         Offset[0.4]}, 
        Offset[0.2]}, "RowsIndexed" -> {}}],
   OutputFormsDump`HeadedColumn],
  Function[BoxForm`e$, 
   TableForm[
   BoxForm`e$, 
    TableHeadings -> {{
      "x = \!\(\*SqrtBox[\(1 - \*FractionBox[\(\[Alpha]\), \(\[Beta]\)]\)]\)",
        "Peierls critical radius[cm]", "Peierls critical Mass [kg]", 
       "Modern value of radius (U235)", "Modern value of mass (U235)", 
       "Modern value of radius (PU239)", "Modern value of mass (PU239)"}, 
      None}, TableSpacing -> {1, 10}]]]], "Output",
 CellChangeTimes->{3.3886630674375*^9, 3.388663347703125*^9, 
  3.40748494815625*^9, 3.407485574984375*^9, 3.439910746390625*^9}]
}, Open  ]],

Cell[TextData[{
 "When Otto Frisch and Rudolf Peierls first used the formula the results were \
astounding: Instead of a few tons of uranium, less than a few tens of kilos \
was needed. In the famous Frisch-Peierls memorandum from March 1940 the \
figure 2.1 cm is mentioned, which was one of the reasons the memorandum \
become a serious concern in the US and the UK during the war. Apart from the \
fact that Peierls' formula underestimates the critical mass,  it turns out \
they used the wrong numerical values for many of variables involved, \
including the density, the fission cross section and the number of neutrons \
per fission and in reality the mass was larger. A modern value for the \
critical mass for a bare sphere of U235 is 56 kg which corresponds to a \
radius of 8.9 cm. In practice the \"device\" (=bomb) the critical mass was \
determined by measuring neutron multiplication from a neutron source \
surrounded by a sphere of the radioactive material. By making the sphere \
successively larger, the multiplication could be extrapolated and the \
critical mass measured (",
 StyleBox["Primer",
  FontSlant->"Italic"],
 ", page 33). "
}], "Text"]
}, Open  ]],

Cell[CellGroupData[{

Cell["References", "Section"],

Cell[TextData[{
 "[1] Rudolf Peierls, 1939. \"",
 StyleBox["Critical conditions in neutron multiplication\". ",
  FontSlant->"Italic"],
 "Proc. Camb. Phil. Soc. 35:610. \n[2] Robert Serber, Edited with an \
Introduction by Richard Rhodes (1992) \"",
 StyleBox["The Los Alamos Primer\".",
  FontSlant->"Italic"],
 " University of California Press, ISBN-0-520-07576-5. \n[3] Richard Rhodes, \
1986. \"",
 StyleBox["The Making of the Atom Bomb\".",
  FontSlant->"Italic"],
 "Touchstone Book, Simon & Schuster New York, ISBN 0-648-81378-5. \n[4] \
Francis Perrin, 1939.\"",
 StyleBox["Calul relatif aux conditions \[EAcute]venetuelles de transmutation \
en ch\[AHat]ine de l'uranium. ",
  FontSlant->"Italic"],
 "Comptes Rendus, 208:1394."
}], "Text"]
}, Open  ]]
}, Open  ]]
},
ScreenStyleEnvironment->"Working",
WindowSize->{900, 856},
WindowMargins->{{82, Automatic}, {Automatic, 11}},
CellLabelAutoDelete->True,
FrontEndVersion->"7.0 for Microsoft Windows (32-bit) (November 10, 2008)",
StyleDefinitions->"Default.nb"
]
(* End of Notebook Content *)

(* Internal cache information *)
(*CellTagsOutline
CellTagsIndex->{}
*)
(*CellTagsIndex
CellTagsIndex->{}
*)
(*NotebookFileOutline
Notebook[{
Cell[CellGroupData[{
Cell[567, 22, 1944, 42, 405, "Title"],
Cell[CellGroupData[{
Cell[2536, 68, 31, 0, 71, "Section"],
Cell[2570, 70, 13397, 352, 921, "Text"],
Cell[CellGroupData[{
Cell[15992, 426, 6117, 176, 859, "Input"],
Cell[22112, 604, 2314, 57, 133, "Output"]
}, Open  ]]
}, Open  ]],
Cell[CellGroupData[{
Cell[24475, 667, 32, 0, 71, "Section"],
Cell[24510, 669, 21896, 715, 1599, "Text"],
Cell[CellGroupData[{
Cell[46431, 1388, 1821, 55, 132, "Input"],
Cell[48255, 1445, 1148, 30, 42, "Output"]
}, Open  ]]
}, Open  ]],
Cell[CellGroupData[{
Cell[49452, 1481, 65, 0, 71, "Section"],
Cell[49520, 1483, 9780, 303, 813, "Text"],
Cell[CellGroupData[{
Cell[59325, 1790, 3100, 92, 337, "Input"],
Cell[62428, 1884, 179, 2, 30, "Output"],
Cell[62610, 1888, 176, 2, 30, "Output"],
Cell[62789, 1892, 8653, 149, 447, "Output"]
}, Open  ]],
Cell[71457, 2044, 2301, 78, 315, "Text"],
Cell[CellGroupData[{
Cell[73783, 2126, 1891, 53, 192, "Input"],
Cell[75677, 2181, 372, 10, 50, "Output"],
Cell[76052, 2193, 1605, 50, 57, "Output"],
Cell[77660, 2245, 368, 9, 30, "Output"],
Cell[78031, 2256, 7239, 126, 447, "Output"]
}, Open  ]],
Cell[85285, 2385, 27, 0, 31, "Input"],
Cell[CellGroupData[{
Cell[85337, 2389, 1548, 47, 165, "Input"],
Cell[86888, 2438, 670, 15, 61, "Output"],
Cell[87561, 2455, 1800, 45, 125, "Output"],
Cell[89364, 2502, 4280, 79, 195, "Output"]
}, Open  ]]
}, Open  ]],
Cell[CellGroupData[{
Cell[93693, 2587, 53, 0, 71, "Section"],
Cell[93749, 2589, 1714, 58, 253, "Text"],
Cell[CellGroupData[{
Cell[95488, 2651, 631, 18, 72, "Input"],
Cell[96122, 2671, 371, 10, 50, "Output"]
}, Open  ]],
Cell[96508, 2684, 193, 4, 29, "Text"],
Cell[CellGroupData[{
Cell[96726, 2692, 394, 11, 31, "Input"],
Cell[97123, 2705, 395, 11, 48, "Output"]
}, Open  ]],
Cell[97533, 2719, 903, 38, 52, "Text"],
Cell[CellGroupData[{
Cell[98461, 2761, 2142, 61, 192, "Input"],
Cell[100606, 2824, 1143, 35, 80, "Output"],
Cell[101752, 2861, 396, 11, 52, "Output"],
Cell[102151, 2874, 3916, 72, 243, "Output"]
}, Open  ]],
Cell[106082, 2949, 2219, 75, 225, "Text"],
Cell[108304, 3026, 1662, 59, 45, "Text"],
Cell[CellGroupData[{
Cell[109991, 3089, 1603, 50, 72, "Input"],
Cell[111597, 3141, 352, 10, 56, "Output"]
}, Open  ]],
Cell[111964, 3154, 4106, 159, 331, "Text"],
Cell[CellGroupData[{
Cell[116095, 3317, 1836, 53, 152, "Input"],
Cell[117934, 3372, 803, 26, 47, "Output"]
}, Open  ]],
Cell[118752, 3401, 2143, 67, 159, "Text"],
Cell[CellGroupData[{
Cell[120920, 3472, 2301, 71, 192, "Input"],
Cell[123224, 3545, 237, 5, 47, "Output"],
Cell[123464, 3552, 349, 10, 49, "Output"],
Cell[123816, 3564, 360, 8, 30, "Output"],
Cell[124179, 3574, 188, 2, 30, "Output"],
Cell[124370, 3578, 186, 2, 30, "Output"],
Cell[124559, 3582, 6337, 109, 243, "Output"],
Cell[130899, 3693, 9607, 163, 229, "Output"]
}, Open  ]],
Cell[140521, 3859, 3314, 118, 280, "Text"],
Cell[CellGroupData[{
Cell[143860, 3981, 1815, 55, 112, "Input"],
Cell[145678, 4038, 425, 10, 33, "Output"],
Cell[146106, 4050, 174, 2, 30, "Output"]
}, Open  ]],
Cell[146295, 4055, 1574, 60, 166, "Text"],
Cell[147872, 4117, 26, 0, 31, "Input"],
Cell[147901, 4119, 333, 5, 47, "Text"],
Cell[CellGroupData[{
Cell[148259, 4128, 1324, 38, 95, "Input"],
Cell[149586, 4168, 2149, 41, 248, "Output"],
Cell[151738, 4211, 2159, 41, 232, "Output"],
Cell[153900, 4254, 4113, 78, 272, "Output"]
}, Open  ]],
Cell[158028, 4335, 4256, 146, 363, "Text"],
Cell[CellGroupData[{
Cell[162309, 4485, 3469, 103, 392, "Input"],
Cell[165781, 4590, 1921, 50, 144, "Output"]
}, Open  ]],
Cell[167717, 4643, 1167, 19, 137, "Text"]
}, Open  ]],
Cell[CellGroupData[{
Cell[168921, 4667, 29, 0, 71, "Section"],
Cell[168953, 4669, 747, 18, 83, "Text"]
}, Open  ]]
}, Open  ]]
}
]
*)

(* End of internal cache information *)