# Assignment 1
Physics 141/241: Computational Physics I

Instructor: Javier Duarte

Spring 2023

Due date for first version: Friday, April 21, 2023 8pm_

Due date for corrected version: Wednesday, April 26, 2023 8pm_

Total points: 110 [141], 150 [241]

## Policies
You are free to collaborate on all the problems, subject to the collaboration policy stated in the syllabus.
You must submit all code used in the homework (if you don't there is a 50% penalty).
You are welcome to use Python, C/C++, or any other language you are comfortable with **within reason**.
Please comment your code such that the TA can follow along and run it without any issues.

The *first version* of the assignment will be graded on effort and completeness and is worth 50% of the total grade.
The *corrected version* of the assignment will be graded on effort and correctness and is worth 50% of the total grade.


## Submission Instructions
**PLEASE NOTE** that there are two steps to submitting your assignment.
Both must be submitted by the deadline.

- Please submit your report as a single `.pdf` file to Gradescope under **Assignment 1 Report** or **Assignment 1 Report Corrections**.
In the report, include any images generated by your code along with your answers to the questions.
  - For instructions specifically pertaining to the Gradescope submission process, see https://www.gradescope.com/get_started#student-submission.
- Please submit your code as a `.zip` archive to Gradescope under **Assignment 1 Code** or **Assignment 1 Code Corrections**. The `.zip` file should contain your code files for all problems. Submit your code either as `.ipynb`, `.py`, or `.c` files (etc.). Follow this naming conventions:
  - For code related to Problem 1(a), include a file named `problem1a.ipynb`, `problem1a.py`, or `problem1a.c` (etc.).


## Introduction
To prepare for your first assignment, you might want to download and modify a sample program.
The source code is available in C ([`leapint.c`](leapint.c)) and in Python ([`leapint.py`](leapint.py)). 
If this link is broken, navigate here: [https://github.com/jmduarte/phys141/blob/main/assignments/assignment1/leapint.py](https://github.com/jmduarte/phys141/blob/main/assignments/assignment1/leapint.py).
To compile the C program, use the command

```bash
gcc -lm -o leapint leapint.c
```

This compilation command creates the executable `leapint`.
The `-lm` links the math libraries during compilation.
To run the compiled program, use the command

```bash
./leapint
```

Alternatively, you can use the python program
```bash
python leapint.py
```

The output consists of four columns, listing (1) time, (2) point index, (3) position `x`, and (4) velocity `v`.
Running the program as supplied yields the trajectory of a point starting at `(X0, V0) = (1, 0)` in the phase space of `(x, v)`, defined by the linear pendulum, or equivalent simple harmonic oscillator, `a(x) = -x`.
The total number of steps taken, number of steps between outputs, and time step used are determined by the parameters `mstep`, `nout`, `dt`, respectively.
These parameters are set in the `main` procedure of the program.

## Problem 1 (141/241) [40 points]
(a) [10 points] Modify the statements which set up initial conditions in the main program to produce trajectories starting from the points `(2, 0)` and `(0, 3)`.
On the `(x, v)` plane, plot these trajectories (e.g. Matplotlib, or your favorite plotting platform) together with the one starting from `(1, 0)`.

(b) [10 points] Replace the linear pendulum in the `accel` routine with the "nonlinear pendulum", `a(x) = -sin(x)`, and again plot trajectories starting from the points `(1, 0)`, `(2, 0)`, and `(0, 3)`. 

(c) [10 points] Animate the phase space motion for the three trajectories simultaneously.

(d) [10 points] Make two-dimensional `(x, v)` plots and three-dimensional `(x, v, t)` plots for (a) and (b).

## Problem 2 (141/241) [20 points]
(a) [10 points] Write leapfrog code for Mercury, Venus, Earth, Mars, Jupiter, Saturn, Uranus, Neptune, Pluto, orbiting around the Sun.
Put the Sun at the origin of your coordinate system.
Do some web-based research for setting up the initial conditions to have realistic scales of the orbits.
Help and links will be provided during lectures and labs.
Ignore the interactions of the planets in the time integration and keep the Sun fixed at the origin.

(b) [10 points] Make a composite plot of all the planetary orbits you calculated.

## Problem 3 (141/241) [30 points]
(a) [10 points] Compare your numerically integrated orbits to the analytic form of the elliptic orbits for Earth and Mars.
Investigate the error behavior over the growing number of revolutions of the planets.
Do you detect a deterioration with the number of revolutions growing?

(b) [10 points] For the two planets, plot the numerically integrated orbit and the analytic elliptic orbit on the same plot.
Each planet should have its own plot.

(c) [10 points] Animate Earth's motion in Matplotlib (or some other application of your choice) over many revolutions.

## Problem 4 (141/241) [20 points]
Prove the exact conservation of angular momentum of the discretized leapfrog-Verlet algorithm on Kepler orbits.

## Problem 5 (241 only) [20 points]
Calculate and plot the orbit of the Halley comet (read about the comet and its orbit in [Wikipedia](https://en.wikipedia.org/wiki/Halley%27s_Comet)).

## Problem 6 (241 only) [20 points]
Put the Halley comet and Earth in the same animation and demonstrate one full Halley year.
See [Wikipedia](https://en.wikipedia.org/wiki/Halley%27s_Comet) and the [comet table](https://en.wikipedia.org/wiki/List_of_numbered_comets) for more information.

## Additional resources
- [Tutorial for Leapfrog Integration](leapfrog.pdf)
- [Leapfrog Python tutorial](http://cvarin.github.io/CSci-Survival-Guide/leapfrog.html)
- [Kepler orbit tutorial](kepler.pdf)
- [The Kepler problem](kepler_py.pdf)
- [Geometric numerical integration illustrated by the Stoermer–Verlet method](stoermer-verlet.pdf)
- [Kepler Problem](http://www.gnm.cl/fgonzalez/pmwiki.php/Programas/KeplerOrbits): Note this code is using Runge-Kutta integration, but your code should use second-order leapfrog integration.