{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 6 - Mapped domains with polar singularity\n",
    "\n",
    "This tutorial provides access to `Struphy domains` which can be defined from analytical formulas or through third-party software (VMEC, GVEC, DESC, etc.).\n",
    "\n",
    "Of particular interest are maps that have a polar singularity (\"magnetic axis\"). Such singularities are challenging numerically; in Struphy two solutions are possible:\n",
    "\n",
    "1. Use `polar splines` when setting up the de Rham sequence (=more expensive)\n",
    "2. Cut out a small hole of the domain around the singularity (=cheap, but problematic for particles).\n",
    "\n",
    "We shall focus on the second possibility in this notebook.\n",
    "\n",
    "## HollowCylinder\n",
    "\n",
    "Let us create the domain `HollowCylinder` with default parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from struphy import domains\n",
    "\n",
    "domain = domains.HollowCylinder()\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A 3D view of the domain can be displayed with"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain.show_3d()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The default parameters of `HollowCylinder` are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for key, val in domain.params.items():\n",
    "    print(key, \"=\", val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use the piece-of-cake parameter `poc` to use only a part of the disk and complete it periodically:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain_poc = domains.HollowCylinder(poc=3)\n",
    "domain_poc.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some relevant domain attributes are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(domain.kind_map)\n",
    "print(domain.pole)\n",
    "print(domain.periodic_eta3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The domain methods are also quite important:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for attr in dir(domain):\n",
    "    if callable(getattr(domain, attr)) and \"__\" not in attr and attr[0] != \"_\":\n",
    "        print(attr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aside from these, the domain object itself is callable: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(domain.__call__)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us change the size of the hole around the pole:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.HollowCylinder(a1=0.05)\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that if we set the inner radius `a1` to zero, the attribute `domain.pole` becomes `True`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.HollowCylinder(a1=0.0)\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(domain.kind_map)\n",
    "print(domain.pole)\n",
    "print(domain.periodic_eta3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## HollowTorus\n",
    "\n",
    "Let us create the domain [HollowTorus](https://struphy.pages.mpcdf.de/struphy/sections/domains.html#struphy.geometry.domains.HollowTorus) with default parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.HollowTorus()\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the attribute `periodic_eta3` is `True` for this mapping:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(domain.kind_map)\n",
    "print(domain.pole)\n",
    "print(domain.periodic_eta3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The default parameters are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for key, val in domain.params.items():\n",
    "    print(key, \"=\", val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us change the size of the hole around the pole, the poloidal angle parametrization (`sfl`)  and the `tor_period`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.HollowTorus(a1=0.05, sfl=True, tor_period=1)\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tokamak\n",
    "\n",
    "[Tokamak](https://struphy.pages.mpcdf.de/struphy/sections/domains.html#struphy.geometry.domains.Tokamak) is the class for mappings for Tokamak MHD equilibria constructed via [field-line tracing](https://struphy.pages.mpcdf.de/struphy/sections/domains.html#struphy.geometry.utilities.field_line_tracing) of a poloidal flux function $\\psi$.\n",
    "\n",
    "Let us create a Tokamak with default parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.Tokamak()\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain.show_3d()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The default parameters are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for key, val in domain.params.items():\n",
    "    if \"cx\" not in key and \"cy\" not in key:\n",
    "        print(key, \"=\", val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The ``Tokamak`` domain is always related to an [AxisymmMHDequilibrium](https://struphy.pages.mpcdf.de/struphy/sections/mhd_equils.html#struphy.fields_background.mhd_equil.base.AxisymmMHDequilibrium), which provides the flux function $\\psi$. In the default parameters this is [AdhocTorus](https://struphy.pages.mpcdf.de/struphy/sections/mhd_equils.html#struphy.fields_background.mhd_equil.equils.AdhocTorus). Instead, we could also look at the default [EQDSKequilibrium](https://struphy.pages.mpcdf.de/struphy/sections/mhd_equils.html#struphy.fields_background.mhd_equil.equils.EQDSKequilibrium):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from struphy import equils\n",
    "\n",
    "mhd_eq = equils.EQDSKequilibrium()\n",
    "\n",
    "domain = domains.Tokamak(equilibrium=mhd_eq)\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us shrink the hole:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.Tokamak(equilibrium=mhd_eq, psi_shifts=[0.2, 2])\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Stellarator mappings\n",
    "\n",
    "Struphy can read data produced by\n",
    "\n",
    "* [GVEC equilibrium code](https://gitlab.mpcdf.mpg.de/gvec-group/gvec)\n",
    "* [DESC equilibrium code](https://desc-docs.readthedocs.io/en/latest/index.html)\n",
    "\n",
    "### GVEC interface\n",
    "\n",
    "The interface is the class [GVECunit](https://struphy.pages.mpcdf.de/struphy/sections/domains.html?highlight=gvec#struphy.geometry.domains.GVECunit).\n",
    "\n",
    "Let us create an instance with default parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain = domains.GVECunit()\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The default parameters are:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for key, val in domain.params.items():\n",
    "    if \"cx\" not in key and \"cy\" not in key and \"cz\" not in key:\n",
    "        print(key, \"=\", val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us put a domain hole around the magnetic axis and use the whole Stellarator (`use_nfp=False`). The parameters must be passed through the [GVECequilibrium](https://struphy.pages.mpcdf.de/struphy/sections/mhd_equils.html#struphy.fields_background.mhd_equil.equils.GVECequilibrium):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gvec_equil = equils.GVECequilibrium(rmin=0.1, use_nfp=False)\n",
    "domain = domains.GVECunit(gvec_equil)\n",
    "domain.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### DESC interface\n",
    "\n",
    "The interface is the class [DESCunit](https://struphy.pages.mpcdf.de/struphy/sections/domains.html?highlight=gvec#struphy.geometry.domains.DESCunit). The parameters must be passed through the [DESCequilibrium](https://struphy.pages.mpcdf.de/struphy/sections/mhd_equils.html#struphy.fields_background.mhd_equil.equils.DESCequilibrium):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "desc_equil = equils.DESCequilibrium(use_nfp=False)\n",
    "domain = domains.DESCunit(desc_equil)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "domain.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "env (3.12.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.12.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}