{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Other features worth mentioning\n",
    "\n",
    "In this and the next chapter on [Useful Packages](09_Useful_Packages.ipynb) I want to give an overview of the functionality and scope Julia can deliver. Details you can find in the referenced links. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Packages we need here:\n",
    "import Pkg; Pkg.add([\"PyCall\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dictionaries and tuples\n",
    "\n",
    "There are plenty of Julia datatypes I have not yet mentioned. This includes dictionaries and tuples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Construct a dictionary, note that dictionaries are typed as well.\n",
    "a = Dict(\"a\" => 1, \"b\" => 2)\n",
    "println(typeof(a))\n",
    "\n",
    "for (l, i) in a  # Loops are over key-value pairs\n",
    "    println(\"key $l maps to $i\")\n",
    "end\n",
    "\n",
    "# Does not work:\n",
    "a[2] = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# A tuple\n",
    "t = (2, \"abc\", )\n",
    "@show t\n",
    "\n",
    "# A named tuple\n",
    "nt = (value=2, description=\"abc\")\n",
    "@show nt\n",
    "\n",
    "println(\"The description is \", nt.description)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/base/collections/\n",
    "- https://docs.julialang.org/en/v1/base/base/#Core.NamedTuple"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Custom types and structs\n",
    "\n",
    "Apart from the built-in types, Julia also offers to declare custom types. Both abstract types and concrete types can be defined:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "abstract type MyAbstract end\n",
    "\n",
    "struct MyPoint <: MyAbstract\n",
    "    x::Float64\n",
    "    y::Float64\n",
    "    description      # Note: Fields do not need to be typed, in this case they are Any\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\"Class methods\" or \"member functions\" as such do not exist in Julia. Instead one uses free functions with type-specific dispatch, e.g."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import Base: length\n",
    "length(p::MyPoint) = sqrt(p.x^2 + p.y^2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "shift_in_x(p::MyPoint, v) = MyPoint(p.x + v, p.y, p.description)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "p = MyPoint(3, 4, \"bla\")\n",
    "length(p)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compared to other languages such structs are used *a lot less* and much functionality in the standard library and third party packages avoids them completely. Still they can be very useful to combine multiple enties, which belong together."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/manual/types/\n",
    "- https://docs.julialang.org/en/v1/manual/constructors/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interacting with the filesystem\n",
    "\n",
    "Julia comes with a built-in library for file system interaction. This makes a lot of use for the `do` syntax we saw before. For example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cd(homedir()) do\n",
    "    println(\"All files and directories in your home are:\")\n",
    "    println(readdir())\n",
    "end\n",
    "# Back to the old directory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Open a file for writing:\n",
    "open(\"testfile\", \"w\") do f\n",
    "    write(f, \"test\")\n",
    "end\n",
    "\n",
    "# Open a file for reading:\n",
    "open(\"testfile\") do f\n",
    "    println(readline(f))\n",
    "end\n",
    "\n",
    "rm(\"testfile\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/base/file/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Unit testing\n",
    "\n",
    "Some functionality for unit testing is part of Julia's standard library. As of now it is pretty basic, but it get's the job done:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using Test\n",
    "@testset \"My first test\" begin\n",
    "    @test 1 == 1  # Should be good, no?\n",
    "    @test_throws BoundsError (1, 2, 3)[5]\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/stdlib/Test/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Macros and metaprogramming\n",
    "\n",
    "Julia has pretty strong metaprogramming support. In fact macros are even an integral part of the core language and their use is frequent in the standard library itself. Common use cases include a temporary change of behavior, to automate the generation of repetitive code or to hide some implementation details under a more convenient syntax. In this course we have already met quite a few macros, for example `@show`, `@time`, `@btime`, `@.`, `@gif`, .... It is extremely uncommon that one needs to add or modify macros for daily work with Julia. Still if you ever come across the need to implement a macro, their syntax is quite readable:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "macro make_addfunction(value)\n",
    "    :(# Start an expression, which is evaluated as code\n",
    "        # $(Symbol(:add, value)) makes a code symbol by concatenating\n",
    "        # two others together\n",
    "        function $(Symbol(:add, value))(x)\n",
    "            x + $value\n",
    "        end\n",
    "    )\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@make_addfunction 10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "add10(4.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/manual/metaprogramming/index.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Checking up on the compiler (@code_warntype and friends)\n",
    "\n",
    "As discussed before type-stability is important in performance-critical parts of the code such that Julia can infer concrete types before JIT compilation and be very specific on the emitted machine code.\n",
    "\n",
    "But introducing type instabilities can come easier as you might think. Consider the following unsuspicious piece of code for example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compute the sum of a column\n",
    "function sumcol(a, n)\n",
    "    @assert ndims(a) == 2\n",
    "    accu = 0\n",
    "    # @inbounds disables bounds check for a[i, n]\n",
    "    @inbounds for i in 1:size(a, 1)  \n",
    "        accu += a[i, n]\n",
    "    end\n",
    "    accu\n",
    "end\n",
    "\n",
    "function average_sumcol(a)\n",
    "    @assert ndims(a) == 2\n",
    "    res = 0\n",
    "    for j in 1:size(a, 2)\n",
    "        res += sumcol(a, j)\n",
    "    end\n",
    "    res / size(a, 2)\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using BenchmarkTools\n",
    "a = randn(ComplexF64, 500, 500)\n",
    "@btime average_sumcol(a);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's reasonably fast, but it actually contains subtle type instabilities. This can be found using `@code_warntype`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@code_warntype average_sumcol(a);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is a rather complicated output, but it roughly speaking tells us with the bold print, that somehow the compiler cannot decide between `Int64` and `Complex{Float64}` for some lines of code and not between `Float64` and `Complex{Float64}` in others. The variables blocks tells us that the `res`  variable declared in line 14 might be a problem. And indeed the issue is that we declare it as `= 0`, which makes it an integer. Therefore the return type of `average_sumcol` is not stable, because it is different for empty arrays (`Int64`) and non-empty arrays (`Complex{Float64}`). A similar spot is line 4. Replacing both instances with `zero(eltype(a))` makes the code completely type stable (and actually a little faster).\n",
    "\n",
    "In this example one might still be able to read and understand the output, but especially for larger functions `@code_warntype` can become long and unclear. Of help is a package with the unrememberable name `Cthulhu` and its macro `@descend_code_warntype`. Unfortunately it does not work in Jupyter notebooks at the moment, so if you want to see it, try it in the REPL."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/stdlib/InteractiveUtils/#InteractiveUtils.@code_warntype\n",
    "- https://github.com/JuliaDebug/Cthulhu.jl\n",
    "- https://docs.julialang.org/en/v1/manual/performance-tips/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Multithreading\n",
    "\n",
    "In this section, we want to very briefly explore one (of many) ways to do multithreading in Julia. For this we first need to increase the number of threads the Julia kernel uses. We check on the CPUs of this machine"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Sys.cpu_summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In my case, I can set the number of threads to `4`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Threads.nthreads() = 4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us recall the computation of the total potential for the multidimensional case in [Dancing_Particles](05_Dancing_Particles.ipynb):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using LinearAlgebra\n",
    "Vmorse(r; re=0.7, α=1.3, D=10) = D*(1 - exp(-α * (r - re)))^2 - D\n",
    "\n",
    "function Vtot(Vpair, x)\n",
    "    n_particles = size(x, 2)\n",
    "    accu = eltype(x)(0)   # Get a zero of the appropriate type  \n",
    "    for i in 1:n_particles-1, j in i+1:n_particles\n",
    "        accu += @views Vpair(norm(x[:, i] .- x[:, j]))\n",
    "    end\n",
    "    accu\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can test it on 5000 particles:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using BenchmarkTools\n",
    "x = randn(2, 5000)\n",
    "@time Vtot(Vmorse, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Puh ... that's getting slow. An obvious place to parallelise is the accumulation over particles. For this we will use the `Threads.@threads` macro. The other required change is to use an atomic accumulator to avoid inconsistencies:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import Base.Threads: @threads, Atomic\n",
    "\n",
    "function pVtot(Vpair, x)\n",
    "    n_particles = size(x, 2)\n",
    "    T = eltype(x)\n",
    "    accu = Atomic{T}(zero(T))    # Atomic accumulator\n",
    "    @threads for i in 1:n_particles-1\n",
    "        local_accu = zero(T)     # Thread-local accumulator\n",
    "        for j in i+1:n_particles\n",
    "            local_accu += @views Vpair(norm(x[:, i] .- x[:, j]))\n",
    "        end\n",
    "        accu[] += local_accu\n",
    "    end\n",
    "    accu[]\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "@time pVtot(Vmorse, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Of course this is still far from perfect, since the load between the threads is not exactly balanced, but we'll leave it here. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/manual/parallel-computing/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Distributed computing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For distributed computing over mulitple nodes, Julia offers `DistributedArrays`, a package for spreading array data across workers and `Distributed` from the standard library for classic building blocks of distributed computing (parallel maps, reductions, channels ...)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### More details\n",
    "- https://juliaparallel.github.io/DistributedArrays.jl/latest/\n",
    "- https://docs.julialang.org/en/v1/stdlib/Distributed/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Mixing Julia and Fortran, C, python, R, Java, ...\n",
    "\n",
    "Mixing Julia and C, C++ or **Fortran** works very well. The underlying reason is because Julia gets compiled via the LLVM compiler infrastructure, which is capable of the other languages out of the box.\n",
    "\n",
    "For example we can directly call the dot product function `ddot` from LAPACK:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Call to lapack to compute dot product\n",
    "v = [1.0, 2.0, 3.0, 4.0]\n",
    "w = [2.0, 3.0, 4.0, 5.0]\n",
    "VdotW = ccall((:ddot_, \"liblapack\"), Float64,\n",
    "              (Ref{Cint}, Ptr{Cdouble}, Ref{Cint}, Ptr{Cdouble}, Ref{Cint}),\n",
    "              length(v), v, 1, w, 1)\n",
    "println(\"v ⋅ w = $VdotW\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(actually ... this should not be done like this in practice, because even if one wants to avoid the builtin julia function `dot`, there are much better LAPACK wrappers in `LinearAlgebra.LAPACK`).\n",
    "\n",
    "As another example, take the following **piece of C code**:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "code = \"\"\"\n",
    "double sum_array(double* array, int n) {\n",
    "    double accu = 0.0;\n",
    "    for (int i = 0; i < n; ++i) {\n",
    "        accu += array[i];\n",
    "    }\n",
    "    return accu;\n",
    "}\n",
    "\"\"\";"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Which we first compile using `cc` into the shared library file `libsums.so`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "open(\"sums.c\", \"w\") do f\n",
    "    write(f, code)\n",
    "end\n",
    "run(`cc -shared -o libsums.so sums.c`)\n",
    "rm(\"sums.c\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using `sum_array` from Julia now similarly only takes an appropriate `ccall`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "v = [1.0, 2.0, 3.0]\n",
    "res = ccall((:sum_array, \"libsums\"), Cdouble,\n",
    "            (Ptr{Cdouble}, Cint), v, length(v))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Of course compiling each time is not so convenient. If you require mixing C++ and Julia a lot, have a look at the `Cxx` package. As a bonus it even adds an interactive C++ shell to Julia's REPL ..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By far the most advanced foreign-language integration has been achieved with **Python**. Here you can easily use *any* python-package from Julia.\n",
    "For example NumPy:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "using PyCall\n",
    "np = pyimport(\"numpy\")\n",
    "\n",
    "np.random.rand(3, 4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A lot of type conversion happens in the background automatically. In this case NumPy arrays are automatically converted to Julia (and thus from row-major to column-major ...). The reverse (i.e. using Julia from python) is also possible, see [Using_Julia_From_Python](08_Using_Julia_From_Python.ipynb) for a quick demonstration."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Other supported languages include **[Java](https://github.com/JuliaInterop/JavaCall.jl) or [R](https://github.com/JuliaInterop/RCall.jl)**, which can be used to import functionality from the respective languages into Julia.\n",
    "\n",
    "##### More details\n",
    "- https://docs.julialang.org/en/v1/manual/calling-c-and-fortran-code/\n",
    "- https://github.com/JuliaPy/PyCall.jl\n",
    "- https://github.com/JuliaInterop/Cxx.jl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Many more features ...\n",
    "\n",
    "... you can find in the [Julia documentation](https://docs.julialang.org) ;)"
   ]
  }
 ],
 "metadata": {
  "@webio": {
   "lastCommId": null,
   "lastKernelId": null
  },
  "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
}