running kharma in tutorial

scripts/tutorial/parfiles/anisotropic_conduction_ideal.par

@@ -68,5 +68,9 @@ single_precision_output = true
 variables = prims
 
 <parthenon/output1>
+file_type = rst
+dt = 5.
+
+<parthenon/output2>
 file_type = hst
-dt = 0.1
+dt = 0.1

scripts/tutorial/tutorial.ipynb

@@ -10,9 +10,7 @@
     "---\n",
     "---\n",
     "\n",
-    "## PREREQUISITES\n",
-    "\n",
-    "### Getting KHARMA\n",
+    "## Getting KHARMA\n",
     "\n",
     "KHARMA can be found on Github [here](https://github.com/AFD-Illinois/kharma).\n",
     "For this tutorial we'll be working with the `tutorial` branch which can be cloned like,\n",
@@ -21,7 +19,7 @@
     "\n",
     "---\n",
     "\n",
-    "### Building pyharm\n",
+    "## Building pyharm\n",
     "\n",
     "Additionally, we'll need `pyharm`, a Python package for analyzing GRMHD data products from HARM-based codes. The code is hosted [here](https://github.com/AFD-Illinois/pyharm), and can be obtained using git,\n",
     "\n",
@@ -49,7 +47,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "### Building KHARMA\n",
+    "## Building KHARMA\n",
     "\n",
     "KHARMA has two main dependencies which must be checked out before compiling it,\n",
     "```\n",
@@ -96,7 +94,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "### Running KHARMA\n",
+    "## Running KHARMA\n",
     "\n",
     "KHARMA can either be run in broadly two ways,\n",
     "- By directly calliing the executable: `./kharma.host -i PROBLEM_PARAMETER_FILE`. This works best when you are running on a personal computer or are on a compute node via an interactive job, and have manually set the module stack and required environment variables correctly.\n",
@@ -111,8 +109,76 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "#### Problem 1: Hot cylinder in pressure equilibrium\n",
-    "\n"
+    "### Problem 1: Pressure equilibrium in ideal MHD\n",
+    "\n",
+    "It is good practice to create a separate directory for each problem run; say `SIMULATION_DIR` is where we'll run the first problem. KHARMA will output data to the directory where the executable is called from. Let's run the following set of commands,\n",
+    "\n",
+    "```\n",
+    "cd SIMULATION_DIR\n",
+    "cp PATH_TO_KHARMA_SOURCE/scripts/tutorial/batch/rusty_gpu.sb ./\n",
+    "cp PATH_TO_KHARMA_SOURCE/scripts/tutorial/parfiles/anisotropic_conduction_ideal.par ./\n",
+    "```\n",
+    "This copies over the job submit script and the parameter file for the problem to `SIMULAITON_DIR`. Before we move on to running the problem, let's take a look at the parameter file.\n",
+    "\n",
+    "KHARMA was designed from ground up with flexibility and extensibility in mind; most parameters are therefore runtime. This means the number of parameters is [long and growing](https://github.com/AFD-Illinois/kharma/wiki/Parameters). However, on a day-to-day basis there are only a limited number of parameters that you'd need to modify depending on the problem and KHARMA sets sensible default values for the rest.\n",
+    "Parameters are divided into blocks and are specified like `block_name/parameter_name=value`. Let's take a look at the parameter file for this problem to see this in practice.\n",
+    "\n",
+    "Now, let's submit the job\n",
+    "\n",
+    "`sbatch rusty_gpu.sb -i anisotropic_conduction_ideal.par`\n",
+    "\n",
+    "In this case we explicitly ask Parthenon (via the `-d` flag in the job submit script) to create the save output (dump files, restart files, etc..) in a separate directory `dumps_kharma`. We could have pointed to the parameter file at `PATH_TO_KHARMA_SOURCE/scripts/tutorial/parfiles/` but it is cleaner this way---each simulation has its own copy of the parameter file and we avoid changing the defaults in the git repository."
+   ]
+  },
+  {
+   "attachments": {},
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "**Congratulations on running your first KHARMA simulation!**\n",
+    "\n",
+    "Let's take a look at the output. The standard output is saved to the `.out` file. When the problem is executing, you'll see a line like this (for explicitly-evolved problems at least),\n",
+    "\n",
+    "```\n",
+    "cycle=50 time=6.0457029249727713e-02 dt=1.4976930529853696e-03 zone-cycles/wsec_step=4.24e+07 wsec_total=4.25e-01 wsec_step=6.18e-03\n",
+    "Max DivB: 1.52656e-16\n",
+    "```\n",
+    "`cycle` denotes the timestep number, `time` is the simulation time in code units (e.g., for a black hole accretion problem it'd be in units of the light-crossing time $GM/c^3$), `dt` is the size of the timestep, and `zone-cycles/wsec_step` is the number of zone cycle per second.\n",
+    "\n",
+    "The simulation output is saved to `dump_kharma`. Files ending with `phdf` are the fluid \"dump\" files and stores the state of the fluid at a given time. Restart files have the `rhdf` extensions.\n",
+    "\n",
+    "It is most convenient to use `pyharm` to work with KHARMA output. At the heart of `pyharm` is the FluidDump object that caches the data in the dump file, and computes derived quantities from the cached data on-the-fly. It acts like a Python dictionary from which various quantities (related to the fluid state, grid, input parameters, etc..) can be accessed by passing the right key. The `load_dump` method is used to obtain the fluid state object. For example,"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import pyharm\n",
+    "import matplotib.pyplot as plt\n",
+    "\n",
+    "dump = pyharm.load_dump(os.path.join(dumpsdir, 'anisotropic_conduction.out0.00000.phdf'))\n",
+    "\n",
+    "fig = plt.figure()\n",
+    "ax  = fig.gca()\n",
+    "ax.pcolormesh(np.squeeze(dump['X1']), np.squeeze(dump['X2']), np.squeeze(dump['Theta']), cmap = 'viridis', shading='gouraud')\n",
+    "ax.set_xlim(0,1)\n",
+    "ax.set_ylim(0,1)\n",
+    "ax.set_xticks([0,0.25,0.5,0.75,1])\n",
+    "ax.set_xticklabels([0,0.25,0.5,0.75,1])\n",
+    "ax.set_yticks([0,0.25,0.5,0.75,1])\n",
+    "ax.set_yticklabels([0,0.25,0.5,0.75,1])\n",
+    "ax.set_xlabel('$x (GM/c^2)$')\n",
+    "ax.set_ylabel('$y (GM/c^2)$')\n",
+    "ax.set_aspect('equal')\n",
+    "divider = make_axes_locatable(ax)\n",
+    "cax = divider.append_axes(\"right\", size=\"5%\", pad=0.05)\n",
+    "cbar = plt.colorbar(temp_plot, cax=cax)\n",
+    "cbar.set_label('P$_{g}/\\\\rho$', rotation=90)\n",
+    "\n",
+    "plt.show()"
    ]
   }
  ],