(***********************************************************************
Mathematica-Compatible Notebook
This notebook can be used on any computer system with Mathematica 3.0,
MathReader 3.0, or any compatible application. The data for the notebook
starts with the line of stars above.
To get the notebook into a Mathematica-compatible application, do one of
the following:
* Save the data starting with the line of stars above into a file
with a name ending in .nb, then open the file inside the application;
* Copy the data starting with the line of stars above to the
clipboard, then use the Paste menu command inside the application.
Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode. Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).
NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing the
word CacheID, otherwise Mathematica-compatible applications may try to
use invalid cache data.
For more information on notebooks and Mathematica-compatible
applications, contact Wolfram Research:
web: http://www.wolfram.com
email: info@wolfram.com
phone: +1-217-398-0700 (U.S.)
Notebook reader applications are available free of charge from
Wolfram Research.
***********************************************************************)
(*CacheID: 232*)
(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[ 15176, 479]*)
(*NotebookOutlinePosition[ 17719, 564]*)
(* CellTagsIndexPosition[ 17405, 549]*)
(*WindowFrame->Normal*)
Notebook[{
Cell[CellGroupData[{
Cell["A Newtonian model for a simple ionic reaction", "Title",
ShowCellTags->False,
CellTags->"title"],
Cell["METRIC Project", "Subtitle"],
Cell["\<\
This problem forms part of your coursework, and should be handed in to your \
tutor by Friday 27th March. Hand in your answers, fully annotated to your \
tutor. Discuss with him/her in what format it is required. From past \
experience you would be well advised to save multiple copies of your work, \
frequently. \
\>", "Subsubtitle",
Evaluatable->False,
CellHorizontalScrolling->False,
TextAlignment->Center],
Cell[CellGroupData[{
Cell["1. Objectives", "Section",
ShowCellTags->False,
CellTags->"objectives"],
Cell[TextData[{
"In this exercise, you are asked to set up and explore a simulation of the \
reaction between a free ",
Cell[BoxData[
\(TraditionalForm\`\(Cl\^-\)\)]],
" ion and a bonded NaCl ion pair. You're going to be using a ",
StyleBox["Newtonian",
FontSlant->"Italic"],
" model: that is to say, the reaction is modelled as an interaction between \
charges moving under forces, and their motion is predicted using Newton's \
second law, ",
StyleBox["F",
FontSlant->"Italic"],
" = ",
StyleBox["ma",
FontSlant->"Italic"],
". For simplicity, the charges will move and lie along a straight line, \
which allows you to treat the problem as a one-dimensional one.\n\nThis \
reaction is a very simple one in chemical terms, but it is still too complex \
to be predicted using traditional pencil-and-paper mathematics, or indeed any \
of ",
StyleBox["Mathematica",
FontSlant->"Italic"],
"'s exact, symbolic functions. Instead, you'll be using a method of ",
StyleBox["numerical approximation",
FontSlant->"Italic"],
", which will of necessity be an imperfect one. One of your objectives is \
to find out how well this numerical approximation seems to model physical \
reality, by finding out whether energy appears to be conserved within the \
model.\n\nDepending on the positions and velocities of the three ions when \
they begin interacting, sometimes a reaction will occur and sometimes it will \
not. A third outcome is also possible: under certain conditions, all three \
ions continue to interact with each other, with no one ion becoming \
dissociated from the other two. Your other main objective is to set the \
initial conditions which will produce these three outcomes and devise a way \
of deciding which has occurred."
}], "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["2. Preliminaries", "Section"],
Cell[TextData[{
"In this section you'll find the ",
StyleBox["Mathematica",
FontSlant->"Italic"],
" code you'll need before you can set up your simulation. However, you \
should note that we have set up the code to make a movie of the reaction, but \
for your report we want you to alter the code to produce plots of position as \
a function of elapsed time. To save you a lot of typing we have placed this \
",
StyleBox["Mathematica",
FontSlant->"Italic"],
" file, ",
StyleBox["MathLabQue2.nb",
FontFamily->"Helvetica"],
", on the hard disk of the Maths Lab computers in the ",
StyleBox["Mathematica",
FontSlant->"Italic"],
" Folder."
}], "Text"],
Cell[CellGroupData[{
Cell["Physical constants", "Subsection",
ShowCellTags->False,
CellTags->"constants"],
Cell[TextData[{
"The first priority is to inform ",
StyleBox["Mathematica",
FontSlant->"Italic"],
" of the values of the various physical constants that the simulation \
requires. We used the Basic Input Palette to produce the suffices, powers and \
Greek letters you see here.\n\nThe permittivity of the vacuum:"
}], "Text"],
Cell[BoxData[
\(\[Epsilon]\_0 = \ 8.85419*^-12\)], "Input"],
Cell["Avogadro's number:", "Text"],
Cell[BoxData[
\(N\_A = \ 6.02205*^23\)], "Input"],
Cell["The relative atomic masses of sodium and chlorine:", "Text"],
Cell[BoxData[{
\(RAM\_Na = \ 22.9898\),
\(RAM\_Cl = 34.9688\)}], "Input"],
Cell["The charge, in Coulombs, on an electron:", "Text"],
Cell[BoxData[
\(e\ = \ 1.60219*^-19\)], "Input"],
Cell["The equilibrium bond length of the NaCl dipole, in metres:", "Text"],
Cell[BoxData[
\(r\_0 = \ 2.5*^-10\)], "Input"],
Cell["\<\
You can now use the relative atomic masses, together with Avogadro's number, \
to calculate the masses of the sodium and chloride ions in kg.\
\>", "Text"],
Cell[BoxData[{
\(m\_Na = \ RAM\_Na\/\(\(10\^3\) N\_A\)\),
\(m\_Cl = \ RAM\_Cl\/\(\(10\^3\) N\_A\)\)}], "Input"]
}, Open ]],
Cell[CellGroupData[{
Cell["Attraction and repulsion forces", "Subsection"],
Cell["\<\
Since, as you'll see, the two chloride ions never approach each other very \
closely, we can model the repulsive force between them using a simple Coulomb \
(electrostatic) repulsion. In newtons, for two negative charges separated by \
a displacement r metres:\
\>", "Text"],
Cell[BoxData[
\(F\_mm[r_]\ =
\(-\(e\^2\/\(4\ \[Pi]\ \(\[Epsilon]\_\(0\ \)\) r\^2\)\)\)\)], "Input"],
Cell["\<\
Between the sodium ion and either one of the sodium ions, on the other hand, \
there is a force that is attractive when the ions are far apart and repulsive \
when they are close together. This force is extremely well modelled by the \
following formula:\
\>", "Text"],
Cell[BoxData[
\(F\_pm[r_]\ =
\(e\^2\/\(4\ \[Pi]\ \[Epsilon]\_0\)\)
\((1\/r\^2 - \((r\_0)\)\^8\/r\^10)\)\)], "Input"]
}, Open ]],
Cell[CellGroupData[{
Cell["Time", "Subsection"],
Cell[TextData[{
"You need to tell ",
StyleBox["Mathematica",
FontSlant->"Italic"],
" the duration of the simulated reaction, in seconds. For different \
reactions under different sets of starting conditions, different durations \
may well be appropriate, but here is a suggestion."
}], "Text"],
Cell[BoxData[
StyleBox[\(T = 1.2*^-12\),
FontColor->GrayLevel[0]]], "Input",
FontColor->RGBColor[1, 0, 0]]
}, Open ]],
Cell[CellGroupData[{
Cell["Initial values of the variables", "Subsection"],
Cell["\<\
Initial positions (in metres) and velocities (in metres per second) of the \
three ions. By changing these values you can observe different behaviours.\
\>", "Text"],
Cell[BoxData[{\(x1\_0 = \(-5.5\) r\_0\), \(x2\_0 = \ 0.0\),
\(x3\_0 = \ 1.1 r\_0\), \(v1\_0 = \ 5000.0\),
StyleBox[\(v2\_0 = \ 0.0\),
FontColor->GrayLevel[0]],
StyleBox[\(v3\_0 = \ 0.0\),
FontColor->GrayLevel[0]]}], "Input",
CellTags->"initvals"],
Cell["\<\
You should save your notebook at this point so that you can retrieve the \
preliminary code without retyping it.\
\>", "Text"]
}, Open ]]
}, Open ]],
Cell[CellGroupData[{
Cell["3. Simulating the reaction", "Section"],
Cell[CellGroupData[{
Cell["Solving the equations of motion", "Subsection"],
Cell[TextData[{
"Each of the three ions moves under the influence of the other two. More \
precisely, the force on each, and therefore the acceleration of each, depends \
on its position relative to the other two. This gives three coupled \
second-order equations of motion, which we solve approximately using ",
StyleBox["Mathematica",
FontSlant->"Italic"],
"'s ",
StyleBox["NDSolve",
FontFamily->"Courier"],
" function, subject to the initial conditions specified above."
}], "Text"],
Cell[BoxData[
\(\(xrules =
NDSolve[{\n\t\t\t
\(\(x1'\)'\)[t] ==
\(1\/m\_Cl\((F\_pm[x2[t] - x1[t]]\ + \ F\_mm[x3[t] - x1[t]])
\)\), \n\t\t\t
\(\(x2'\)'\)[t] ==
\(1\/m\_Na\)
\((\ \(-F\_pm[x2[t] - x1[t]]\)\ + \ F\_pm[x3[t] - x2[t]])\), \n
\t\t\t\(\(x3'\)'\)[t] ==
\(1\/m\_Cl\)
\((\(-F\_mm[x3[t] - x1[t]]\)\ - \ F\_pm[x3[t] - x2[t]])\), \n
\t\t\tx1[0] == x1\_0, x2[0] == x2\_0, x3[0] == x3\_0, \n\t\t\t
\(x1'\)[0] == v1\_0, \(x2'\)[0] == v2\_0, \(x3'\)[0] == v3\_0}, \n
\t\t{x1, x2, x3}, {t, 0, T}, MaxSteps -> 2000]; \)\)], "Input",
ShowCellTags->False,
CellTags->"diffeq"],
Cell[TextData[{
"The solutions are given in the form of three \"hidden\" functions of ",
StyleBox["t",
FontSlant->"Italic"],
", defined only for 0 \[LessEqual] ",
StyleBox["t",
FontSlant->"Italic"],
" \[LessEqual] ",
StyleBox["T",
FontSlant->"Italic"],
"."
}], "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["Positions and velocities in physical space", "Subsection"],
Cell[TextData[{
"We need a way of getting actual values out of the three \"hidden\" \
functions that came out of the call to ",
StyleBox["NDSolve",
FontFamily->"Courier"],
"."
}], "Text"],
Cell[BoxData[{
\(pos1[t_]\ := \ \((x1[t] /. xrules)\)[\([1]\)]\),
\(pos2[t_]\ := \ \((x2[t] /. xrules)\)[\([1]\)]\),
\(pos3[t_]\ := \ \((x3[t] /. xrules)\)[\([1]\)]\),
\(vel1[t_]\ := \ \((\(x1'\)[t] /. xrules)\)[\([1]\)]\),
\(vel2[t_]\ := \ \((\(x2'\)[t] /. xrules)\)[\([1]\)]\),
\(vel3[t_]\ := \ \((\(x3'\)[t] /. xrules)\)[\([1]\)]\)}], "Input"],
Cell["\<\
You should satisfy yourself that these definitions have worked, by reading \
values with commands like\
\>", "Text"],
Cell[BoxData[{
\(pos1[0.5 T]\),
\(vel3[0.3 T]\)}], "Input"],
Cell["and so on.", "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["Creating an animation", "Subsection"],
Cell["\<\
You are now in a position to create a 401-frame animation of the motion of \
the three ions:\
\>", "Text"],
Cell[BoxData[
\(threeIonFrame[t_]\ :=
Show[Graphics[{AbsolutePointSize[15], \n
\t\t\t\t{RGBColor[0.8, 0.8, 0], Point[{pos1[t], 0}]}, {
RGBColor[0.4, 0.4, 0.4], Point[{pos2[t], 0}]}, {
RGBColor[0.8, 0.8, 0], Point[{pos3[t], 0}]}}], \n\t\t
PlotRange -> {{\(-15\) r\_0, 15 r\_0}, Automatic}]\)], "Input"],
Cell[BoxData[
\(\(Table[threeIonFrame[t], {t, 0, T, T/100}]; \)\)], "Input",
CellTags->"animation"],
Cell[TextData[{
"To view this animation, click on the cell bracket that groups all the \
frames together and choose ",
StyleBox["Animate Selected Graphics",
FontFamily->"Helvetica",
FontWeight->"Bold"],
" from the ",
StyleBox["Cell",
FontFamily->"Helvetica",
FontWeight->"Bold"],
" menu."
}], "Text"]
}, Open ]]
}, Open ]],
Cell[CellGroupData[{
Cell["4. The problem", "Section"],
Cell[CellGroupData[{
Cell["(i) Recoding", "Subsection",
CellTags->"prob1"],
Cell["\<\
Your first task is to re-code the example program so that you can plot the \
positions of each ion as a function of time. Make sure that your plots are \
clearly and unambiguously labelled.\
\>", "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["(ii) The three possible outcomes", "Subsection",
CellTags->"probii"],
Cell[TextData[{
"Now find initial conditions for which a reaction appears to occur; for \
which it does ",
StyleBox["not",
FontSlant->"Italic"],
" appear to occur; and a third set for which neither chloride ion appears \
to dissociate, and the three-ion interaction appears to continue \
indefinitely.\n\nYou should feel free to change the values of ",
Cell[BoxData[
StyleBox[\(x1\_0\),
FontWeight->"Bold"]]],
", ",
Cell[BoxData[
StyleBox[\(x2\_0\),
FontWeight->"Bold"]]],
", ",
Cell[BoxData[
StyleBox[\(x3\_0\),
FontWeight->"Bold"]]],
" and ",
Cell[BoxData[
StyleBox[\(v1\_0\),
FontWeight->"Bold"]]],
", but you should probably leave ",
Cell[BoxData[
StyleBox[\(v2\_0\),
FontWeight->"Bold"]]],
" and ",
Cell[BoxData[
StyleBox[\(v3\_0\),
FontWeight->"Bold"]]],
" set at zero in the first instance."
}], "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["(iii) Energy conservation", "Subsection",
CellTags->"prob3"],
Cell[TextData[{
"The total mechanical energy of the system is formed of the ",
StyleBox["kinetic",
FontSlant->"Italic"],
" energies of each of the three particles (which you can calculate using \
the velocity functions you have defined) and the ",
StyleBox["potential",
FontSlant->"Italic"],
" energy that arises from their relative positions in the force field \
(which you can calculate using the position functions you have defined, based \
on what you know about potentials).\n\nThe simulation is based on a numerical \
approximation. In a fundamental sense, then, it's not to be trusted, since no \
numerical method is exact. If we found, for example, that energy appeared not \
be be conserved, at least approximately, in our model, then we should reject \
it.\n\nInvestigate how well energy appears to be conserved within this model. \
Report on the performance of the simulation in this regard."
}], "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell["(iv) A test for dissociation", "Subsection",
CellTags->"prob4"],
Cell[TextData[{
"In part (ii) of this problem, you uncovered three different types of \
apparent behaviour of the system. However, there's a problem. How can we \
tell, if an ion ",
StyleBox["appears",
FontSlant->"Italic"],
" to be dissociated, whether it really is?\n\nDevise a test for deciding \
whether, when dissociation appears to have occurred, it really has. Your test \
need not be perfect but it should take you further than mere observation. \
Apply your test to your \"reaction\" and \"no reaction\" examples from part \
(ii), and comment on your findings"
}], "Text"]
}, Open ]]
}, Open ]]
}, Open ]]
},
FrontEndVersion->"Microsoft Windows 3.0",
ScreenRectangle->{{0, 640}, {0, 424}},
ScreenStyleEnvironment->"Presentation",
WindowSize->{489, 333},
WindowMargins->{{Automatic, 4}, {5, Automatic}},
PrintingCopies->1,
PrintingPageRange->{1, Automatic},
PageHeaders->{{Cell[
TextData[ {
ValueBox[ "FileName"]}], "Header"], Inherited, Cell[
TextData[ {
CounterBox[ "Page"]}], "PageNumber"]}, {Cell[
TextData[ {
ValueBox[ "FileName"]}], "Header"], Inherited, Cell[
TextData[ {
CounterBox[ "Page"]}], "PageNumber"]}},
PageFooters->{{Cell[
TextData[ "rt"], "Footer"], Cell[
TextData[ "v1.1"], "Footer"], Cell[
TextData[ "10/97"], "Footer"]}, {Cell[
TextData[ "rt"], "Footer"], Cell[
TextData[ "v1.1"], "Footer"], Cell[
TextData[ "10/97"], "Footer"]}},
PageFooterLines->{True, True}
]
(***********************************************************************
Cached data follows. If you edit this Notebook file directly, not using
Mathematica, you must remove the line containing CacheID at the top of
the file. The cache data will then be recreated when you save this file
from within Mathematica.
***********************************************************************)
(*CellTagsOutline
CellTagsIndex->{
"title"->{
Cell[1731, 51, 106, 2, 213, "Title",
CellTags->"title"]},
"objectives"->{
Cell[2328, 70, 81, 2, 72, "Section",
CellTags->"objectives"]},
"constants"->{
Cell[4987, 139, 88, 2, 68, "Subsection",
CellTags->"constants"]},
"initvals"->{
Cell[7935, 242, 274, 6, 153, "Input",
CellTags->"initvals"]},
"diffeq"->{
Cell[9029, 278, 729, 16, 459, "Input",
CellTags->"diffeq"]},
"animation"->{
Cell[11538, 359, 105, 2, 61, "Input",
CellTags->"animation"]},
"prob1"->{
Cell[12079, 384, 55, 1, 68, "Subsection",
CellTags->"prob1"]},
"probii"->{
Cell[12387, 396, 76, 1, 68, "Subsection",
CellTags->"probii"]},
"prob3"->{
Cell[13431, 436, 68, 1, 68, "Subsection",
CellTags->"prob3"]},
"prob4"->{
Cell[14472, 460, 71, 1, 68, "Subsection",
CellTags->"prob4"]}
}
*)
(*CellTagsIndex
CellTagsIndex->{
{"title", 16519, 515},
{"objectives", 16605, 518},
{"constants", 16695, 521},
{"initvals", 16787, 524},
{"diffeq", 16873, 527},
{"animation", 16961, 530},
{"prob1", 17047, 533},
{"probii", 17134, 536},
{"prob3", 17221, 539},
{"prob4", 17307, 542}
}
*)
(*NotebookFileOutline
Notebook[{
Cell[CellGroupData[{
Cell[1731, 51, 106, 2, 213, "Title",
CellTags->"title"],
Cell[1840, 55, 34, 0, 76, "Subtitle"],
Cell[1877, 57, 426, 9, 254, "Subsubtitle"],
Cell[CellGroupData[{
Cell[2328, 70, 81, 2, 72, "Section",
CellTags->"objectives"],
Cell[2412, 74, 1793, 36, 620, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[4242, 115, 35, 0, 72, "Section"],
Cell[4280, 117, 682, 18, 195, "Text"],
Cell[CellGroupData[{
Cell[4987, 139, 88, 2, 68, "Subsection",
CellTags->"constants"],
Cell[5078, 143, 335, 7, 170, "Text"],
Cell[5416, 152, 63, 1, 38, "Input"],
Cell[5482, 155, 34, 0, 45, "Text"],
Cell[5519, 157, 53, 1, 38, "Input"],
Cell[5575, 160, 66, 0, 45, "Text"],
Cell[5644, 162, 82, 2, 61, "Input"],
Cell[5729, 166, 56, 0, 45, "Text"],
Cell[5788, 168, 53, 1, 38, "Input"],
Cell[5844, 171, 74, 0, 45, "Text"],
Cell[5921, 173, 50, 1, 38, "Input"],
Cell[5974, 176, 165, 3, 70, "Text"],
Cell[6142, 181, 120, 2, 114, "Input"]
}, Open ]],
Cell[CellGroupData[{
Cell[6299, 188, 53, 0, 68, "Subsection"],
Cell[6355, 190, 284, 5, 120, "Text"],
Cell[6642, 197, 111, 2, 66, "Input"],
Cell[6756, 201, 278, 5, 120, "Text"],
Cell[7037, 208, 137, 3, 68, "Input"]
}, Open ]],
Cell[CellGroupData[{
Cell[7211, 216, 26, 0, 68, "Subsection"],
Cell[7240, 218, 303, 7, 120, "Text"],
Cell[7546, 227, 118, 3, 38, "Input"]
}, Open ]],
Cell[CellGroupData[{
Cell[7701, 235, 53, 0, 68, "Subsection"],
Cell[7757, 237, 175, 3, 70, "Text"],
Cell[7935, 242, 274, 6, 153, "Input",
CellTags->"initvals"],
Cell[8212, 250, 136, 3, 70, "Text"]
}, Open ]]
}, Open ]],
Cell[CellGroupData[{
Cell[8397, 259, 45, 0, 72, "Section"],
Cell[CellGroupData[{
Cell[8467, 263, 53, 0, 68, "Subsection"],
Cell[8523, 265, 503, 11, 170, "Text"],
Cell[9029, 278, 729, 16, 459, "Input",
CellTags->"diffeq"],
Cell[9761, 296, 294, 11, 70, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[10092, 312, 64, 0, 68, "Subsection"],
Cell[10159, 314, 197, 6, 70, "Text"],
Cell[10359, 322, 381, 6, 153, "Input"],
Cell[10743, 330, 126, 3, 70, "Text"],
Cell[10872, 335, 70, 2, 61, "Input"],
Cell[10945, 339, 26, 0, 45, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[11008, 344, 43, 0, 68, "Subsection"],
Cell[11054, 346, 116, 3, 70, "Text"],
Cell[11173, 351, 362, 6, 245, "Input"],
Cell[11538, 359, 105, 2, 61, "Input",
CellTags->"animation"],
Cell[11646, 363, 326, 11, 95, "Text"]
}, Open ]]
}, Open ]],
Cell[CellGroupData[{
Cell[12021, 380, 33, 0, 72, "Section"],
Cell[CellGroupData[{
Cell[12079, 384, 55, 1, 68, "Subsection",
CellTags->"prob1"],
Cell[12137, 387, 213, 4, 95, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[12387, 396, 76, 1, 68, "Subsection",
CellTags->"probii"],
Cell[12466, 399, 928, 32, 220, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[13431, 436, 68, 1, 68, "Subsection",
CellTags->"prob3"],
Cell[13502, 439, 933, 16, 395, "Text"]
}, Open ]],
Cell[CellGroupData[{
Cell[14472, 460, 71, 1, 68, "Subsection",
CellTags->"prob4"],
Cell[14546, 463, 590, 11, 220, "Text"]
}, Open ]]
}, Open ]]
}, Open ]]
}
]
*)
(***********************************************************************
End of Mathematica Notebook file.
***********************************************************************)