01_Basics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 1: Creating magnetic structure and calculating 1D mPDF patterns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Contents of this tutorial"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Class structure of the mpdf package\n",
    "- Simple example: Calculating the mPDF from a spin dimer\n",
    "- Creating a magnetic structure from a CIF file using diffpy.Structure\n",
    "- Creating a magnetic structure from an MCIF file\n",
    "- Creating a magnetic structure by defining your own unit cell\n",
    "- Creating a magnetic structure with multiple species: Simple ferrimagnet\n",
    "- Creating more complex magnetic structures: One-dimensional spin helix\n",
    "- Playing with some calculation parameters\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Class structure of the mpdf package"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The core classes in the mpdf package the MPDFcalculator class, the MagStructure class, and the optional (but often very useful) MagSpecies class. There are also numerous other functions and classes to aid in mPDF calculations, many of which we will introduce here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### MagStructure class"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The main job of a MagStructure object is to contain a numpy array of atomic positions (MagStructure.atoms) and a corresponding numpy array of magnetic moment vectors (MagStructure.spins), which are the only two absolutely required ingredients for calculating the mPDF. Additional optional information such as the Lande g-factor or magnetic form factor can be stored in a MagStructure object if it is needed (as it often is for performing mPDF refinements). A diffpy.structure object corresponding to the atomic structure of the material can also be stored for convenience when generating the atomic positions and spin vectors. Multiple types of magnetic species (e.g. spins on different ionic sites that may have different moment sizes) can be stored inside a single MagStructure object."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### MagSpecies class"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The MagSpecies class is intended to provide an easy way to generate the atomic positions and spin vectors of a particular type of magnetic species in a structure. It takes a diffpy.structure object and packages additional information pertaining to the magnetic structure, such as the magnetic propagation vector(s), magnetic basis vector(s), spin/orbital quantum numbers, and magnetic form factor. Alternatively, the user can define a unit cell populated with positions and spin orientations of the magnetic atoms, and this will be used rather than a diffpy.structure object. The structural information is used in class methods that automatically generate the atomic positions and corresponding spin vectors, which can then be stored in a MagStructure object. There are no limits on the number of propagation and basis vectors, allowing for arbitrarily complex magnetic structures. MagSpecies objects can be loaded directly into a MagStructure object or created from inside the MagStructure object. Although a MagSpecies object is not required to calculate the mPDF, it will usually provide the easiest way to populate the atom and spin arrays that are required for the mPDF. The reason for having separate MagStructure and MagSpecies classes is that a single magnetic structure may contain multiple magnetic species."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### MPDFcalculator class"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The MPDFcalculator class handles the details of the calculation and contains information that is not directly related to the magnetic structure, such as the real-space calculation range, instrumental resolution effects, scale factors, and experimental parameters like q<sub>min</sub> and q<sub>max</sub>. It requires a MagStructure object as input, from which it extracts the atomic positions and spin vectors to be used in calculating the mPDF. The mPDFcalculator class has options for calculating both the properly normalized and the unnormalized mPDF (the unnormalized mPDF is frequently what we obtain experimentally)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Simple example: Calculating the mPDF from a pair of spins"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will now create a very simple magnetic structure consisting of just two spins and then calculate the corresponding mPDF."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Import the necessary libraries\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from diffpy.mpdf import *\n",
    "from diffpy.structure import loadStructure\n",
    "\n",
    "### Set all plots to be inline\n",
    "%matplotlib notebook\n",
    "\n",
    "# Note: I often have to run this cell twice to get the\n",
    "# matplotlib plots to display in interactive mode"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Create a MagStructure object\n",
    "mstruc = MagStructure()\n",
    "\n",
    "### Create two atoms in the structure at the positions (0,0,0) and (4,0,0)\n",
    "mstruc.atoms = np.array([[0,0,0],[4,0,0]])\n",
    "\n",
    "### Create two spin vectors corresponding to the atoms. Let's make them antiferromagnetic.\n",
    "mstruc.spins = np.array([[0,0,1],[0,0,-1]])\n",
    "\n",
    "### Create the MPDFcalculator object and load the magnetic structure into it\n",
    "mc = MPDFcalculator(mstruc)\n",
    "mc.rmin = 0\n",
    "mc.rmax = 15\n",
    "### Calculate and plot the mPDF!\n",
    "r, fr = mc.calc() # Use calc() if you want to extract the numerical results of the calculation\n",
    "mc.plot() # Use plot() if you just want to plot the mPDF without extracting the numerical arrays."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Comment on the sign and position of the peak in the mPDF."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's make a copy of this structure but make the dimer ferromagnetic instead of antiferromagnetic."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### make a copy of the magnetic structure\n",
    "mstruc2 = mstruc.copy()\n",
    "\n",
    "### Set the spins to be ferromagnetic\n",
    "mstruc2.spins = ????\n",
    "\n",
    "### Create another mPDF calculator\n",
    "mc2 = MPDFcalculator(mstruc2)\n",
    "\n",
    "### Calculate the mPDF from the new calculator\n",
    "r2, fr2 = ????\n",
    "\n",
    "### Compare the antiferromagnetic and ferromagnetic mPDFs\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r,fr,'b-',label='Antiferro')\n",
    "ax.plot(r2,fr2,'r-',label='Ferro')\n",
    "\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'f ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend()\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise\n",
    "Calculate the mPDF for a 1D chain of 10 atoms with alternating magnetization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Creating a magnetic structure from a CIF file using diffpy.structure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This example will show how to use the MagSpecies class and additional features in the MagStructure class to quickly generate the atomic positions and spin vectors from a diffpy.Structure object. It will also show how to calculate the unnormalized mPDF, which is useful for comparison to and refinement against neutron total scattering data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load the atomic structure from the CIF file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now let's create a diffpy.structure object from a CIF file for MnO. This has all the atomic information,\n",
    "### but none of the magnetic information. Note that a .stru file generated from PDFgui or diffpy.structure\n",
    "### would also work just as well.\n",
    "mno = loadStructure('files/MnO_cubic.cif')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### As a refresher, let's take a look at the mno structure object.\n",
    "print(mno)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create the MagSpecies object"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### We will now create a MagSpecies object to extend the mno structure object.\n",
    "mspec = MagSpecies()\n",
    "\n",
    "### Load the mno structure and give the magnetic species a label\n",
    "mspec.struc = mno\n",
    "mspec.label = 'Mn2+' # you can choose anything for the label\n",
    "\n",
    "### Now we need to tell it which atoms in MnO are magnetic. Consult the printed\n",
    "### list of atoms from the previous cell and provide the (zero-based) indices of the Mn atoms.\n",
    "mspec.strucIdxs = [????]\n",
    "\n",
    "### Now we provide the magnetic propagation and basis vectors, which are known from previous neutron\n",
    "### diffraction studies of MnO.\n",
    "k = np.array([0.5,0.5,0.5])\n",
    "s = np.array([1,-1,0]) \n",
    "mspec.kvecs = np.array([k])\n",
    "mspec.basisvecs = np.array([s])\n",
    "\n",
    "### Below, we specify the origin relative to which the propagation vectors are defined.\n",
    "### It should always correspond to the Cartesian coordinates of the position of one of\n",
    "### the spins in the structure. By default it is zero, so the line below isn't necessary\n",
    "### in this particular case, but we show it just for demonstrative purposes.\n",
    "mspec.origin = np.array([0,0,0])\n",
    "\n",
    "### Now we provide information about the magnetic form factor. We tell the MagSpecies object the type of magnetic\n",
    "### ion, and it looks up magnetic form factor in a table.\n",
    "mspec.ffparamkey = 'Mn2'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create the MagStructure object"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now we can create a MagStructure object and load mspec into it.\n",
    "mstruc = MagStructure()\n",
    "mstruc.loadSpecies(mspec)\n",
    "### Now we will generate the atomic positions and spins. It is important to do the atoms first, since the\n",
    "### spins are generated by applying the propagation and basis vectors to the atomic positions. These methods\n",
    "### use the information contained in the atomic and magnetic structures to generate arrays of atomic positions\n",
    "### and spin vectors.\n",
    "mstruc.makeAtoms()\n",
    "mstruc.makeSpins()\n",
    "\n",
    "### And we make the magnetic form factor:\n",
    "mstruc.makeFF()\n",
    "\n",
    "### Alternative approach: run mstruc.makeAll() and it takes care of the three previous commands."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### View the magnetic structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### diffpy.mpdf includes a very rudimentary magnetic structure viewer\n",
    "### which can be useful at times. Let's visualize all the spins found\n",
    "### within 6 Angstroms of the origin.\n",
    "distances = np.apply_along_axis(np.linalg.norm, 1, mstruc.atoms)\n",
    "viewmask = distances < 6\n",
    "mstruc.visualize(mstruc.atoms[viewmask], mstruc.spins[viewmask], showcrystalaxes=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now we can create the mPDF calculator and load the magnetic structure\n",
    "mc = MPDFcalculator()\n",
    "mc.magstruc = mstruc\n",
    "\n",
    "### Plot the mPDF\n",
    "????"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Creating a magnetic structure from an MCIF file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you have an MCIF file available (see, for example, the database at http://webbdcrista1.ehu.es/magndata/ ), then you can create a MagStructure object with just one line of code."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mcif = 'files/1.31_MnO.mcif'\n",
    "mstruc = create_from_mcif(mcif, ffparamkey='Mn2')\n",
    "mstruc.makeAll()\n",
    "\n",
    "mc = MPDFcalculator(mstruc)\n",
    "\n",
    "mc.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise\n",
    "Go to the MAGNDATA website and download the MCIF file for KMnF4. Plot the mPDF for this compound."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "????"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Creating a magnetic structure by defining your own unit cell"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This example will show you how to create a magnetic structure from a unit cell that you define yourself. You will have to provide the lattice vectors, the positions of the magnetic atoms in the unit cell, and the magnetic moments corresponding to those atoms."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Create the magnetic species and turn off the diffpy.structure option.\n",
    "mspec = MagSpecies(useDiffpyStruc=False)\n",
    "\n",
    "### Define the lattice vectors of the unit cell. Let's make a cubic unit cell.\n",
    "a = 4.0\n",
    "mspec.latVecs = np.array([[a,0,0],????])\n",
    "\n",
    "### Define the positions of the magnetic atoms in the unit cell (in fractional coordinates). Let's make\n",
    "### We'll make a body-centered structure.\n",
    "mspec.atomBasis = np.array([[0,0,0], [0.5,0.5,0.5]])\n",
    "\n",
    "### Define the magnetic moments in the same order as the list of atoms.\n",
    "### Let's make an antiferromagnet.\n",
    "mspec.spinBasis = np.array([[1,1,1], [????]])\n",
    "\n",
    "\n",
    "### Create the magnetic structure object and load mspec.\n",
    "mstruc = MagStructure()\n",
    "mstruc.loadSpecies(mspec)\n",
    "mstruc.makeAtoms()\n",
    "mstruc.makeSpins()\n",
    "\n",
    "### Let's visualize the first unit cell to make sure we have what we expect.\n",
    "visAtoms = np.array([[0,0,0],[a,0,0],[0,a,0],[0,0,a],[a,a,0],\n",
    "                     [a,0,a],[0,a,a],[a,a,a],[0.5*a,0.5*a,0.5*a]])\n",
    "visSpins = mstruc.spinsFromAtoms(visAtoms, fractional=False)\n",
    "mstruc.visualize(visAtoms, visSpins)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now we can set up the MPDFcalculator and plot the mPDF.\n",
    "mc = MPDFcalculator(mstruc)\n",
    "mc.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Creating a magnetic structure with multiple species: Simple ferrimagnet"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will create a ferrimagnetic structure to illustrate the use of multiple magnetic species within a single magnetic structure. Let's build another antiferromagnetic body-centered cubic structure but with two different spin species, one with a large moment and one with a small moment."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Make a magnetic species with a large moment"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Create the first magnetic species and turn off the diffpy.structure option.\n",
    "mspec1 = MagSpecies(????)\n",
    "\n",
    "### Define the lattice vectors of the unit cell. Let's make a cubic unit cell.\n",
    "a = 4.0\n",
    "mspec1.latVecs = np.array([[a,0,0],[0,a,0],[0,0,a]])\n",
    "\n",
    "### Define the atomic position and magnetic moment.\n",
    "mspec1.atomBasis = np.array([0,0,0])\n",
    "mspec1.spinBasis = np.array([0,0,1])\n",
    "### if you have multiple species, they must have distinct labels. A random hex\n",
    "### string is generated by default, but you can also set your own label.\n",
    "mspec1.label = 'big'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Make a magnetic species with a small moment"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now make the other species, starting with mspec1 as a template\n",
    "mspec2 = mspec1.copy()\n",
    "mspec2.atomBasis = np.array([0.5,0.5,0.5])\n",
    "mspec2.spinBasis = np.array([????])\n",
    "mspec2.label = 'small'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create and view the magnetic structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Create the magnetic structure object and load the species.\n",
    "mstruc = MagStructure()\n",
    "mstruc.loadSpecies(????)\n",
    "mstruc.loadSpecies(????)\n",
    "mstruc.makeAll()\n",
    "### Note: when multiple MagSpecies are loaded into a MagStructure using the makeAll()\n",
    "### command, the calcIdxs attribute of the MagStructure is automatically updated to\n",
    "### include one spin from each distinct species.\n",
    "\n",
    "### Again, let's visualize the first unit cell to make sure we have what we expect.\n",
    "visAtoms = np.array([[0,0,0],[a,0,0],[0,a,0],[0,0,a],[a,a,0],[a,0,a],[0,a,a],[a,a,a],[0.5*a,0.5*a,0.5*a]])\n",
    "visSpins = mstruc.spinsFromAtoms(visAtoms, fractional=False)\n",
    "mstruc.visualize(visAtoms, visSpins, showcrystalaxes=True, axesorigin=np.array([-1,-1,-1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculate the mPDF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now we can set up the MPDFcalculator and plot the mPDF.\n",
    "mc = MPDFcalculator(????)\n",
    "mc.plot(linearTermMethod='autoslope')\n",
    "# Note: we use the linearTermMethod argument for magnetic structures\n",
    "# with a net magnetization, such as this ferrimagnet. For more info,\n",
    "# see the ferromagnet.ipynb notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Creating more complex magnetic structures: One-dimensional spin helix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This example shows how to use multiple magnetic propagation vectors to create a non-collinear magnetic structure--in this case, a one-dimensional spin helix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### We will load in a structure from a CIF file and modify it to simulate a 1-D material.\n",
    "astruc = loadStructure(\"files/MnO_cubic.cif\")\n",
    "astruc.lattice.a = 3.0\n",
    "astruc.lattice.b = 150.0 # we make the b and c lattice constants unphysically large\n",
    "astruc.lattice.c = 150.0\n",
    "\n",
    "### Create the magnetic species object.\n",
    "helix = MagSpecies(astruc)\n",
    "\n",
    "### Set up the magnetic propagation and basis vectors for a helical spin configuration.\n",
    "k = np.array([np.sqrt(2)/10,0,0]) # make the period incommensurate with the lattice\n",
    "helix.kvecs=np.array([k, -k])\n",
    "\n",
    "Sk = 0.5*(np.array([0,0,1])+0.5j*np.array([0,1,0])) # j is the numpy symbol for the imaginary unit\n",
    "helix.basisvecs = np.array([Sk,Sk.conj()])\n",
    "\n",
    "### Populate with atoms and spins.\n",
    "helix.rmaxAtoms = 70.0 # we will generate spins within a sphere of radius 70 A\n",
    "helix.makeAtoms()\n",
    "helix.makeSpins()\n",
    "helix.label = 'helix'\n",
    "\n",
    "### Create the magnetic structure object.\n",
    "mstruc = MagStructure()\n",
    "mstruc.loadSpecies(helix)\n",
    "\n",
    "#### Visualize the spins.\n",
    "x,y,z = mstruc.atoms.transpose()\n",
    "mask = np.logical_and(z==0,np.logical_and(y==0,np.abs(x)<30))\n",
    "visatoms = mstruc.atoms[mask]\n",
    "visspins = spinsFromAtoms(mstruc, visatoms, fractional=False)\n",
    "mstruc.visualize(visatoms, visspins)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Now plot the mPDF.\n",
    "mc = MPDFcalculator(mstruc)\n",
    "mc.rmax=70.0\n",
    "\n",
    "mc.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exploring some of the additional parameters in the mPDF"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we will take a look at some of the various calculation options for the mPDF, including:\n",
    "- Whether the normalized mPDF, unnormalized mPDF, or both are calculated (see Frandsen & Billinge, Acta Cryst. A71 (2015), 325-334; \n",
    "https://doi.org/10.1107/S205327331500306X\n",
    ").\n",
    "- Qmax to simulate experimental conditions.\n",
    "- Parameters that dampen and/or broaden the mPDF, simulating the effects of thermal motion and instrumental resolution.\n",
    "- \"Ordered\" scale factor corresponding to the magnitude of the locally ordered moment, and \"paramagnetic\" scale factor which depends only on the magnitude and spatial extent of a single localized moment, not any correlations between moments.\n",
    "- The r-range for the calculation.\n",
    "\n",
    "For this, we will create the antiferromagnetic MnO structure again."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create the MnO magnetic structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Create the MagStructure object and populate the atom/spin arrays.\n",
    "### Feel free to copy and paste from before."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculate the normalized and unnormalized mPDFs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mc = MPDFcalculator(????)\n",
    "r, fr = mc.calc() ### by default, just the normalized mPDF is calculated\n",
    "r, dr = mc.calc(normalized=False) ### setting the normalized arg to False calculates the unnormalized mPDF\n",
    "\n",
    "### Let's compare the two quantities\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, fr/fr.max(), 'b-', label='Normalized')\n",
    "ax.plot(r, dr/dr.max(), 'r-', label='Unnormalized')\n",
    "\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'f ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "\n",
    "### Or you can calculate them both together:\n",
    "r, fr, dr = mc.calc(both=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Model the effect of finite Qmax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Default mPDF calculator\n",
    "mc = MPDFcalculator(mstruc)\n",
    "r, fr = mc.calc()\n",
    "\n",
    "### Adjust qmax to simulate termination ripples\n",
    "mcA = mc.copy()\n",
    "mcA.qmax = 20.0\n",
    "rA,frA = mcA.calc()\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, fr, 'b-', label='Default')\n",
    "ax.plot(rA, frA, 'r-', label='With qmax')\n",
    "ax.plot(r, frA-fr-35, 'g-', label='Difference')\n",
    "\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'f ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Change the intrinisic mPDF peak width"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Change the intrinsic mPDF peak width to simulate thermal motion\n",
    "mcB = mc.copy()\n",
    "mstrucB = mstruc.copy()\n",
    "mstrucB.Uiso = 0.015 # atomic displacement parameter; default is 0.01. Try different values to see the effect.\n",
    "mcB.magstruc = mstrucB\n",
    "rB, frB = mcB.calc()\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, fr, 'b-', label='Default')\n",
    "ax.plot(rB, frB, 'r-', label='Broadened')\n",
    "ax.plot(r,frB - fr - 35, 'g-', label='Difference')\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'f ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Implement a finite correlation length"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Include a damping factor to simulate a finite correlation length\n",
    "### and/or instrumental resolution effects.\n",
    "mcC = mc.copy()\n",
    "mstrucC = mstruc.copy()\n",
    "mstrucC.corrLength = 10 # correlation length in Angstroms\n",
    "mcC.magstruc = mstrucC\n",
    "rC, frC = mcC.calc()\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, fr, 'b-', label='Default')\n",
    "ax.plot(rC, frC, 'r-', label='Finite correlation length')\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'f ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Change the ordered scale factor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Reduce the ordered scale factor, simulating a reduced ordered moment.\n",
    "mcD = mc.copy()\n",
    "mcD.ordScale = 0.5*mc.ordScale\n",
    "rD, frD, drD = mcD.calc(both=True) # both the normalized and unnormalized mPDF\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, dr, 'b-', label='Default')\n",
    "ax.plot(rD, drD, 'r-', label='Reduced scale factor')\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'd ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Change the calculation range"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Since the r-range is longer, we also need to re-generate the atoms and spins\n",
    "### to fill a sphere of at least radius 40 Angstroms. The code below does this.\n",
    "### See what happens if you set rmaxAtoms to something too small, like 10 or 15.\n",
    "mcif = 'files/1.31_MnO.mcif'\n",
    "mstruc = create_from_mcif(mcif, ffparamkey='Mn2', rmaxAtoms=40)\n",
    "mstruc.makeAll()\n",
    "\n",
    "### Change the r-range for the calculation from 20 Angstroms (default) to 40.\n",
    "mcE = MPDFcalculator(mstruc)\n",
    "mcE.rmax = 40 # in Angstroms\n",
    "\n",
    "\n",
    "rE,frE,drE = mcE.calc(both=True)\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.add_subplot(111)\n",
    "\n",
    "ax.plot(r, dr, 'bo', label='Default')\n",
    "ax.plot(rE, drE, 'r-', label='Longer r-grid',lw=2)\n",
    "ax.set_xlabel(r'r ($\\AA$)')\n",
    "ax.set_ylabel(r'd ($\\AA^{-2}$)')\n",
    "\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close('all')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercise\n",
    "Make a MagStructure object for MnTe (you can find it on MAGNDATA) and set the correlation length to 15 A. Calculate the mPDF out to a maximum r-value such that the mPDF signal is more than 90% suppressed.\n",
    "\n",
    "Hint: If you use the create_from_mcif function, you can pass an optional argument for rmaxAtoms. Otherwise, you will need to make sure each MagSpecies instance has rmaxAtoms set to an appropriate value."
   ]
  },
  {
   "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.7.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}