Mathematica solving problems
����������� In this assignment, we will be learning how to navigate the Wolfram language. The best software pack- age to begin on is Mathematica. One can access Mathematica in two ways:
1) Go to the dungeons (basement) of AP&M, and any computer lab should have Mathematica already installed. Or for a complete list of computer labs with Mathematica installed, check out here: http://insci14.ucsd.edu/acs_sql/scripts/lablookup/ And go to Specialty Software → Mathematica → All
2) Follow the instructions on this page: https://acms.ucsd.edu/services/software/available-software/mathematica.html to install Mathematica on your personal machine (for free! [or well possibly included as part of your tuition money, so you might as well use it]). Some students have already had some issues with the method, so post on piazza if you get stuck!
����������������������������� The first thing we need to be able to do is have Wolfram remember things for us. This is done by storing information into variables. Simply type a name you want your variable to be, put an equals sign next to it, and then whatever information you want to store. Note this equals sign is actually an “assignment” operator, since it assigns the right hand side to the left hand side. After you’re done typing, execute the cell by pressing Shift + Enter (lines and groups of lines in wolfram come together in “cells”, and you can click on different cells with the right brackets appearing on the righthand side of the notebook.)
� = �
�
� = �
�
� + �
�
Note you can always re-execute a cell in the Wolfram language by clicking onto that cell and hitting shift + enter. Further note that you can prevent the output by appending a semi-colon to the end of a line of code. You can also put multiple lines of code into a single cell. Compare:
� = ��
� = ��
� / �
��
��
�
�
and:
� = ���
� = ���
� / �
�
�
These very basic calculations are actually extremely helpful, see the following: For example, the force between two points charges in physics:
�� = �� * ��-��
�� = �� * ��-��
� = � * ��-��
ϵ� = ���� * ��-���
�
� π ϵ�
�� ��
��
�������
Now I don’t have to worry about punching a wrong number into a calculator, and having to redo the whole calculation! I can also make sure the formula is correct, and edit it without having to redo the whole calculation!
Notice how I used an “epsilon naught”. Special characters are easy to type into the Wolfram language, just hit Escape + letter + Escape to give the greek equivalent to that letter.
escape + a + escape = α escape + b + escape = β escape + c + escape = χ
Woops, sometimes you have to just hit the first letter that sounds like the greek letter. We were looking for gamma, so let’s hit escape g escape = γ. Wolfram will likely give you suggestions for characters, I personally find those annoying after while, so I turn it off. Capitals:
escape + G + escape = Γ escape + X + escape = Ξ escape + O + escape = Ω
Also notice how the numbers I typed previously had exponents for scientific notation. And furthermore the formula had fractions which looked quite similar to how we would write it by hand. Mathematica
2��� XC_Assignment.nb
excels at making formulas look “nice”. To type an exponent: either use the carrot: ^, or type (at least on a mac): control 6
����
����
���
����
To type in fractions, use the forward slash: /, or type (on a mac): control /
�� / �
��
�
��
�
��
�
There are many ways to beautify your formulas, often we will just learn them as we need/want them.
That’s it for variables, let’s now learn about functions. To define a function, say f(x) = x2, use the follow- ing syntax:
�[�_] = ���
In wolfram language, the arguments of functions go into “square brackets”: [ ]. For initially defining functions you put an underscore after the variable in the bracket. Then use assignment operator, =. Then type your function.
You can then later on call your function with an argument, as follows: A generic argument:
�[�]
��
�[��������������������]
���������������������
A numeric argument:
�[�]
��
�[-�]
�
Wolfram has a huge library of built-in functions: For example, to plot something, just type Plot. Note that all built-in Wolfram functions begin with a Capital Letter.
XC_Assignment.nb ���3
����[�[�]� {�� -�� �}]
-1.0 -0.5 0.5 1.0
0.2
0.4
0.6
0.8
1.0
Finding Wolfram functions to do what you want can sometimes be a pain, if you have a hunch for what you want to do, search the documentation. To open the documentation, go to Help → Wolfram Documen- tation, or (on a mac) command shift f, or just hit f1. Personally, I type a word, highlight it, and hit com- mand shift f.
The documentation for the Wolfram language is quite impressive compared to some languages. Here we can learn how to use the function, and many examples are provided. In programming, learning by example is one of the best ways to learn.
One very useful function is the integrate function:
���������[���[�]� {�� �� �}]
��������� ������ � ������� ����������� �������� �� ����� (�) ���{��� �� �}�
�
�
���[��] ⅆ��
Here wolfram remembers that I defined x=10 earlier, and it’s giving me an error. Sometimes if you want to reuse variables that were previously defined, you need to clear them:
�����[�� �]
Now x and y have been given a fresh new start:
���������[���[�]� {�� �� �}]
-� + ���[�]
If you want to convert this into a decimal expression, use the function N or NIntegrate:
�[-� + ���[�]]
��������
����������[���[�]� {�� �� �}]
��������
NIntegrate has the advantage that it can integrate functions that might not have anti-derivatives:
4��� XC_Assignment.nb
����������� ��
� {�� �� �}
�
�
�� ��
ⅆ�
Wolfram can’t figure it out the antiderivative so it spit it back out at us, but Wolfram can compute it numerically:
������������ ��
� {�� �� �}
�������
You may have noticed the integral sign in the above outputs. We can actually perform integrals in Wolfram in a beautiful way with escape characters: use escape int escape to make an integral sign, control - to do the lower limit, and while your cursor is on the lower limit: control 5 for the upper limit, and finally use escape dd escape to give it the differential.
θ�
θ�
���[θ] ⅆθ
���[θ�] - ���[θ�]
Derivatives, can by done with the escape pd escape, pd for partial derivative, and the subscript with ctrl -
∂� � �
� �
∂� (���[�] ���[�])
���[�]
�
Since we are learning about calculus in multiple dimensions, let’s move onto vectors. Vectors in Wol- fram are called Lists, and written with “curly brackets” or simply named, “braces” : { }.
� = {�� �}�
� = {-�� �}�
Suppose we want to make a line that goes through these two vectors. We can write parametric equa- tions for this straight forwardly with a convex combination:
�����[�]�
�[�_] = � � + (� - �) �
{-� + � �� � + �}
Given parametric equations, we can visual the line with a *parametric plot*:
XC_Assignment.nb ���5
��������������[�[�]� {�� -��� ��}]
-20 -10 10 20
-10
-5
5
10
Notice that we had to specify the range of t. The full line includes all points t ∈ ℝ, but typically, we only care about visualizing a subset of values for t.
When plotting it is useful to be able to have control over the plot range, we can encode exactly the plot range by adding an option to ParametricPlot[], (check the help documentation for other fun options):
�������� = ��������������[�[�]� {�� -��� ��}� ��������� → {{-�� �}� {-�� �}}]
-3 -2 -1 1 2 3
-3
-2
-1
1
2
3
Let’s verify that the line is going through the points we specified earlier, A and B. To visual vectors as arrows in Mathematica, we need to refer to the graphics library for arrows. You can read the help documentation for more info:
6��� XC_Assignment.nb
���� = ��������[�����[{{�� �}� �}]]
���� = ��������[�����[{{�� �}� �}]]
The above makes arrows that go from the origin to the point A or B. Now it would be great if we could show everything on the same plot. To do this use the function called Show[].
XC_Assignment.nb ���7
����[��������� ����� ����]
-3 -2 -1 1 2 3
-3
-2
-1
1
2
3
Great! Now let me explain a thing about convex combinations. Convex combinations are sums of
objects where the coefficients themselves add up to one, for example here’s convex combination: α A →
+
β B →
, with α + β =1. The constraint α + β =1 can be manipulated to eliminate β: β = 1-α, so we can re-
write the convex combination without β: α A →
+ (1-α) B →
.
Can you think of what the significance of coefficients adding up to one is? Hint: probabilities add up to one! So in fact, a convex combination is just a *weighted average*. Of course, probabilities can only be positive. Watch what happens when we restrict the coefficients t and (1-t) in the previous plot to be positive:
8��� XC_Assignment.nb
�������� = ��������������[� � + (� - �) �� {�� -��� ��}� ��������� → {{-�� �}� {-�� �}}�
�������������� → ��������[{�� �� �}� � > � �� (� - �) > �]]�
����[��������� ����� ����]
-3 -2 -1 1 2 3
-3
-2
-1
1
2
3
The constraints, t > 0 && (1-t) > 0, are equivalent to saying t ∈ (0,1)
So the line is really just a weighted average between two points. But this only gives us the *interpolation* between the points. We can *extrapolate* by letting the coefficients go negative, i.e., lifting the constraint that t∈(0,1) to t∈ℝ.
XC_Assignment.nb ���9
�������� =
��������������[� � + (� - �) �� {�� -��� ��}� ��������� → {{-�� �}� {-�� �}}]�
����[��������� ����� ����]
-3 -2 -1 1 2 3
-3
-2
-1
1
2
3
We can convert parametric equations into implicit equations by eliminating the parameter(s). In general parameter elimination is not easy, but there are algorithms to do it for us. Just as you don’t need to know how an engine works in order to drive a car, you don’t necessarily need to know how to parameter elimination works to make a computer do it for you. The function to do parameter elimination is called GroebnerBasis[], and it is used in the following way: GroebnerBasis[ equations, variables to keep, variables to eliminate]. (Note: in a class on linear algebra, you learn how to eliminate *linear* parameters, but Groebner basis will eliminate *non-linear* parameters - and you might only learn about this in a special topics course! )
�������������[� � + (� - �) � ⩵ {�� �}� {�� �}� {�}]
{� + � - � �}
Note that when typing in entire equations as an input, you must use *DOUBLE EQUALS*, ==, because a single equals, =, in computer science is not actually an equality, but an assignment.
Anyway, as you can see, we wanted to keep variables {x,y} and eliminate the parameter {t}. Note, we usually assume the output of the GroebnerBasis command is equal to zero.
Remember: to plot implicit equations, we use a ContourPlot!
10��� XC_Assignment.nb
���� = �����������[� + � - � � ⩵ �� {�� -�� �}� {�� -�� �}]�
����[����� ����� ����]
-3 -2 -1 0 1 2 3
-3
-2
-1
0
1
2
3
Same exact thing as before (well the axes have changed a little, but the line is the same)
The same exact discussion works for planes in 3-dimensions of course: (however, we have to use the 3D versions of the above functions, e.g., ParametricPlot3D instead of ParametricPlot)
XC_Assignment.nb ���11
� = {�� �� �}�
� = {�� �� �}�
� = {�� �� �}�
���� = ����������[�����[����[{{�� �� �}� �}]]]�
���� = ����������[�����[����[{{�� �� �}� �}]]]�
���� = ����������[�����[����[{{�� �� �}� �}]]]�
��������� = ����������������[� � + � � + (� - � - �) ��
{�� -�� �}� {�� -�� �}� ��������� → {{-�� �}� {-�� �}� {-�� �}}�
����� → ������ ���������� → {�� �� �}� ���� → �����]�
����[���������� ����� ����� ����]
Note, you can click on the graph to move around the camera angle, are you convinced that the plane is going through the 3 points?
Notice what happens when we restrict the values of coefficients s, t, (1-s-t) to be *positive*:
12��� XC_Assignment.nb
��������� = ����������������[� � + � � + (� - � - �) ��
{�� -�� �}� {�� -�� �}� ��������� → {{-�� �}� {-�� �}� {-�� �}}�
����� → ������ ���������� → {�� �� �}� ���� → ������
�������������� → ��������[{�� �� �� �� �}� � > � �� � > � �� (� - � - �) > �]]�
����[���������� ����� ����� ����]
The positive coefficient convex combination gives us an *interpolation*. We can *extrapolate* by letting the coefficients be positive or negative, to get an entire plane
For fun: note what happens when we don’t use a convex combination:, say let’s change (1 - s - t) to
1 - s2 - t):
XC_Assignment.nb ���13
��������� = ����������������� � + � � + � - �� - � ��
{�� -�� �}� {�� -�� �}� ��������� → {{-�� �}� {-�� �}� {-�� �}}�
����� → ������ ���������� → {�� �� �}� ���� → ������
����[���������� ����� ����� ����]
Not even close to a plane. If the coefficient don’t add up to one, you don’t get the linear interpolation between the points! And if you don’t get the interpolation correct, you certainly aren’t going to get the extrapolation correct either. Now, instead, let’s change (1 - s - t) to (1 - t):
14��� XC_Assignment.nb
��������� = ����������������[� � + � � + (� - �) ��
{�� -�� �}� {�� -�� �}� ��������� → {{-�� �}� {-�� �}� {-�� �}}�
����� → ������ ���������� → {�� �� �}� ���� → �����]�
����[���������� ����� ����� ����]
The surface is still linear, but it doesn’t interpolate between the points we wanted.
Anyway, moving on, computing the implicit equation for the plane is easy with a computer, no geometri- cal considerations necessary!:
�������������[� � + � � + (� - � - �) � == {�� �� �}� {�� �� �}� {�� �}]
{-�� + � � + � � + � �}
Again, implicit equations can be visualized with *contour plots*, in this case, ContourPlot3D:
XC_Assignment.nb ���15
��������� = �������������[-�� + � � + � � + � � ⩵ �� {�� -�� �}� {�� -�� �}�
{�� -�� �}� ����� → ������ ���������� → {�� �� �}� ���� → �����]�
����[���������� ����� ����� ����]
Same plot as before!
Okay, let’s take a moment to let our hippocampus store this into long-term memory:
Parametric equations → parametric plot Implicit equations → contour plot But both give the same thing — which method to use depends on what you’re given.
Let’s do a more sophisticated example now, here’s parametric equations for *some* surface in 3D.
�[�_� �_] = � � � - �� - �� � � � � - �� - �� � � � ��
16��� XC_Assignment.nb
�����������������[�� �]� {�� -�� �}� {�� -�� �}�
�������������� → ��������{�� �� �� �� �}� � - �� - �� > ��
���� → ������ ���������� → ��� ����� → ������ ���� → �����
Note that I added the RegionFunction to avoid taking square roots of negative numbers. It’s kind of a pretty surface in a way...
Eliminating the parameters to get the implicit equation seems a little daunting if we had to do it by hand. Thankfully we now know how to drive the Mathematica Mobile, and GrobnerBasis[] will eliminate the parameters for us (please note, Groebner basis only eliminates parameters from algebraic expressions, it cannot always help you if you use transcendental equations by including trigs or logarithms ):
�������������� � � - �� - �� � � � � - �� - �� � � � � ⩵ {�� �� �}� {�� �� �}� {�� �}
�� �� - � � � � + �� �� + �� ��
Nice! An implicit equation, but does it really correspond to the same surface as before?
XC_Assignment.nb ���17
��������������� �� - � � � � + �� �� + �� �� ⩵ �� {�� -�� �}� {�� -�� �}� {�� -�� �}�
���� → ������ ���������� → ��� ����� → ������ ���� → ������ ���������� → ��
Phew, it turns out I wasn’t lying about anything! (You might notice the quality of the plot is slightly worse in the contour plot compared to the parametric plot - this is because contour plots require more computa- tion. You can increase the quality by increasing the PlotPoints option, but don’t go too high or it’ll take a really long time! )
���������� Save your current Mathematica file, called a notebook, into some directory on the machine you’re working on. Once you do that, type the command:
������������[�����������������[]]
This is where the plots will go to (that we’ll export later on)
���������
The following parametric equations
r → (t) = ±
1-t2
1+t2 , t 1-t2
1+t2 for t∈ (-1,1),
have a parametric plot as the following ( I need two plots for each of the ± sign):
18��� XC_Assignment.nb
����� = �������������� � - ��
� + �� �
� � - ��
� + �� � {�� -�� �}�
����� = ��������������- � - ��
� + �� �
� � - ��
� + �� � {�� -�� �}�
����� = ���������� ������ ��������� → {{-�� �}� {-�� �}}�
��������� → �� →
(�) = ± � - ��
� + �� �
� � - ��
� + �� �
-1.0 -0.5 0.5 1.0
-1.0
-0.5
0.5
1.0
r → (t) = ±
1 - t2
1 + t2 , t 1 - t2
1 + t2
(By the way, square roots can by typed into Mathematica on a mac with control shift 2, or just type Sqrt[] ) This curve is called a lemniscate. Find the implicit equations for this lemniscate, and create a contour plot, and label the plot with the implicit equation you find. Export the plot in the following way:
������[��������������� �����]
������������
but be sure to change “plot3” to the whatever variable name you store the ContourPlot into! Upload the image to TritonEd when you complete the assignment.
���������
Recall from physics that Newton’s law is an equation for unknown parametric equations which describe the position of an object. Let’s use Mathematica to solve Newton’s laws for the parametric equations
XC_Assignment.nb ���19
that describe position of an electrically charged particle (e.g., an electron) in an external electric field in 2-dimensions, and then *animate* the motion of the particle. Here’s Newton’s law in this context:
d
dt2 r → (t) =
q
m E →
where r → (t) is the (parametric) position of the particle, E
→
is an external electric field, q is the charge of the particle, m is the mass of the particle. Don’t worry if you’ve never heard of these, the important thing is that we are solving for parametric equations, just think of everything else as a given function.
Let’s make the external electric field with some randomness so everyone’s plot will be different.
{��� ��} = ����������[{-�� �}� {�}]�
{��� ��} = ����������[{-�� �}� {�}]�
{��� ��} = ����������[{-�� �}� {�}]�
����[�_� �_] �=
-�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��}
+ -�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��}
+ -�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��}
Physically, this corresponds to 3 *fixed* attractive charged particles at positions (x1,y1), (x2,y2), and (x3,y3). Our free particle is going to wiggle and move around them! The electric field is a vector-valued function, and we can visualize vector valued functions with a VectorPlot[].
20��� XC_Assignment.nb
����������[����[�� �]� {�� -�� �}� {�� -�� �}]
-1.0 -0.5 0.0 0.5 1.0
-1.0
-0.5
0.0
0.5
1.0
Note that the electric field gets weak *very fast* with distance, so the arrows quickly shrink and disap- pear. We can actually force all the arrows to be the same size:
�������� =
����������[����[�� �]� {�� -�� �}� {�� -�� �}� ����������� → {����� �� ����}]
-1.0 -0.5 0.0 0.5 1.0
-1.0
-0.5
0.0
0.5
1.0
Note that a VectorPlot of a vector valued function E →
(x, y) will simply put the tail of the vector E →
(x0, y0),
XC_Assignment.nb ���21
at the point (x0,y0).
Okay, let’s solve Newton’s laws! First we need to pick an initial position and initial velocity of our parti- cle. Let’s let the position be random, and the initial velocity be zero, and let’s compute the position for an initial parameter ti=0 to the final parameter value tf=1, so that t∈(ti, tf).
{��� ��} = ����������[{-�� �}� {�}]�
{���� ���} = {�� �}�
{��� ��} = {�� �}�
Newton’s law is a differential equation, we can solve a differential equation with the command NDSolve[]:
��� = �������[{{���[�]� ���[�]} ⩵ �� ����[�[�]� �[�]]� �[�] ⩵ ���
�[�] ⩵ ��� ��[�] ⩵ ���� ��[�] ⩵ ���}� {�� �}� {�� ��� ��}][[�]]�
{����� ����} = {�� �} /� ����
Note that the solution, i.e., the parametric equations, in general do not have a closed form expression. The answers are instead stored as numerical data. So you’ll just have to accept the syntax of the above code for now because that’s just how you manipulate numerical data in this case. Also note that I set q=50 and m=1, (there’s a 50 as a coefficient of the electric field), you can modify that number if you want.
������� = ��������������[{����[�]� ����[�]}� {�� �� �}� ��������� → ���]�
����[��������� �������]
-1.0 -0.5 0.0 0.5 1.0
-1.0
-0.5
0.0
0.5
1.0
Er, wow, no wonder physics is so hard, that’s the path a charged particle follows?! (Notice how the particle tries to go in the direction of the arrows but its big fat inertia makes it hard to turn in time!) How the heck are we supposed to do *any* physics without a computer?!
Since we are using some randomness, occasionally we would like to run our code again quickly to regenerate the plot using new numbers. We might as well put all the necessary code into a single cell so we only have to run one cell (if you do put everything into one cell, MAKE SURE EVERY INDIVID-
22��� XC_Assignment.nb
UAL PIECE OF THE CODE WORKS BEFORE YOU SHOVE IT ALL INTO A SINGLE CELL, OR ELSE IT IS A NIGHTMARE TO DEBUG):
{��� ��} = ����������[{-�� �}� {�}]�
{��� ��} = ����������[{-�� �}� {�}]�
{��� ��} = ����������[{-�� �}� {�}]�
{��� ��} = ����������[{-�� �}� {�}]�
{���� ���} = {�� �}�
{��� ��} = {�� �}�
����[�_� �_] �=
-�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��}
+ -�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��}
+ -�
(� - ��)� + (� - ��)� �/�
{� - ��� � - ��} �
��� = �������[{{���[�]� ���[�]} ⩵ �� ����[�[�]� �[�]]� �[�] ⩵ ���
�[�] ⩵ ��� ��[�] ⩵ ���� ��[�] ⩵ ���}� {�� �}� {�� ��� ��}][[�]]�
{����� ����} = {�� �} /� ����
������� = ��������������[{����[�]� ����[�]}� {�� �� ��}� ��������� → ���]�
�������� =
����������[����[�� �]� {�� -�� �}� {�� -�� �}� ����������� → {����� �� ����}]�
����[��������� �������]
-1.0 -0.5 0.0 0.5 1.0
-1.0
-0.5
0.0
0.5
1.0
You could just run the above cell now until you get an output you are happy with.
XC_Assignment.nb ���23
If you get a warning message like:
��������������� ������� � ���� �
�� � == ������������������� � ���� ���� �� ����������� ����� ����������� �� ����� ������ ��������� � �
� � ���� � �� � -������������������� � �����
Just set, tf = 0.5 for the remainder of the problem, or just try running the code again until no errors appear. (Sometimes two charges will collide and the computer gets mad at us)
Now let’s animate it. An animation just plays a bunch of frames one by one with slight variations. So let’s make a collection of frames. We can do so with the Table[] command.
������ = ���
����� = �����
������� = ��������������[{����[�]� ����[�]}� {�� ��� ��}� ��������� → ���]�
����[��������� �������]�
��� �� + ����� ��� �� - ��
������ �
This will generate frames that each will plot t ∈(ti , tp) for increasing values of tp until tf, to give about 60 total frames. You can increase the number of frames to get a smoother plot, but it takes longer to compute. Now export as an animated gif:
������[��������������� �����]
������������
You can view the animation if you open the saved .gif file in a webbrowser, such as firefox/tor/safari You can also save as other video formats, but I want you to upload a .gif to TritonEd as your finished work.
You can also view the animation within Mathematica:
24��� XC_Assignment.nb
�����������[�����]
-1.0 -0.5 0.0 0.5 1.0
-1.0
-0.5
0.0
0.5
1.0
Note that only in *exceedingly* rare and *exceedingly* special circumstances can you turn the paramet- ric equations from Newton’s law into implicit equations. One special case is two planets that orbit each other - in which the implicit equations correspond to an ellipse, cool huh? But in general, don’t even think about finding the implicit equations, since after all the parameter now gets the meaning of time, so why would you want to eliminate that?
������� Please upload your two .gif files to TritonEd to receive credit
XC_Assignment.nb ���25