{ "cells": [ { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "view-in-github" }, "source": [ "\"Open" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "K5K3W5LX0s6i" }, "source": [ "# Post processing of mumax3 output data in python" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "pXovm2LT0s6k" }, "source": [ "When mumax3 is used for actual research, then most likely the RAW output data will need some further processing. \n", "This jupyter notebook demonstrates how mumax3 output data can be processed, analyzed, and visualized using Python. After defining mumax3 helper functions to run mumax3 and read the output data in Python, we will study three different use cases.\n", "\n", "Table of contents:\n", "1. [getting ready](#gettingready)\n", "2. [mumax3 helper functions](#helperfunctions)\n", "3. [standard problem 4 revisited](#std4)\n", "4. [skyrmion excitation (temporal fft)](#skyrmionexcitation)\n", "5. [spin wave dispersion (temporal + spatial fft)](#spinwaves)" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "hnoBIsYw0s6l" }, "source": [ " \n", "## 1. Getting ready" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "Gv71NUjz0s6m" }, "source": [ "If you want to run the examples in this notebook, then mumax3 must be is installed on the system that is running this notebook. Go to [mumax.github.io/download](http://mumax.github.io/download) for more information on how to install mumax3 on your system. If you are running this notebook in a google colaboratory session, then it suffices to run the cell below to install mumax3. " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 366 }, "colab_type": "code", "id": "9dDyrCd90s6o", "outputId": "610ede0e-3c92-4f60-b400-0935f0ea56c9" }, "outputs": [], "source": [ "try:\n", " import google.colab\n", "except ImportError:\n", " pass\n", "else:\n", " !wget https://mumax.ugent.be/mumax3-binaries/mumax3.10_linux_cuda10.1.tar.gz\n", " !tar -xvf mumax3.10_linux_cuda10.1.tar.gz\n", " !rm mumax3.10_linux_cuda10.1.tar.gz\n", " !rm -rf mumax3.10 && mv mumax3.10_linux_cuda10.1 mumax3.10\n", " import os\n", " os.environ['PATH'] += \":/content/mumax3.10\"" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "1Si5CtCG0s6w" }, "source": [ "In the examples presented in this notebook, we will use the numpy and pandas libraries for post processing mumax3 output data, and matplotlib to visualize this data. So let's import these libraries." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "rmod7y4x0s6y" }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import pandas as pd" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "bZIHecFr0s63" }, "source": [ " \n", "## 2. mumax3 helper functions" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "gyY3fg800s64" }, "source": [ "Let's start by writing a function which reads a mumax3 output table and puts the date in a Pandas dataframe." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "W7Nn_7Gl0s66" }, "outputs": [], "source": [ "def read_mumax3_table(filename):\n", " \"\"\"Puts the mumax3 output table in a pandas dataframe\"\"\"\n", "\n", " from pandas import read_table\n", " \n", " table = read_table(filename)\n", " table.columns = ' '.join(table.columns).split()[1::2]\n", " \n", " return table" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "3HkwQV5f0s6-" }, "source": [ "Mumax3 does not only write output data in a table, it can also write field data to ovf files. The function below converts all ovf files in the output directory to numpy files. These files are then loaded using the numpy.load function. The result is a python dictionary of the field data (a numpy array) with the ovf filename as key." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "apKunwIp0s6_" }, "outputs": [], "source": [ "def read_mumax3_ovffiles(outputdir):\n", " \"\"\"Load all ovffiles in outputdir into a dictionary of numpy arrays \n", " with the ovffilename (without extension) as key\"\"\"\n", " \n", " from subprocess import run, PIPE, STDOUT\n", " from glob import glob\n", " from os import path\n", " from numpy import load\n", "\n", " # convert all ovf files in the output directory to numpy files\n", " p = run([\"mumax3-convert\",\"-numpy\",outputdir+\"/*.ovf\"], stdout=PIPE, stderr=STDOUT)\n", " if p.returncode != 0:\n", " print(p.stdout.decode('UTF-8'))\n", "\n", " # read the numpy files (the converted ovf files)\n", " fields = {}\n", " for npyfile in glob(outputdir+\"/*.npy\"):\n", " key = path.splitext(path.basename(npyfile))[0]\n", " fields[key] = load(npyfile)\n", " \n", " return fields" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "HFgCSbwf0s7E" }, "source": [ "The function below executes a mumax3 script and returns the data of the output table in a pandas dataframe, and the saved fields as numpy arrays." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "wKe6laCH0s7F" }, "outputs": [], "source": [ "def run_mumax3(script, name, verbose=False):\n", " \"\"\" Executes a mumax3 script and convert ovf files to numpy files\n", " \n", " Parameters\n", " ----------\n", " script: string containing the mumax3 input script\n", " name: name of the simulation (this will be the name of the script and output dir)\n", " verbose: print stdout of mumax3 when it is finished\n", " \"\"\"\n", " \n", " from subprocess import run, PIPE, STDOUT\n", " from os import path\n", "\n", " scriptfile = name + \".txt\" \n", " outputdir = name + \".out\"\n", "\n", " # write the input script in scriptfile\n", " with open(scriptfile, 'w' ) as f:\n", " f.write(script)\n", " \n", " # call mumax3 to execute this script\n", " p = run([\"mumax3\",\"-f\",scriptfile], stdout=PIPE, stderr=STDOUT)\n", " if verbose or p.returncode != 0:\n", " print(p.stdout.decode('UTF-8'))\n", " \n", " if path.exists(outputdir + \"/table.txt\"):\n", " table = read_mumax3_table(outputdir + \"/table.txt\")\n", " else:\n", " table = None\n", " \n", " fields = read_mumax3_ovffiles(outputdir)\n", " \n", " return table, fields" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "GlPfPFFV0s7K" }, "source": [ " \n", "## 3. Standard problem 4 revisited" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "vy2OXKAh0s7L" }, "source": [ "Let's start by putting the mumax3 input script for standard problem 4 in the string variable `script`." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "nX50Gun30s7L" }, "outputs": [], "source": [ "script=\"\"\"\n", "SetGridsize(128, 32, 1)\n", "SetCellsize(500e-9/128, 125e-9/32, 3e-9)\n", "\n", "Msat = 800e3\n", "Aex = 13e-12\n", "alpha = 0.02\n", "\n", "m = uniform(1, .1, 0)\n", "relax()\n", "\n", "autosave(m, 200e-12)\n", "tableadd(e_total)\n", "tableautosave(10e-12)\n", "\n", "B_ext = vector(-24.6E-3, 4.3E-3, 0)\n", "run(1e-9)\n", "\"\"\"" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "hCWtYmCe0s7P" }, "source": [ "Now we can execute this mumax3 script using the `run_mumax3` helper function. This function returns the output table in a pandas dataframe and a dictionary of the saved fields. If `verbose=True`, then the log output of mumax3 will be printed out when mumax3 finishes. Note that on the left side of the cell, you can check if the simulation is finished or still running." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 54 }, "colab_type": "code", "id": "Hcv86Mr90s7Q", "outputId": "ae7dc897-916c-43ca-ce60-e94710700204" }, "outputs": [], "source": [ "table, fields = run_mumax3( script, name=\"standardproblem4\", verbose=False )" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "FKEfolRG0s7U" }, "source": [ "The table data is put in a pandas dataframe. This makes it very easy to plot and analyze the table data." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 282 }, "colab_type": "code", "id": "okLIFMBU0s7U", "outputId": "25408939-82da-4503-e26e-a235c5cbbf71" }, "outputs": [], "source": [ "print(table)\n", "\n", "plt.figure()\n", "\n", "nanosecond = 1e-9\n", "plt.plot( table[\"t\"]/nanosecond, table[\"mx\"])\n", "plt.plot( table[\"t\"]/nanosecond, table[\"my\"])\n", "plt.plot( table[\"t\"]/nanosecond, table[\"mz\"])\n", "\n", "plt.xlabel(\"Time (ns)\")\n", "plt.ylabel(\"Magnetization\")\n", "\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "ObqhouM80s7Y" }, "source": [ "Let's print out the keys of the field dictionary to get an idea about the content of this dictionary:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print(fields.keys())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As you can see, the dictionary contains 6 fields with keys corresponding to the ovf filenames (without extension). Now, let's focus on the first field 'm000000'." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "m = fields[\"m000000\"]\n", "\n", "print(\"type =\", type(m))\n", "print(\"shape =\", m.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that the field data is stored in a 4D numpy array: `m[ic,iz,iy,ix]` with `ic` the index for the magnetization component (X=0,Y=1,Z=2), and with `(ix,iy,iz)` the cell index. We can now use the numpy and matplotlib libraries to process and visualize the field date. E.g. the cell below demonstrates how one can easily create an intensity plot of the absolute value of $m_y$ " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "vgFynTqy0s7Z" }, "outputs": [], "source": [ "def show_abs_my(m):\n", " my_abs = np.abs( m[1,0,:,:] )\n", " plt.figure()\n", " plt.imshow(my_abs, vmin=0, vmax=1, cmap=\"afmhot\")\n", " plt.show()\n", "\n", "show_abs_my(fields[\"m000001\"])" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "PSxMSAIT0s7c" }, "source": [ " \n", "## 4. Skyrmion excitation\n", "\n", "In this section, we will study the excitation spectrum of a skyrmion in a nanodisc. This example is largely based on [this PRB paper](https://journals.aps.org/prb/pdf/10.1103/PhysRevB.90.064410) by Joo-Von Kim and Felipe Garcia-Sanchez.\n", "\n", "The mumax3 script belows excites the skyrmion with a sinc pulse with maximum frequency of 50GHz. This means that we will excite modes with frequencies below 50GHz. " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "g3_uu7am0s7d" }, "outputs": [], "source": [ "# NUMERICAL PARAMETERS RELEVANT FOR THE SPECTRUM ANALYSIS\n", "fmax = 50e9 # maximum frequency (in Hz) of the sinc pulse\n", "T = 2e-9 # simulation time (longer -> better frequency resolution)\n", "dt = 1/(2*fmax) # the sample time (Nyquist theorem taken into account)\n", "\n", "# Note that this is a format string, this means that the statements inside the \n", "# curly brackets get evaluated by python. In this way, we insert the values of \n", "# the variables above in the script.\n", "script=f\"\"\"\n", "diam := 100e-9\n", "nx := 32\n", "setgridsize(nx,nx,1)\n", "setcellsize(diam/nx,diam/nx,1e-9)\n", "setgeom(circle(diam))\n", "\n", "Msat = 1e6\n", "Aex = 15e-12\n", "Dind = 3.0e-3\n", "Ku1 = 1e6\n", "AnisU = vector(0,0,1)\n", "alpha = 0.001\n", "\n", "t0 := 1/({fmax})\n", "B_ext = vector(0, 0, 0.0005 * sinc( 2*pi*{fmax}*(t-t0)))\n", "m = neelskyrmion(-1,1)\n", "minimize()\n", "\n", "autosave(m,{dt})\n", "tableautosave({dt})\n", "run({T})\n", "\"\"\"\n", " \n", "table, fields = run_mumax3(script,\"skyrmion\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The excitation spectrum of the skyrmion, can be found by taking the Fast Fourier Transform (FFT) of the spatial average of the magnetization deviation $m_z(t)-m_z(0)$." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "HFsrRoMD0s7i" }, "outputs": [], "source": [ "# FAST FOURIER TRANSFORM\n", "dm = table[\"mz\"] - table[\"mz\"][0] # average magnetization deviaton\n", "spectr = np.abs(np.fft.fft(dm)) # the absolute value of the FFT of dm\n", "freq = np.linspace(0, 1/dt, len(dm)) # the frequencies for this FFT\n", "\n", "# PLOT THE SPECTRUM\n", "plt.plot(freq/1e9, spectr)\n", "plt.xlim(0,fmax/1e9)\n", "plt.ylabel(\"Spectrum (a.u.)\")\n", "plt.xlabel(\"Frequency (GHz)\")\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are clearly two peaks in the excitation spectrum. One can find the corresponding frequencies using a peak finding algorithm. Here, we will just determine the position of the peaks visually." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# MODE IDENTIFICATION\n", "mode1_idx = 7 # determined visually\n", "mode2_idx = 58 # determined visually\n", "mode1_freq = freq[mode1_idx]\n", "mode2_freq = freq[mode2_idx]\n", "print(\"Mode 1 frequency: %.2f GHz\"%(mode1_freq/1e9))\n", "print(\"Mode 2 frequency: %.2f GHz\"%(mode2_freq/1e9))\n", "\n", "# PLOT THE SPECTRUM\n", "plt.plot(freq/1e9, spectr)\n", "plt.axvline(mode1_freq/1e9, lw=1, ls='--', c='gray')\n", "plt.axvline(mode2_freq/1e9, lw=1, ls='--', c='gray')\n", "plt.xlim(0,fmax/1e9)\n", "plt.ylabel(\"Spectrum (a.u.)\")\n", "plt.xlabel(\"Frequency (GHz)\")\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We have found the resonant frequencies by applying a FFT on the average magnetization. In addition, we can apply the FFT on the magnetization in every cell to get more information on these resonant modes." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "A81Do_8B0s7l" }, "outputs": [], "source": [ "# Stack all snapshots (4D arrays) of the magnetization on top of each other\n", "# The results in a single 5D array (first index is the snapshot index)\n", "m = np.stack([fields[key] for key in sorted(fields.keys())])\n", "\n", "# Select the z component and the (only) layer z=0\n", "mz = m[:,2,0,:,:]\n", "\n", "# Apply the FFT for every cell\n", "mz_fft = np.fft.fft(mz, axis=0)\n", "\n", "# Select the the two modes\n", "mode1 = mz_fft[mode1_idx]\n", "mode2 = mz_fft[mode2_idx]\n", "\n", "# Plot the result\n", "plt.figure(figsize=(10,4))\n", "plt.subplot(1,3,1)\n", "plt.title(\"$m_z$\")\n", "plt.imshow(mz[0])\n", "plt.subplot(1,3,2)\n", "plt.title(\"Mode 1\")\n", "plt.imshow(np.abs(mode1)**2)\n", "plt.subplot(1,3,3)\n", "plt.title(\"Mode 2\")\n", "plt.imshow(np.abs(mode2)**2)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "wXG9W0oG0s7p" }, "source": [ " \n", "## 5. Ferromagnetic spinwave dispersion relation" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "IenALN_20s7q" }, "source": [ "Disregarding dipole-dipole interactions, the dispersion relation $f(k)$ in a ferromagnetic wire is given by\n", "\n", "$$ f(k) = \\frac{\\gamma}{2\\pi}\\left[ \\frac{2A}{M_s} k^2 + B \\right] $$\n", "\n", "with $A$ the exchange stiffness, $M_s$ the saturation magnetization, and $B$ an externally applied field. In this section, we will try to reproduce this dispersion relation numerically using mumax3.\n", "\n", "The mumax3 script below simulates a uniformly magnetized nanowire with an applied field perpendicular to the wire. Spin waves are excited by a sinc pulse with a maximum frequency of 20GHz at the center of the simulation box. Note that the demagnetization is disabled in this simulation." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "cgSybvYj0s7u" }, "outputs": [], "source": [ "# NUMERICAL PARAMETERS\n", "fmax = 20e9 # maximum frequency (in Hz) of the sinc pulse\n", "T = 1e-8 # simulation time (longer -> better frequency resolution)\n", "dt = 1/(2*fmax) # the sample time\n", "dx = 4e-9 # cellsize\n", "nx = 1024 # number of cells\n", "\n", "# MATERIAL/SYSTEM PARAMETERS\n", "Bz = 0.2 # Bias field along the z direction\n", "A = 13e-12 # exchange constant\n", "Ms = 800e3 # saturation magnetization\n", "alpha = 0.05 # damping parameter\n", "gamma = 1.76e11 # gyromagnetic ratio\n", "\n", "# Note that this is a format string, this means that the statements inside the \n", "# curly brackets get evaluated by python. In this way, we insert the values of \n", "# the variables above in the script.\n", "script=f\"\"\"\n", "setgridsize({nx},1,1)\n", "setcellsize({dx},{dx},{dx})\n", "\n", "enabledemag = false\n", "Aex = {A}\n", "Msat = {Ms}\n", "alpha = {alpha}\n", "\n", "Bz := {Bz}\n", "B_ext = vector(0,0,{Bz})\n", "defregion(1,rect(2*{dx},inf))\n", "B_ext.setregion(1, vector(0.01 * sinc( 2*pi*{fmax}*(t-{T}/2)), 0, {Bz}))\n", "\n", "m = uniform(0,0,1) \n", "autosave(m,{dt})\n", "run({T})\n", "\"\"\"\n", " \n", "table, fields = run_mumax3(script,\"spinwaves\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Applying the 2D Fast Fourier Transform on the magnetization snapshots yield the dispersion relation of the ferromagnetic nanowire." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": {}, "colab_type": "code", "id": "Pi8qUMmO0s7w" }, "outputs": [], "source": [ "# Stack all snapshots of the magnetization on top of each other\n", "m = np.stack([fields[key] for key in sorted(fields.keys())])\n", "\n", "# Select the x component\n", "mx = m[:,0,0,0,:]\n", "\n", "# Apply the two dimensional FFT\n", "mx_fft = np.fft.fft2(mx)\n", "mx_fft = np.fft.fftshift(mx_fft)\n", "\n", "plt.figure(figsize=(10,6))\n", "\n", "# Show the intensity plot of the 2D FFT\n", "extent = [ -(2*np.pi)/(2*dx), (2*np.pi)/(2*dx), -1/(2*dt), 1/(2*dt)] # extent of k values and frequencies\n", "plt.imshow(np.abs(mx_fft)**2, extent=extent, aspect='auto', origin='lower', cmap=\"inferno\")\n", "\n", "# Plot the analytical derived dispersion relation \n", "k = np.linspace(-2e8,2e8,1000)\n", "freq_theory = A*gamma*k**2 /(np.pi*Ms) + gamma*Bz /(2*np.pi)\n", "plt.plot(k,freq_theory,'r--',lw=1)\n", "plt.axhline(gamma*Bz/(2*np.pi),c='g',ls='--',lw=1)\n", "\n", "plt.xlim([-2e8,2e8])\n", "plt.ylim([0,fmax])\n", "plt.ylabel(\"$f$ (Hz)\")\n", "plt.xlabel(\"$k$ (1/m)\")\n", "\n", "plt.show()" ] } ], "metadata": { "@webio": { "lastCommId": null, "lastKernelId": null }, "accelerator": "GPU", "colab": { "include_colab_link": true, "name": "postprocessing.ipynb", "provenance": [] }, "kernelspec": { "display_name": "Python 3", "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.8.3" } }, "nbformat": 4, "nbformat_minor": 1 }