{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# IPython Notebook for turning in solutions to the problems in the Essentials of Paleomagnetism Textbook by L. Tauxe"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problems in Chapter 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem 1:  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Given that:\n",
    "\n",
    "$$\n",
    "\\nabla V_m = - \\bigl(\n",
    "{ {\\partial}\\over {\\partial r} }\n",
    "{ {m \\cos \\theta} \\over {4 \\pi r^2}} + \n",
    "{ {1\\over r} }\n",
    "{ {\\partial}\\over {\\partial \\theta} }\n",
    " { { m\\cos \\theta}\\over { 4 \\pi r^2} }\n",
    " \\bigr)\n",
    "$$\n",
    "\n",
    "it follows that:\n",
    "\n",
    "$$\n",
    "H_r = -{ {\\partial}\\over {\\partial r} }\n",
    "{ {m \\cos \\theta} \\over {4 \\pi r^2}} = { {2m\\cos \\theta } \\over {4 \\pi r^3}} = { {m\\cos \\theta } \\over {2 \\pi r^3} }\n",
    "$$\n",
    "and \n",
    "$$\n",
    "H_{\\theta}=\n",
    "{ {1\\over r} }\n",
    "{ {\\partial}\\over {\\partial \\theta} }\n",
    " { { m\\cos \\theta}\\over { 4 \\pi r^2} } \n",
    "={{ m\\sin \\theta} \\over {4 \\pi r^3}}.\n",
    "$$\n",
    "Plug in the values for $m, \\theta, r$:\n",
    "\n",
    "$$\n",
    "H_r = { {8 \\times 10^{22} \\cos(45)} \\over {2 \\pi \\cdot (6 \\times 10^6)^3} } [ {{Am^2}\\over{m^3} } ] = 41.7 A/m\n",
    "$$\n",
    "$$\n",
    "H_{\\theta} = \n",
    "{ {8 \\times 10^{22} \\sin(45)} \\over {4 \\pi \\cdot (6 \\times 10^6)^3} } [ {{Am^2}\\over{m^3} } ] = 20.8 A/m\n",
    "$$\n",
    "Because $B = \\mu_o H$, and $\\mu_o = 4 \\pi \\times 10^{-7} [ {{kg m^2} \\over {A^2s^2}}\\cdot {1\\over m} ]$,\n",
    "$B_r = \\mu_o H_r \n",
    "= 52 \\mu$T and $B_{\\theta} = 26 \\mu$T.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H_r=     41.7 H_theta=     20.8\n"
     ]
    }
   ],
   "source": [
    "# code to calculate H_r and H_theta\n",
    "import numpy as np\n",
    "m=80e21 # magnetic moment of the dipole in Am^2\n",
    "theta=np.radians(45) # 45 degrees converted to radians\n",
    "r= 6e6 # approximately the radius of the Earth in meters\n",
    "H_r=(m*np.cos(theta))/(2.*np.pi * r**3)\n",
    "H_theta = (m*np.sin(theta))/(4.*np.pi *r**3)\n",
    "print 'H_r= ','%7.1f'%(H_r), 'H_theta= ', '%7.1f'%(H_theta)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "B_r =  52 uT\n",
      "B_theta = 26 uT\n"
     ]
    }
   ],
   "source": [
    "mu_o=(4.*np.pi)*1e-7 # mu_o in Henry's per meter\n",
    "B_r=mu_o*H_r\n",
    "B_theta=mu_o*H_theta\n",
    "print 'B_r = ','%i'%(B_r*1e6), 'uT' # B_r in microtesla\n",
    "print 'B_theta =','%i'%(B_theta*1e6),'uT' # B_theta in microtesla"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem 2:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here is an interactive script:  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "B_SI:  35.0\n",
      "m_SI:  2.78e-23\n",
      "H_SI:  10185.9163579\n",
      "B_cgs:  3.500e+05\n",
      "m_cgs:  2.7800e-20\n",
      "H_cgs:  128.0\n"
     ]
    }
   ],
   "source": [
    "def convert_units(B,m,H,units):\n",
    "    \"\"\"\n",
    "    This function converts input values of B, m, and H in\n",
    "    units of either cgs or SI to SI and cgs respectively.\n",
    "    \"\"\"\n",
    "    if units==\"cgs\": # cgs=> SI\n",
    "        B_out = B*1e-4 # convert gauss to tesla\n",
    "        m_out= m*1e-3 # convert emu to Am^2\n",
    "        H_out= H/(4*np.pi*1e-3) # convert oe to A/m\n",
    "    elif units==\"SI\": # SI=> cgs\n",
    "        B_out= B*1e4 # convert  tesla to gauss\n",
    "        m_out=m*1e3 # convert Am^2 to emu\n",
    "        H_out=H*(4*np.pi*1e-3) # convert A/m to oe\n",
    "    else: # mistake\n",
    "        print 'mistake in your units, try again '   \n",
    "    return B_out,m_out,H_out # returns the output variables\n",
    "B,m,H,units=3.5e5,2.78e-20,128.,'cgs' # define some values\n",
    "B_SI,m_SI,H_SI = convert_units(B,m,H,units) # call the function\n",
    "print 'B_SI: ',B_SI # print out the converted values\n",
    "print 'm_SI: ',m_SI\n",
    "print 'H_SI: ',H_SI\n",
    "B_cgs,m_cgs,H_cgs = convert_units(B_SI,m_SI,H_SI,'SI') # do the revers\n",
    "print 'B_cgs: ','%5.3e'%(B_cgs)\n",
    "print 'm_cgs: ','%8.4e'%(m_cgs)\n",
    "print 'H_cgs: ',H_cgs\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Problem 3:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "a) This problem boils down to finding the value for ${\\bf m}$ in Equation 1.8 in Chapter 1 that would give rise to a radial field of 10$\\mu$T at a depth  of 2890 km (radius of the Earth minus radius of the dipole source).      Because $\\theta$ is 0$^{\\circ}$, we do not have to worry about the tangential component yet, only the radial field, which was given by:\n",
    "\n",
    "$$\n",
    "H_r = {1\\over {2\\pi}}{{m \\cos\\theta}\\over {r^3}}.\n",
    "$$\n",
    "\n",
    "We were given the field strength in units of tesla, so we have to convert to $H$ by the relation $H=B/\\mu_o $.   Taking  $\\cos \\theta$ as unity here, we solve for  $m$ and get: \n",
    "\n",
    "$$\n",
    "{m}={{2\\pi H_r r^3}} = {{2\\pi B_r r^3}\\over{\\mu_o}}\n",
    "$$\n",
    "\n",
    "Plugging  in the values, we get:\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " m =      1.2  ZAm^2\n"
     ]
    }
   ],
   "source": [
    "B_r=10e-6 # radial field in microtesla\n",
    "r=2890*1e3 # radius converted to meters\n",
    "m=2.*np.pi*B_r*r**3/mu_o\n",
    "print 'm = ', '%7.1f'%(m*1e-21) ,  ' ZAm^2' # m in ZAm^2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "b) To compare 10 $\\mu$T with the field produced by an axial dipole of 80 ZAm$^2$, we need the second part of Equation 1.8 in the text:\n",
    "\n",
    "$$\n",
    "H_{\\theta} =  {{m \\sin\n",
    "\\theta}\\over{4\\pi r^3}}.\n",
    "$$\n",
    "\n",
    "Realizing the $\\theta$ is 90 minus the latitude, we have $\\theta=30$ and \n",
    "$$\n",
    "{H_r} ={ {80 \\times 10^{21}\\times \\cos 30  }\\over {2 \\pi (6370 \\times 10^3)^3 }}\n",
    "$$\n",
    "and \n",
    "$$\n",
    "{H_{\\theta}} = { {80 \\times 10^{21}\\times  \\sin 30  }\\over {4 \\pi (6370 \\times 10^3)^3 }}\n",
    "$$\n",
    "which are  42.7   and  12.3 Am$^{-1}$ respectively.  The total field would be $\\sqrt{H_r^2+H_{\\theta}^2}$ or 44.4 Am$^{-1}$ .  To compare this with the value of 10$\\mu$T, we must convert this to units of $B$, so\n",
    "multiplying 44.4 by $\\mu_o$ we have 55.8 $\\mu$T.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H_r =     42.6  A/m\n",
      "H_theta =     12.3  A/m\n",
      "total field =    44.4  A/m\n",
      "in uT, this is =    55.8\n"
     ]
    }
   ],
   "source": [
    "mu_o=(4.*np.pi)*1e-7 # mu_o in Henry's per meter\n",
    "m_dipole=80e21 # axial dipole of the Earth's magnetic field in ZAm^2\n",
    "theta=90.-60. # 60 degree latitude expressed as a colatitude\n",
    "theta=np.radians(theta) # convert to radians\n",
    "r_earth=6371*1e3 # radius of the Earth in meters\n",
    "H_r= (m_dipole * np.cos(theta))/(2.*np.pi*r_earth**3)\n",
    "H_theta= (m_dipole * np.sin(theta))/(4.*np.pi*r_earth**3)\n",
    "print 'H_r = ','%7.1f'%(H_r) , ' A/m'\n",
    "print 'H_theta = ','%7.1f'%(H_theta), ' A/m'\n",
    "print 'total field =','%7.1f'%(np.sqrt(H_r**2+H_theta**2)), ' A/m'\n",
    "print 'in uT, this is =','%7.1f'%(mu_o*np.sqrt(H_r**2+H_theta**2)*1e6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Problem 4:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Knowing that $B = \\mu_o H$, work out the fundamental units of $\\mu_o$ in SI units.    Fundamental units are the units for time (seconds), distance (meters), mass (kilograms), electrical potential (volts), electrical current (amps) and so on.   Units work just like other variables in algebra; you can multiply and divide them.  So, \n",
    "\n",
    "$B$ is in tesla which boil down to  [kg A$^{-1}$ s $^{-2}$]  in fundamental units\n",
    "\n",
    "$H $ is in [A/m] which already is in fundamental units.  \n",
    "\n",
    "We know that $\\mu_o = {B\\over H}$, so the units of $\\mu_o$ are the units for $B$ divided by the units for $H$:\n",
    "\n",
    "i.e., \n",
    "\n",
    "$\n",
    "{\\mu_o  [} {{ kg \\cdot m} \\over {A^2 s^2} }].\n",
    "$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}