{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "single_cuboid-0",
   "metadata": {
    "tags": [
     "readme"
    ]
   },
   "source": [
    "# Single cuboid\n",
    "\n",
    "Calculate the field and its analytical gradient outside a uniformly magnetized cube. All lengths here are in nm, magnetization is in A/m, fields are in T, and gradients are in T/nm.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "single_cuboid-1",
   "metadata": {
    "tags": [
     "readme"
    ]
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "from microcubed import Magnet\n",
    "\n",
    "cube = Magnet(size=[100, 100, 100], center=[0, 0, 0], magnetization=[0, 0, 8e5])\n",
    "points = np.array([[0, 0, -150], [80, 0, -150]]).T\n",
    "field = cube.Bfield(points)  # (3, N): Bx, By, Bz\n",
    "gradient = cube.dBfield(points)  # (3, 3, N): derivative axis, field axis, point\n",
    "print(field)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "single_cuboid-2",
   "metadata": {},
   "source": [
    "## Verify the result\n",
    "\n",
    "These assertions also run in GitHub Actions.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "single_cuboid-3",
   "metadata": {},
   "outputs": [],
   "source": [
    "assert field.shape == (3, 2)\n",
    "assert gradient.shape == (3, 3, 2)\n",
    "assert np.isfinite(field).all() and np.isfinite(gradient).all()\n",
    "# Reflection symmetry makes transverse components vanish on the z axis.\n",
    "np.testing.assert_allclose(field[:2, 0], 0, atol=1e-14)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "single_cuboid-4",
   "metadata": {},
   "source": [
    "## Visualize the field and source geometry\n",
    "\n",
    "The outlines are XY projections of the source cuboids onto the field map,\n",
    "not intersections with the sampling plane. The white line is the\n",
    "projected material boundary. The cube is above the sampling plane.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "single_cuboid-5",
   "metadata": {
    "tags": [
     "readme"
    ]
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig, ax = cube.plot_2d(x=(-250, 250, 61), y=(-250, 250, 61), z=-150, component=\"z\")\n",
    "for boundary in cube.union_boundary(\"xy\"):\n",
    "    ax.plot(*boundary, \"w-\", linewidth=2, label=\"Material boundary\")\n",
    "ax.set(xlabel=\"x (nm)\", ylabel=\"y (nm)\", title=\"Bz (T) at z = -150 nm\", aspect=\"equal\")\n",
    "ax.legend(loc=\"upper right\", fontsize=8, facecolor=\"#555555\", labelcolor=\"white\", framealpha=0.95)\n",
    "fig.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "single_cuboid-view-6",
   "metadata": {},
   "source": [
    "## Line profiles\n",
    "\n",
    "Sample the same exterior line at `y = 0`, `z = -150 nm`.\n",
    "The upper panel shows all three field components. The lower panel shows analytical derivatives along the line; these use `dBfield()`, not finite differences."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "single_cuboid-view-7",
   "metadata": {},
   "outputs": [],
   "source": [
    "line_x = np.linspace(-250, 250, 201)\n",
    "line_points = np.vstack([line_x, np.zeros_like(line_x), np.full_like(line_x, -150)])\n",
    "line_field = cube.Bfield(line_points)\n",
    "fig, axes = plt.subplots(2, 1, figsize=(8, 6), sharex=True, layout=\"constrained\")\n",
    "for component, label in enumerate(\"xyz\"):\n",
    "    axes[0].plot(line_x, line_field[component], label=f\"B{label}\")\n",
    "axes[0].set(ylabel=\"B (T)\", title=\"Field along y = 0, z = -150 nm\")\n",
    "line_gradient = cube.dBfield(line_points)\n",
    "for component, label in enumerate(\"xyz\"):\n",
    "    axes[1].plot(line_x, line_gradient[0, component], label=f\"dB{label}/dx\")\n",
    "axes[1].set(xlabel=\"x (nm)\", ylabel=\"Gradient (T/nm)\")\n",
    "assert np.isfinite(line_gradient).all()\n",
    "assert np.isfinite(line_field).all()\n",
    "for ax in axes:\n",
    "    ax.grid(alpha=0.2)\n",
    "    ax.legend(fontsize=8, ncol=2)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "single_cuboid-view-8",
   "metadata": {},
   "source": [
    "## 3D field and source geometry\n",
    "\n",
    "The arrows sample only the exterior region below the cuboids. Arrow direction\n",
    "shows the field direction; color represents its magnitude in tesla. Arrow lengths\n",
    "are normalized for readability, so they do not encode strength. Translucent blue\n",
    "surfaces show the cuboid faces at their actual position."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "single_cuboid-view-9",
   "metadata": {},
   "outputs": [],
   "source": [
    "from mpl_toolkits.mplot3d.art3d import Poly3DCollection\n",
    "\n",
    "fig = plt.figure(figsize=(9, 6))\n",
    "ax = fig.add_subplot(projection=\"3d\")\n",
    "cube.plot_3d(\n",
    "    x=(-250, 250, 11),\n",
    "    y=(-250, 250, 9),\n",
    "    z=(-250, -75, 4),\n",
    "    max_points=396,\n",
    "    normalize=True,\n",
    "    length=25,\n",
    "    ax=ax,\n",
    ")\n",
    "for magnet in [cube]:\n",
    "    faces = magnet.corners.T[[[0, 1, 3, 2], [4, 5, 7, 6], [0, 1, 5, 4], [2, 3, 7, 6], [0, 2, 6, 4], [1, 3, 7, 5]]]\n",
    "    ax.add_collection3d(\n",
    "        Poly3DCollection(\n",
    "            faces,\n",
    "            facecolor=\"#60a5fa\",\n",
    "            edgecolor=\"#1e3a8a\",\n",
    "            linewidth=0.6,\n",
    "            alpha=0.18,\n",
    "        )\n",
    "    )\n",
    "ax.set(\n",
    "    xlim=(-250, 250),\n",
    "    ylim=(-250, 250),\n",
    "    zlim=(-250, 70),\n",
    "    xlabel=\"x (nm)\",\n",
    "    ylabel=\"y (nm)\",\n",
    "    zlabel=\"z (nm)\",\n",
    "    title=\"Exterior field and source cuboids\",\n",
    ")\n",
    "ax.set_box_aspect((500, 500, 320))\n",
    "ax.view_init(elev=22, azim=-60)\n",
    "ax.set_position([0.01, 0.06, 0.7, 0.88])\n",
    "fig.axes[-1].set_position([0.88, 0.22, 0.025, 0.56])\n",
    "fig.axes[-1].set_ylabel(\"|B| (T)\")\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}