{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# [01_Variables_Control_Packages](01_Variables_Control_Packages.ipynb)\n", "\n", "Compute $2^{100}$:\n", "- As fast as possible\n", "- As exact as possible" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# This will be wrong:\n", "2^100" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# An exact way to compute the answer is using arbitrary-precision integers\n", "BigInt(2)^100" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Usually faster than arbitrary-precision computations are\n", "# floating point operations\n", "2.0^100" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise\n", "Compute\n", "$$ 15! \\qquad 100! \\qquad \\left(\\begin{array}{c} 100 \\\\ 15 \\end{array}\\right) $$\n", "with the Julia you know so far." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "res15 = 1\n", "for i in 1:15\n", " res15 = res15 * i\n", "end\n", "\n", "res100 = BigInt(1)\n", "for i in 1:100\n", " res100 = res100 * i\n", "end\n", "\n", "res85 = BigInt(1)\n", "for i in 1:(100 - 15)\n", " res85 = res85 * i\n", "end\n", "\n", "res100 / res85 / res15" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# [02_Functions_Types_Dispatch](02_Functions_Types_Dispatch.ipynb)\n", "\n", "### Exercise\n", "- Take timesteps of $\\Delta t = 0.1$ and start at $x_n = 0$ and $v_n = 1$. Propagate the dynamics for 5 steps in the harmonic potential $V = x^2$.\n", "- Now do the same thing using the confining potential:\n", " $$ V_\\alpha(x) = \\left\\{ \\begin{array}{ll} (|x| - 1)^\\alpha & |x| > 1 \\\\ 0 & \\text{else}\\end{array}\\right.$$\n", " for $\\alpha = 4$." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "function acceleration(x)\n", " if x > 1 || x < -1\n", " -4 * sign(x) * (abs(x) - 1)^3\n", " else\n", " 0\n", " end\n", "end\n", "\n", "x = 0.\n", "v = 1.\n", "for i in 1:5\n", " x, v = euler(acceleration, 0.1, x, v)\n", "end\n", "println(x, \" \", v)\n", "\n", "#\n", "# A shorter, but equivalent way using the \"do\" Syntax\n", "#\n", "\n", "x = 0.\n", "v = 1.\n", "for i in 1:5\n", " x, v = euler(0.1, x, v) do x\n", " if x > 1 || x < -1\n", " -4 * sign(x) * (abs(x) - 1)^3\n", " else\n", " 0\n", " end\n", " end\n", "end\n", "println(x, \" \", v)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise\n", "Which of the following type are subtypes of another?\n", "Try to guess first and then verify by using the operator `<:`.\n", "\n", "```julia\n", "Float64 AbstractFloat Integer\n", "Number AbstractArray Complex\n", "Real Any Nothing\n", "```" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The following type chains (subtype ``<:`` supertype) are true:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Float64 <: AbstractFloat <: Real <: Number <: Any" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Integer <: Real <: Number <: Any" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Complex <: Number <: Any" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Nothing <: Any" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "AbstractArray <: Any" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# [03_Arrays](03_Arrays.ipynb)\n", "\n", "### Exercise\n", "Create the following arrays using Julia code:\n", "$$\\left(\\begin{array}{ccccc}\n", " 2&2&2&2&2 \\\\\n", " 2&2&2&2&2 \\\\\n", " 2&2&2&2&2 \\\\\n", " \\end{array}\\right) \\qquad\n", " \\left(\\begin{array}{cccc}\n", " 0.1&0.5&0.9&1.3\\\\\n", " 0.2&0.6&1.0&1.4\\\\\n", " 0.3&0.7&1.1&1.5\\\\\n", " 0.4&0.8&1.2&1.6\\\\\n", " \\end{array}\\right)\n", " $$" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "ones(Int, 3, 5) + ones(Int, 3, 5)\n", "\n", "# or \n", "\n", "2ones(Int, 3, 5)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "reshape(collect(1:16), 4, 4) / 10" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise: Propagating a particle in 2D\n", "We want improve our capabilities to perform simulations to two dimensions using the 2D harmonic potential\n", "$$ V\\left( \\begin{array}{c} x_1 \\\\ x_2 \\end{array}\\right) = x_1^2 + x_2^2 $$\n", "and the acceleration map\n", "$$ \\vec{A}_V = - \\nabla V. $$\n", "The adapted forward-Euler scheme is\n", "$$\\left\\{\\begin{array}{l}\n", "\\vec{v}^{(n+1)} = \\vec{v}^{(n)} + \\vec{A}_V(\\vec{x}^{(n)}) \\Delta t\\\\\n", "\\vec{x}^{(n+1)} = \\vec{x}^{(n)} + \\vec{v}^{(n)} \\Delta t\\\\\n", "\\end{array}\\right. .$$\n", "\n", "Change the `euler` function we introduced in [02_Functions_Types_Dispatch.ipynb](02_Functions_Types_Dispatch.ipynb)\n", "accordingly (Hint: Suprisingly few changes are needed) and run the dynamics for five steps using\n", "$$x_n = \\left( \\begin{array}{c} 0 \\\\ 0 \\end{array}\\right) \\qquad v_n = \\left( \\begin{array}{c} 1 \\\\ 0 \\end{array}\\right) .$$\n", "Check against your previous implementation." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# The previous euler implementation can be used without any additional change:\n", "euler(A, Δt, xn, vn) = xn + vn * Δt, vn + A(xn) * Δt\n", "\n", "# Assuming x to be the vector (x_1 \\\\ x_2), the derivative of V can be written as:\n", "A2d(x) = -2x\n", "\n", "# Therefore propagating 5 steps is equal to:\n", "x = [0, 0]\n", "v = [1, 0]\n", "for i in 1:5\n", " x, v = euler(A2d, 0.1, x, v)\n", "end\n", "println(x)\n", "println(v)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# [05_Dancing_Particles](05_Dancing_Particles.ipynb)\n", "\n", "**Exercise:**\n", "\n", "A much more stable integrator than the `euler` we used so far is the verlocity Verlet:\n", "\n", "$$\\left\\{\\begin{array}{l}\n", "x^{(n+1)} = x^{(n)} + v^{(n)} \\Delta t + \\frac{A_V(x^{(n)})}{2} \\Delta t^2\\\\\n", "v^{(n+1)} = v^{(n)} + \\frac{A_V(x^{(n))} + A_V(x^{(n+1)})}{2} \\Delta t\\\\\n", "\\end{array}\\right. $$\n", "\n", "- Program this algorithm, taking care that it supports multi-dimensional positions and velocities as well. In practice we would like to avoid recomputing $A_V(x)$ as much as possible, since this is usually the expensive step \n", " of the dynamics. For our purposes there is no need to keep an eye on that.\n", "- How does the previous dynamics look like in this example. Does this algorithm conserve energy (phase-space plot)?\n", "- Also look at the Morse potential" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# One way to code the verlet function is:\n", "function verlet(A, Δt, xn, vn)\n", " An = A(xn)\n", " x_next = xn + vn * Δt + An/2 * Δt^2\n", " v_next = vn + (An + A(x_next)) / 2 * Δt\n", " x_next, v_next\n", "end" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise\n", "Program the total potential function for a matrix $\\textbf{x}$. A useful function is `norm` from the `LinearAlgebra` package." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# One solution:\n", "function Vtot(Vpair, x)\n", " n_particles = size(x, 2)\n", " accu = zero(eltype(x)) # Get a zero of the appropriate type\n", " for i in 1:n_particles, j in i+1:n_particles\n", " accu += Vpair(norm(x[:, i] .- x[:, j]))\n", " end\n", " accu\n", "end" ] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.3.0", "language": "julia", "name": "julia-1.3" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.3.0" } }, "nbformat": 4, "nbformat_minor": 2 }