Skip to content

Refactor colibritd-pde tutorial: Add "CFD_EULER" use case #5323

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
JulienCalistoTD wants to merge 5 commits into
base: main
Choose a base branch
from
Open

Refactor colibritd-pde tutorial: Add "CFD_EULER" use case #5323

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
d029eff
Refactor colibritd-pde tutorial: Add "CFD_EULER" use case and Updated…
JulienCalistoTD Jul 2, 2026
73df2f2
Apply suggestions from code review
abbycross Jul 7, 2026
78a98f6
spelling
abbycross Jul 7, 2026
fefb452
tox -e fix
abbycross Jul 7, 2026
a9e06c8
add version-info tag
abbycross Jul 7, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
174 changes: 130 additions & 44 deletions docs/guides/colibritd-pde.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@
"---\n",
"\n",
"\n",
"{/* cspell:ignore CMAES, Hypoelastic, edgecolor, royalblue, rstride, cstride, colibritd, xlabel, ylabel, zlabel, Jaffali */}"
"{/* cspell:ignore CMAES, Hypoelastic, edgecolor, royalblue, rstride, cstride, colibritd, xlabel, ylabel, zlabel, Jaffali, Pressureless, Colorplot, viridis, fontsize, fontweight */}"
]
},
{
Expand All @@ -26,7 +26,7 @@
},
{
"cell_type": "markdown",
"id": "b9852687",
"id": "9cd91354",
"metadata": {
"tags": [
"version-info"
Expand Down Expand Up @@ -64,36 +64,57 @@
},
{
"cell_type": "markdown",
"id": "e2817b13",
"id": "be28949c",
"metadata": {},
"source": [
"### Computational fluid dynamics\n",
"\n",
"The inviscid Burgers' equation, models flowing non-viscous fluids as follows:\n",
"**The inviscid Burgers' equation** models flowing non-viscous fluids as follows:\n",
"\n",
"$$\\frac{\\partial u}{\\partial t} + u\\frac{\\partial u}{\\partial x} = 0,$$\n",
"\n",
"$u$ represents the fluid speed field. This use-case has a temporal boundary condition: you can select the initial condition and then allow the system to relax. Currently, the only accepted initial conditions are linear functions: $ax + b$.\n",
"$u$ represents the fluid speed field. This use-case has a temporal boundary condition: you can select the initial condition and then allow the system to relax. Currently, the only accepted initial conditions are linear functions: $ax + b$. The analytical solution is:\n",
"\n",
"$$u(t, x) = \\frac{ax + b}{at + 1}.$$\n",
"\n",
"**The pressureless Euler equations** model compressible inviscid fluid flow with damping as follows:\n",
"\n",
"$$\\frac{\\partial g}{\\partial t} + u\\frac{\\partial g}{\\partial x} + g\\frac{\\partial u}{\\partial x} = 0,$$\n",
"\n",
"$$u\\frac{\\partial g}{\\partial t} + g\\frac{\\partial u}{\\partial t} + u^2\\frac{\\partial g}{\\partial x} + 2gu\\frac{\\partial u}{\\partial x} + \\frac{\\mu}{(1+t)^{\\lambda}} g u = 0,$$\n",
"\n",
"$g$ represents the density field, $u$ the velocity field and $\\mu$ a damping coefficient. In our formulation, we pose $\\lambda = 1$, so it will not be used as a parameter in the following. This use-case has temporal boundary conditions: $g(0, x) = e^{-x}$ and $u(0, x) = x$. The analytical solution is:\n",
"\n",
"$$g(t, x) = \\frac{1-\\mu}{(1+t)^{1-\\mu} - \\mu} \\exp\\!\\left(\\frac{(\\mu-1)\\, x}{(1+t)^{1-\\mu} - \\mu}\\right),$$\n",
"\n",
"$$u(t, x) = \\frac{(1-\\mu)\\, x}{\\left((1+t)^{1-\\mu} - \\mu\\right)(1+t)^{\\mu}}.$$\n",
"\n",
"The arguments for CFD's differential equations are on a fixed grid, as follows:\n",
"\n",
"- $t$ is between 0 and 0.95 with 30 sample points. $x$ is between 0 and 0.95 with a step size of 0.2375.\n",
"- $t$ is between 0 and 0.95 with 41 sample points. $x$ is between 0 and 0.95 with 41 sample points.\n",
"\n",
"### Material deformation\n",
"\n",
"This use case focuses on hypoelastic deformation with the one-dimensional tensile test, in which a bar fixed in space is pulled at its other extremity. We describe the problem as follows:\n",
"This use case focuses on **hypoelastic deformation with the one-dimensional tensile test**, in which a bar fixed in space is pulled at its other extremity. We describe the problem as follows:\n",
"\n",
"$$u' - \\frac{\\sigma}{3K} - \\frac{2}{\\sqrt{3}}\\epsilon_0\\left(\\frac{\\sigma'}{\\sigma_0\\sqrt{3}}\\right)^n = 0$$\n",
"$$u' - \\frac{\\sigma}{3K} - \\frac{2}{\\sqrt{3}}\\epsilon_0\\left(\\frac{\\sigma'}{\\sigma_0\\sqrt{3}}\\right)^n = 0,$$\n",
"\n",
"$$\\sigma' - b = 0,$$\n",
"\n",
"$K$ represents the bulk modulus of the material being stretched, $n$ the exponent of a power law, $b$ the force per unit mass, $\\epsilon_0$ the proportional stress limit, $\\sigma_0$ the proportional strain limit, $u$ the stress function, and $\\sigma$ the strain function.\n",
"$K$ represents the bulk modulus of the material being stretched, $n$ the exponent of a power law, $b$ the force per unit mass, $\\epsilon_0$ the proportional stress limit, $\\sigma_0$ the proportional strain limit, $u$ the stress function, and $\\sigma$ the strain function. The analytical solution is:\n",
"\n",
"$$\\sigma(x) = \\sigma_0 - bx,$$\n",
"\n",
"$$u(x) = -\\frac{3^{-(3+n)/2}}{2bK(1+n)\\,\\sigma_0^{n}}\\Biggl[3^{(1+n)/2}b^2\\sigma_0^n(1+n)x^2 - 2\\cdot 3^{(1+n)/2}b\\sigma_0^n(1+n)\\sigma_0 x - 12\\epsilon_0 K\\sigma_0^{1+n}$$\n",
"$$- 12b\\epsilon_0 K\\sigma_0^n x\\left(\\frac{-bx+\\sigma_0}{\\sigma_0}\\right)^n + 12\\epsilon_0 K\\sigma_0^{n+1}\\left(\\frac{-bx+\\sigma_0}{\\sigma_0}\\right)^n - 12\\epsilon_0 K \\sigma_0^{1+n}\\Biggr],$$\n",
"\n",
"where $\\sigma_0 = g(0)$ is the boundary condition on the strain at $x=0$.\n",
"\n",
"The considered bar is of unitary length. This use-case has a boundary condition for surface stress $t$, or the amount of work needed to stretch the bar.\n",
"\n",
"The arguments for MD's differential equations are on a fixed grid, as follows:\n",
"\n",
"- $x$ is between 0 and 1 with a step size of 0.04."
"- $x$ is between 0 and 1 with 30 sample points."
]
},
{
Expand All @@ -105,10 +126,11 @@
"\n",
"The following table presents statistics on various runs of our function.\n",
"\n",
"| Example | Number of qubits | Initialization | Error | Total time (min) | Runtime usage (min) |\n",
"| ---------------------------- | ---------------- | --------------------- | --------- | ---------------- | ------------------- |\n",
"| Inviscid Burgers' equation | 50 | `PHYSICALLY_INFORMED` | $10^{-2}$ | 66 | 25 |\n",
"| Hypoelastic 1D tensile test | 18 | `RANDOM` | $10^{-2}$ | 123 | 100 |"
"| Example | Number of qubits | Initialization | Error | Total time (min) | Runtime usage (min) |\n",
"| ------------------------------ | ---------------- | --------------------- | --------- | ---------------- | ------------------- |\n",
"| Inviscid Burgers' equation | 50 | `PHYSICALLY_INFORMED` | $10^{-2}$ | 66 | 25 |\n",
"| Pressureless Euler's equations | 70 | `PHYSICALLY_INFORMED` | $10^{-2}$ | 48 | 34 |\n",
"| Hypoelastic 1D tensile test | 18 | `RANDOM` | $10^{-2}$ | 123 | 100 |"
]
},
{
Expand All @@ -124,12 +146,18 @@
{
"cell_type": "code",
"execution_count": null,
"id": "c15f6a96-8c21-4a9b-944b-8ec81a4b3696",
"id": "95a715d2",
"metadata": {},
"outputs": [],
"source": [
"from qiskit_ibm_catalog import QiskitFunctionsCatalog\n",
"\n",
"catalog = QiskitFunctionsCatalog(\n",
" channel=\"ibm_cloud / ibm_quantum_platform\",\n",
" instance=\"USER_CRN / HGP\",\n",
" token=\"USER_API_KEY / IQP_API_TOKEN\",\n",
")\n",
"\n",
"catalog = QiskitFunctionsCatalog(channel=\"ibm_quantum_platform\")\n",
"\n",
"# Verify that you have access to the function\n",
Expand All @@ -139,7 +167,7 @@
{
"cell_type": "code",
"execution_count": null,
"id": "95a715d2",
"id": "4ec04623",
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -161,9 +189,9 @@
"id": "446ac943",
"metadata": {},
"source": [
"### Computational Fluid Dynamics (CFD)\n",
"### Inviscid Burgers' equation (CFD)\n",
"\n",
"When initial conditions are set to $u(0,x) = x$, the results are as follows:"
"For Burgers' equation, when initial conditions are set to $u(0,x) = x$, the results are as follows:"
]
},
{
Expand All @@ -174,7 +202,9 @@
"outputs": [],
"source": [
"# launch the simulation with initial conditions u(0,x) = a*x + b\n",
"job = quick.run(use_case=\"cfd\", physical_parameters={\"a\": 1.0, \"b\": 0.0})"
"job = quick.run(\n",
" use_case=\"CFD_BURGER\", physical_parameters={\"a\": 1.0, \"b\": 0.0}\n",
")"
]
},
{
Expand Down Expand Up @@ -208,25 +238,83 @@
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"_ = plt.figure()\n",
"ax = plt.axes(projection=\"3d\")\n",
"\n",
"# plot the solution using the 3d plotting capabilities of pyplot\n",
"t, x = np.meshgrid(solution[\"samples\"][\"t\"], solution[\"samples\"][\"x\"])\n",
"ax.plot_surface(\n",
" t,\n",
" x,\n",
" solution[\"functions\"][\"u\"],\n",
" edgecolor=\"royalblue\",\n",
" lw=0.25,\n",
" rstride=26,\n",
" cstride=26,\n",
" alpha=0.3,\n",
")\n",
"ax.scatter(t, x, solution[\"functions\"][\"u\"], marker=\".\")\n",
"ax.set(xlabel=\"t\", ylabel=\"x\", zlabel=\"u(t,x)\")\n",
"\n",
"plt.show()"
"def plot_result_3d(result):\n",
" fig = plt.figure()\n",
" ax = fig.add_subplot(projection=\"3d\")\n",
"\n",
" t, x = np.meshgrid(result[\"samples\"][\"t\"], result[\"samples\"][\"x\"])\n",
"\n",
" ax.plot_surface(\n",
" t,\n",
" x,\n",
" result[\"functions\"][\"u\"],\n",
" edgecolor=\"royalblue\",\n",
" lw=0.25,\n",
" rstride=26,\n",
" cstride=26,\n",
" alpha=0.3,\n",
" )\n",
" ax.scatter(t, x, result[\"functions\"][\"u\"], marker=\".\")\n",
" ax.set(xlabel=\"t\", ylabel=\"x\", zlabel=\"u(t,x)\")\n",
"\n",
" plt.show()\n",
"\n",
"\n",
"# Call\n",
"plot_result_3d(solution)"
]
},
{
"cell_type": "markdown",
"id": "4408d52e",
"metadata": {},
"source": [
"### Pressureless Euler's equation (CFD)\n",
"\n",
"For Euler's equation, when initial conditions are set to $g(0, x) = e^{-x}$ and $u(0, x) = x$, for a given $\\mu$ (here, $\\mu = 0.1$) and $\\lambda = 1$, the results are as follows:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e0412513",
"metadata": {},
"outputs": [],
"source": [
"# Launches the solving for an arbitrary mu\n",
"job = quick.run(use_case=\"CFD_EULER\", physical_parameters={\"mu\": 0.1})\n",
"\n",
"solution = job.result()\n",
"\n",
"\n",
"# Colorplot function\n",
"def plot_result_2d(result):\n",
" fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n",
"\n",
" configs = {\n",
" \"g\": {\"cmap\": \"viridis\", \"title\": \"g(t, x)\"},\n",
" \"u\": {\"cmap\": \"plasma\", \"title\": \"u(t, x)\"},\n",
" }\n",
"\n",
" t = result[\"samples\"][\"t\"]\n",
" x = result[\"samples\"][\"x\"]\n",
"\n",
" for ax, (field, cfg) in zip(axes, configs.items()):\n",
" v = result[\"functions\"][field]\n",
"\n",
" im = ax.contourf(t, x, v, levels=50, cmap=cfg[\"cmap\"])\n",
" fig.colorbar(im, ax=ax, label=cfg[\"title\"])\n",
"\n",
" ax.set_xlabel(\"t\")\n",
" ax.set_ylabel(\"x\")\n",
" ax.set_title(cfg[\"title\"], fontsize=13, fontweight=\"bold\")\n",
"\n",
" plt.tight_layout()\n",
" plt.show()\n",
"\n",
"\n",
"plot_result_2d(solution)"
]
},
{
Expand All @@ -246,11 +334,9 @@
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"\n",
"# select the properties of your material\n",
"# Select the properties of your material\n",
"job = quick.run(\n",
" use_case=\"md\",\n",
" use_case=\"MD\",\n",
" physical_parameters={\n",
" \"t\": 12.0,\n",
" \"K\": 100.0,\n",
Expand All @@ -261,7 +347,7 @@
" },\n",
")\n",
"\n",
"# plot the result\n",
"# Plot the result\n",
"solution = job.result()\n",
"\n",
"_ = plt.figure()\n",
Expand Down Expand Up @@ -311,12 +397,12 @@
"metadata": {},
"outputs": [],
"source": [
"job = quick.run(use_case=\"mdf\", physical_params={})\n",
"job = quick.run(use_case=\"MD\", physical_params={})\n",
"\n",
"print(job.error_message())\n",
"\n",
"\n",
"# or write a wrapper around it for a more human readable version\n",
"# A wrapper can also be used for a more human readable version\n",
"def pprint_error(job):\n",
" print(\"\".join(eval(job.error_message())[\"error\"]))\n",
"\n",
Expand Down
23 changes: 14 additions & 9 deletions docs/tutorials/colibritd-pde.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -120,7 +120,7 @@
},
{
"cell_type": "markdown",
"id": "642a4fb3",
"id": "939789bb",
"metadata": {},
"source": [
"## Step 1: Set properties of the problem to solve\n",
Expand All @@ -145,13 +145,13 @@
"representing the spatial dimension. The general form of the equation is called\n",
"the viscous Burgers' equation and reads:\n",
"\n",
"$\\frac{\\partial u}{\\partial t} + u \\frac{\\partial u}{\\partial x} = \\nu \\frac{\\partial^2 u}{\\partial^2 x},$\n",
"$$\\frac{\\partial u}{\\partial t} + u \\frac{\\partial u}{\\partial x} = \\nu \\frac{\\partial^2 u}{\\partial x^2},$$\n",
"\n",
"where $u(x,t)$ is the fluid speed field at a given position $x$ and time $t$, and $\\nu$\n",
"is the viscosity of the fluid. Viscosity is an important property of a fluid\n",
"that measures its rate-dependent resistance to motion or deformation, and thus it\n",
"plays a crucial role in the determination of the dynamics of a fluid. When the\n",
"viscosity of the fluid is null ($\\nu$ = 0), the equation becomes a conservation\n",
"viscosity of the fluid is null ($\\nu = 0$), the equation becomes a conservation\n",
"equation that can develop discontinuities (shock waves), due to the lack of its\n",
"internal resistance. In this case, the equation is called the inviscid Burgers'\n",
"equation and is a special case of nonlinear wave equation.\n",
Expand All @@ -162,15 +162,20 @@
"70% of aerodynamic theory deals with inviscid flows.\n",
"\n",
"This tutorial uses the inviscid Burgers' equation as a CFD example to solve\n",
"on IBM® QPUs using QUICK-PDE, given by the equation:\n",
"on IBM® QPUs using QUICK-PDE, given by:\n",
"\n",
"$\\frac{\\partial u}{\\partial t} + u\\frac{\\partial u}{\\partial x} = 0 $\n",
"$$\\frac{\\partial u}{\\partial t} + u\\frac{\\partial u}{\\partial x} = 0.$$\n",
"\n",
"The initial condition for this problem is set to a linear function:\n",
"$u(t=0,x) = ax + b,\\text{ with }a,b\\in\\mathbb{R}$\n",
"$u(t=0,x) = ax + b,\\text{ with }a,b\\in\\mathbb{R},$\n",
"where $a$ and $b$ are arbitrary constants influencing the shape of the\n",
"solution. You can adjust $a$ and $b$ and see how they influence the\n",
"solving process and the solution."
"solving process and the solution. For this choice of initial condition,\n",
"the equation admits the following analytical solution:\n",
"\n",
"$$u(t, x) = \\frac{ax + b}{at + 1},$$\n",
"\n",
"which holds as long as no shock has formed, as in, for $t < -1/a$ when $a < 0$."
]
},
{
Expand All @@ -181,7 +186,7 @@
"outputs": [],
"source": [
"job = quick.run(\n",
" use_case=\"cfd\",\n",
" use_case=\"CFD_BURGER\",\n",
" physical_parameters={\"a\": 1.0, \"b\": 1.0},\n",
")"
]
Expand Down Expand Up @@ -270,7 +275,7 @@
"outputs": [],
"source": [
"job_2 = quick.run(\n",
" use_case=\"cfd\",\n",
" use_case=\"CFD_BURGER\",\n",
" physical_parameters={\"a\": 0.5, \"b\": 0.25},\n",
" nb_qubits={\"u\": {\"t\": 2, \"x\": 1}},\n",
" depth={\"u\": 3},\n",
Expand Down