02 - Particles.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "AMUSE tutorial on particle sets\n",
    "====================\n",
    "\n",
    "\n",
    "AMUSE particle sets are a handy tool for storing data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture _preamble\n",
    "!pip show amuse-framework || pip install amuse-framework\n",
    "# Load the amuse units module\n",
    "from amuse.units import units\n",
    "from amuse.lab import Particles"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Declare a single particle\n",
    "sun_and_earth = Particles(2)\n",
    "sun = sun_and_earth[0]\n",
    "sun.mass = 1 | units.MSun\n",
    "sun.position = (0, 0, 0) | units.au\n",
    "sun.velocity = (0, 0, 0) | units.kms\n",
    "print(f\"Sun={sun}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Now declare the Earth\n",
    "from amuse.units.constants import G\n",
    "\n",
    "def relative_orbital_velocity(mass, distance):\n",
    "    return (G * mass / distance).sqrt()\n",
    "\n",
    "earth = sun_and_earth[1]\n",
    "earth.mass = 1 | units.MEarth\n",
    "earth.position = (1, 0, 0) | units.au\n",
    "\n",
    "vorb = relative_orbital_velocity(\n",
    "    sun_and_earth.total_mass(), \n",
    "    earth.position.length()\n",
    ")\n",
    "earth.velocity = (0, 1, 0) * vorb\n",
    "print(f\"Earth={earth}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sun_and_earth.move_to_center()\n",
    "print(sun_and_earth)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, let's give the particles a specific name (or other attribute)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "setattr(sun_and_earth, \"name\", \"\")\n",
    "sun_and_earth.name = [\"sun\", \"earth\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "How we have declared the particles and moved them to the center of mass.\n",
    "We can also search for a specific particle.\n",
    "For example, the one with the \"earth\" in the attribute \"name\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "earth = sun_and_earth[sun_and_earth.name==\"earth\"]\n",
    "print(f\"Earth={earth}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can add a moon in orbit around the earth"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "moon = Particles(1)\n",
    "moon.name = \"moon\"\n",
    "moon.mass = 7.34767309e+22 | units.kg\n",
    "moon.position = (384400, 0, 0) | units.km\n",
    "vorb = relative_orbital_velocity(\n",
    "    earth.mass + moon.mass, \n",
    "    moon.position.length()\n",
    ")\n",
    "moon.velocity = (0, 1, 0) * vorb\n",
    "print(f\"moon={moon}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Moon, however, is now somewhere inside the Sun with zero velocity - which is not good.\n",
    "We will have to relocate the Moon to make it orbit around the Earth.\n",
    "We do that by simply adding the positions and velocity of Earth to the Moon's."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "moon.position += earth.position\n",
    "moon.velocity += earth.velocity"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can add the moon to the Sun and Earth system"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sun_and_earth.add_particle(moon)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note the use of the singular here, because we only add a single particle to the particle set sun_and_earth.\n",
    "It is probably better to rename the sun_and_earth now."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solar_system = sun_and_earth"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is important now to recenter the entire system, because by adding the moon we shifted the center of mass."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solar_system.move_to_center()\n",
    "print(f\"Solar system: {solar_system}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now manipulate the planetary system, or query it.\n",
    "for example by querying the masses."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"mass=\", solar_system.mass.in_(units.MEarth))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This gives us a list of the masses of all objects, in units the the earth's mass. In fact, each of the particle's attributes is a simple numpy array: it can be assigned and manipulated as such.\n",
    "\n",
    "Another way to ecquire the same information could be done as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"mass={solar_system.mass/solar_system[1].mass}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In some (hopefully rare) cases you may want to use the particle set or its attributes as simple numpy arrays, without the units.\n",
    "This is easily achieved by stripping the unit from the array. This can be realied by explicitely querying the selected parameter with that specific unit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solar_system.position.value_in(units.parsec)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "which, in this case, gives you a 2-dimensional *numpy.array* of the positions of star, planet and moon in units of a parsec."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, you may want to query the particle set solar_system.\n",
    "For example by asking what attributes it has.\n",
    "This can be done as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dir(solar_system)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Or get some general help on the underlying particle class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(solar_system)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You have performed some rudimentary operations on a particle set.\n",
    "It is now time to experiment a little for yourself.\n",
    "\n",
    "\n",
    "Assignments and questions:\n",
    "---------------\n",
    "\n",
    "### Assignment 1:\n",
    "Add the planet Jupiter (see [Wikipedia](https://en.wikipedia.org/wiki/Jupiter)) to your small planetary system.\n",
    "\n",
    "### Assignment 2:\n",
    "Your planetary system is notoriously planar, and initialy the earth and moon are positioned along the Cartesian x-axis with the velocity vector in the Cartesian y-direction. \n",
    "\n",
    "Make the Sun-Earth-Moon system more realistic by introducing a small inclination to the Earth's and Moon's orbits and by giving them a random mean anomaly.\n",
    "\n",
    "The [Orbital element module](https://github.com/amusecode/amuse/blob/main/src/amuse/ext/orbital_elements.py) of AMUSE could come in handy. \n",
    "\n",
    "### Assignment 3:\n",
    "Calculate the total  gravitational binding energy of solarsystem.\n",
    "\n",
    "Now displace the entire particle set by 100 parsec and give it a linear velocity of 100km/s in the z-direction.\n",
    "Then calculate the binding energy of the system again.\n",
    "\n",
    "Did the binding energy of the Solar system change by this translation?\n",
    "\n",
    "### Question 1:\n",
    "Particle sets have the attribute *get_binaries()*.\n",
    "If you use this function to check the binaries in your system you will find that (without Jupiter) you have 3 binaries. Explain why the Sun is in a binary with the Moon. You may want to take a look at the [source code](https://github.com/amusecode/amuse/blob/main/src/amuse/datamodel/particle_attributes.py).\n",
    "\n",
    "### Assignment 4:\n",
    "Generate another particle set with a 2 solar-mass star and two planets of 10 and 100 Earth masses in circular orbits at 0.1 and 0.6 au.\n",
    "Place this second planetary system at apocenter around your Solar system (true anomaly of 180 degrees) at a semimajor axis of 60 au with an eccentricity of 0.6.\n",
    "Then move the entire system to the center of mass.\n",
    "\n",
    "### Question 2:\n",
    "Which of the orbits of the binary star with planets from *Assignment 4* has the highest binding energy?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}