NGSolve24

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# VTK Output for 3D (and 2D) Visualization\n",
    "\n",
    "Besides the convenient visualizations of the Netgen GUI and the NGSolve webgui, we can also use the VTK output for 3D visualization. \n",
    "\n",
    "This is particularly useful for polished visualizations of 3D simulations where you want to apply several filters to the visualization, e.g. to show the solution on a slice of the domain, or to show the solution as a surface plot and possibly combine these things. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup dependencies\n",
    "\n",
    "In this unit we will use `pyvista` to visualize the VTK output and we will need `vtk` as well. You can install these packages via pip in your terminal (or a virtual env.):\n",
    "\n",
    "```bash\n",
    "pip3 install pyvista vtk ipywidgets trame trame-vtk trame-vuetify --break-system-packages --upgrade\n",
    "```\n",
    "or execute the following cell:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!PIP_BREAK_SYSTEM_PACKAGES=1 pip3 install pyvista vtk ipywidgets trame trame-vtk trame-vuetify --upgrade"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Functionality on dummy data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to, you can skip to the next section for the visualization of a flow field."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generating the VTK output"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ngsolve import *\n",
    "from netgen.occ import unit_cube"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#VTKOutput?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh = Mesh (unit_cube.GenerateMesh(maxh=0.4))\n",
    "gfu = GridFunction (L2(mesh=mesh, order=3))\n",
    "gfu.Set(sin(4*pi*x*y*z))\n",
    "gfv = GridFunction (L2(mesh=mesh, order=3))\n",
    "gfv.Set(cos(x**2+y**2-z**2))\n",
    "\n",
    "#define VTK output for visualization:\n",
    "vtk = VTKOutput(mesh,coefs=[gfu,gfv],names=[\"gfu\",\"gfv\"],filename=\"vtkout\",subdivision=3, floatsize=\"single\",legacy=False)\n",
    "#execute the VTK output:\n",
    "vtk.Do()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Reading the VTK data to `pyvista`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pyvista\n",
    "pyvista.set_jupyter_backend('html')\n",
    "visobj = pyvista.read('vtkout.vtu')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Simple plotting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Option 1: Call `plot`: (will give no result for jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "visobj.plot(scalars=\"gfu\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Option 1b: Call `plot` with several options:  (will give no result for jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "visobj.plot(scalars=\"gfu\", show_scalar_bar=True, show_edges=False, \n",
    "            edge_color='black', cmap='jet', lighting=True, \n",
    "            n_colors=16, clim=[-1,1], \n",
    "            scalar_bar_args={'title':'u', 'height':1.5, 'width':0.5, \n",
    "                             'position_x':0.1, 'position_y':0.1})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Option 2: Use a visualization object (`pyvista.Plotter`) and add the data (possibly step by step) to it:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot = pyvista.Plotter()\n",
    "plot.add_mesh(visobj, scalars=\"gfu\", cmap=\"terrain\")\n",
    "plot.show()\n",
    "plot.export_html(\"plot1.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For jupyter book take a look [here](plot1.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Some filters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Slices: (will give no result for jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "orthoslices = visobj.slice_orthogonal()\n",
    "orthoslices.plot(scalars=\"gfu\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Iso-surfaces: (will give no result for jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contour = visobj.contour([0.25], scalars=\"gfu\", rng=[-1,1])\n",
    "print(contour)\n",
    "contour.plot(color=\"lightblue\", smooth_shading=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Combining stuff: (will give no result for jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "def translation_matrix(x, y, z):\n",
    "    return np.array([\n",
    "        [1, 0, 0, x],\n",
    "        [0, 1, 0, 0],\n",
    "        [0, 0, 1, 0],\n",
    "        [0, 0, 0, 1]])\n",
    "box_moved_left = visobj.copy().transform(translation_matrix(-1,0,0))\n",
    "orthoslices_moved_right = orthoslices.combine().transform(translation_matrix(1,0,0))\n",
    "\n",
    "plot = pyvista.Plotter(shape=(1,2))\n",
    "plot.subplot(0,0)\n",
    "plot.add_text(\"subplot1\", font_size=24)\n",
    "plot.add_mesh(box_moved_left, scalars=\"gfu\", cmap=\"jet\")\n",
    "plot.add_mesh(contour, color=\"lightblue\", opacity=0.7)\n",
    "plot.add_mesh(orthoslices_moved_right, scalars=\"gfu\", cmap=\"jet\")\n",
    "plot.subplot(0,1)\n",
    "plot.add_text(\"subplot2\", font_size=24)\n",
    "plot.add_mesh(visobj, scalars=\"gfv\")\n",
    "\n",
    "plot.show()\n",
    "plot.export_html(\"plot2.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For jupyter book take a look [here](plot2.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualizing for the flow problem\n",
    "\n",
    "We assume that you ran the [`navierstokes_hdg_chorin.ipynb` example](navierstokes_hdg_chorin.ipynb) successfully (in parallel) and have the VTK output files available."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pyvista\n",
    "pyvista.set_jupyter_backend('html')\n",
    "visobj = pyvista.read('parstokes.pvd')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "visobj = visobj.combine()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(next cell generates no output in jupyter book)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "orthoslices = visobj.slice_orthogonal().combine()\n",
    "orthoslices.plot(scalars=\"velocity\",cmap=\"jet\",clim=[0,2.25])\n",
    "orthoslices.plot(scalars=\"pressure\",cmap=\"coolwarm\",clim=[-0.06,0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot = pyvista.Plotter()\n",
    "arrows = orthoslices.glyph(scale=\"velocity\", orient=\"velocity\", tolerance=0.03)\n",
    "\n",
    "seed = pyvista.Box(bounds=(-2,1,-1,1,-1,1), level=5)\n",
    "\n",
    "strl = visobj.streamlines_from_source(\n",
    "    seed,\n",
    "    vectors=\"velocity\",\n",
    "    max_step_length=1,\n",
    "    min_step_length=0.1,\n",
    "    max_time=500,\n",
    "    max_steps=5000,\n",
    "    max_error=1e-6,\n",
    "    integration_direction=\"forward\",\n",
    ")\n",
    "\n",
    "plot.add_mesh(orthoslices, scalars=\"velocity\", cmap=\"jet\")\n",
    "#plot.add_mesh(arrows, color=\"black\")\n",
    "plot.add_mesh(strl, color=\"black\")\n",
    "#plot.add_mesh(pyvista.Cylinder(center=(0.5,0.2,0.2), direction=(0,0,1), height=0.41, radius=0.05), color=\"white\")\n",
    "plot.show()\n",
    "plot.export_html(\"plot5.html\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When considered in html export (jupyter book), take a look [here](plot5.html)"
   ]
  }
 ],
 "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.12.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}