Gumawa ng modelo ng umaagos na walang-lagkit na likido gamit ang QUICK-PDE

Paalala

Ang mga Qiskit Functions ay isang eksperimentong feature na available lamang sa mga IBM Quantum® Premium Plan, Flex Plan, at On-Prem (sa pamamagitan ng IBM Quantum Platform API) Plan users. Nasa preview release status ang mga ito at maaaring magbago.

Tinatayang paggamit: 50 minuto sa isang Heron r2 processor. (TANDAAN: Ito ay isang tantiya lamang. Maaaring mag-iba ang inyong runtime.)

Tandaan na ang execution time ng function na ito ay karaniwang higit sa 20 minuto, kaya maaaring nais ninyong hatiin ang tutorial na ito sa dalawang seksyon: ang una, kung saan babasahin ninyo ito at ilulunsad ang mga jobs, at ang pangalawa ilang oras mamaya (na nagbibigay ng sapat na oras para makumpleto ang mga jobs), upang magtrabaho sa mga resulta ng mga jobs.

Pinagmulan

Ang tutorial na ito ay naglalayong magturo sa introductory level kung paano gamitin ang QUICK-PDE function upang malutas ang mga kumplikadong multi-physics problems sa 156Q Heron R2 QPUs sa pamamagitan ng paggamit ng ColibriTD's H-DES (Hybrid Differential Equation Solver). Ang pinagbabatayan na algorithm ay inilalarawan sa H-DES paper. Tandaan na ang solver na ito ay maaari ring makalutas ng non-linear equations.

Ang mga multi-physics problems - kabilang ang fluids dynamics, heat diffusion, at material deformation, upang pangalanan ang ilan - ay maaaring pangkalahatang ilarawan sa pamamagitan ng Partial Differential Equations (PDEs).

Ang mga problemang tulad nito ay lubhang may kaugnayan para sa iba't ibang industriya at bumubuo ng isang mahalagang sangay ng applied mathematics. Gayunpaman, ang paglutas ng non-linear multivariate coupled PDEs gamit ang mga klasikal na tools ay nananatiling hamon dahil sa pangangailangan ng exponentially malaking dami ng resources.

Ang function na ito ay angkop para sa mga equation na may tumataas na komplikasyon at variables, at ito ay ang unang hakbang upang buksan ang mga posibilidad na dati ay itinuturing na hindi malulutas. Upang ganap na ilarawan ang isang problema na nimodelo ng mga PDEs, kinakailangan na malaman ang initial at boundary conditions. Maaaring lubhang baguhin ng mga ito ang solusyon ng PDE at ang landas upang mahanap ang kanilang solusyon.

Ituturo sa inyo ng tutorial na ito kung paano:

Tukuyin ang mga parameter ng initial condition function.
Isaayos ang qubit number (na ginagamit upang i-encode ang function ng differential equation), depth, at shot number.
Patakbuhin ang QUICK-PDE upang malutas ang pinagbabatayan na differential equation.

Mga Kinakailangan

Bago simulan ang tutorial na ito, siguraduhing mayroon kayo ng sumusunod na naka-install:

Qiskit SDK v2.0 o mas bago (pip install qiskit)
Qiskit Functions Catalog (pip install qiskit-ibm-catalog)
Matplotlib (pip install matplotlib)
Access sa QUICK-PDE function. Punan ang form upang humiling ng access.

Pag-setup

Mag-authenticate gamit ang inyong API key at piliin ang function tulad ng sumusunod:

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit-ibm-catalog

import numpy as np
import matplotlib.pyplot as plt
from qiskit_ibm_catalog import QiskitFunctionsCatalog

catalog = QiskitFunctionsCatalog(
    channel="ibm_quantum_platform",
    instance="INSTANCE_CRN",
    token="YOUR_API_KEY",  # Use the 44-character API_KEY you created and saved from the IBM Quantum Platform Home dashboard
)

quick = catalog.load("colibritd/quick-pde")

Hakbang 1: Itakda ang mga properties ng problemang sosolusyunan

Sinasaklaw ng tutorial na ito ang user experience mula sa dalawang pananaw: ang physical problem na natutukoy ng mga initial conditions, at ang algorithmic component sa paglutas ng isang fluid dynamics example sa isang quantum computer.

Ang Computational Fluid Dynamics (CFD) ay may malawak na hanay ng mga aplikasyon, at samakatuwid ay mahalaga na pag-aralan at lutasin ang pinagbabatayan na mga PDEs. Ang isang mahalagang pamilya ng mga PDEs ay ang mga Navier-Stokes equations, na isang sistema ng nonlinear partial differential equations na naglalarawan ng paggalaw ng mga likido. Lubhang may kaugnayan ang mga ito para sa mga scientific problems at engineering applications.

Sa ilalim ng ilang kundisyon, ang mga Navier-Stokes equations ay nabawasan sa Burgers' equation, isang convection–diffusion equation na naglalarawan ng mga phenomena na nangyayari sa fluid dynamics, gas dynamics, at nonlinear acoustics, upang pangalanan ang ilan, sa pamamagitan ng pagmodelo ng dissipative systems.

Ang one-dimensional version ng equation ay umaasa sa dalawang variables: $t \in \mathbb{R}_{\geq 0}$ na nagmo-modelo ng temporal dimension, $x \in \mathbb{R}$ na kumakatawan sa spatial dimension. Ang pangkalahatang anyo ng equation ay tinatawag na viscous Burgers' equation at nagsasaad:

$\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = \nu \frac{\partial^2 u}{\partial^2 x},$

kung saan ang $u(x,t)$ ay ang fluid speed field sa isang ibinigay na posisyon $x$ at oras $t$ , at ang $\nu$ ay ang viscosity ng likido. Ang viscosity ay isang mahalagang property ng isang likido na sumusukat sa rate-dependent resistance nito sa paggalaw o deformation, at samakatuwid ay gumaganap ng mahalagang papel sa pagtukoy ng dynamics ng isang likido. Kapag ang viscosity ng likido ay null ( $\nu$ = 0), ang equation ay nagiging conservation equation na maaaring makabuo ng mga discontinuities (shock waves), dahil sa kakulangan ng internal resistance nito. Sa kasong ito, ang equation ay tinatawag na inviscid Burgers' equation at isang special case ng nonlinear wave equation.

Mahigpit na sinasalita, ang mga inviscid flows ay hindi nangyayari sa kalikasan, ngunit kapag nagmo-modelo ng aerodynamic flow, dahil sa infinitesimal effect ng transport, ang paggamit ng isang inviscid description ng problema ay maaaring maging kapaki-pakinabang. Nakakagulat na higit sa 70% ng aerodynamic theory ay tumatalakay sa mga inviscid flows.

Ginagamit ng tutorial na ito ang inviscid Burgers' equation bilang isang CFD example upang lutasin sa mga IBM® QPUs gamit ang QUICK-PDE, na ibinibigay ng equation:

$\frac{\partial u}{\partial t} + u\frac{\partial u}{\partial x} = 0$

Ang initial condition para sa problemang ito ay nakatakda sa isang linear function: $u(t=0,x) = ax + b,\text{ with }a,b\in\mathbb{R}$ kung saan ang $a$ at $b$ ay arbitrary constants na nakakaimpluwensya sa hugis ng solusyon. Maaari ninyong isaayos ang $a$ at $b$ at tingnan kung paano nila naimpluwensyahan ang proseso ng paglutas at ang solusyon.

job = quick.run(
    use_case="cfd",
    physical_parameters={"a": 1.0, "b": 1.0},
)

print(job.result())

{'functions': {'u': array([[1.        , 0.96112378, 0.9230742 , 0.88616096, 0.85058445,
        0.81644741, 0.78376878, 0.75249908, 0.72253689, 0.69374562,
        0.66597013, 0.63905258, 0.61284684, 0.58723093, 0.56211691,
        0.53745752, 0.51324915, 0.48953036, 0.46637547, 0.44388257,
        0.4221554 , 0.40127848, 0.38128488, 0.36211604, 0.34357308,
        0.32525895, 0.30651089, 0.28632252, 0.26325504, 0.23533692],
       [1.2375    , 1.19267729, 1.14850734, 1.10544526, 1.06382155,
        1.02385326, 0.98565757, 0.94926734, 0.91464784, 0.88171402,
        0.85034771, 0.82041411, 0.79177677, 0.76431068, 0.73791248,
        0.71250742, 0.68805224, 0.66453346, 0.64196021, 0.62035121,
        0.59971506, 0.5800232 , 0.56117499, 0.54295419, 0.52497612,
        0.50662498, 0.48698059, 0.4647339 , 0.43809065, 0.40466247],
       [1.475     , 1.4242308 , 1.37394048, 1.32472956, 1.27705866,
        1.23125911, 1.18754636, 1.1460356 , 1.10675879, 1.06968242,
        1.03472529, 1.00177563, 0.9707067 , 0.94139043, 0.91370806,
        0.88755732, 0.86285533, 0.83953655, 0.81754494, 0.79681986,
        0.77727473, 0.75876792, 0.74106511, 0.72379234, 0.70637915,
        0.687991  , 0.66745028, 0.64314527, 0.61292625, 0.57398802],
       [1.7125    , 1.65578431, 1.59937362, 1.54401386, 1.49029576,
        1.43866495, 1.38943515, 1.34280386, 1.29886974, 1.25765082,
        1.21910288, 1.18313715, 1.14963664, 1.11847019, 1.08950364,
        1.06260722, 1.03765842, 1.01453964, 0.99312968, 0.97328851,
        0.95483439, 0.93751264, 0.92095522, 0.90463049, 0.88778219,
        0.86935702, 0.84791997, 0.82155665, 0.78776186, 0.74331358],
       [1.95      , 1.88733782, 1.82480676, 1.76329816, 1.70353287,
        1.6460708 , 1.59132394, 1.53957212, 1.49098069, 1.44561922,
        1.40348046, 1.36449867, 1.32856657, 1.29554994, 1.26529921,
        1.23765712, 1.21246152, 1.18954273, 1.16871442, 1.14975716,
        1.13239406, 1.11625736, 1.10084533, 1.08546864, 1.06918523,
        1.05072304, 1.02838966, 0.99996803, 0.96259746, 0.91263913]])}, 'samples': {'t': array([0.        , 0.03275862, 0.06551724, 0.09827586, 0.13103448,
       0.1637931 , 0.19655172, 0.22931034, 0.26206897, 0.29482759,
       0.32758621, 0.36034483, 0.39310345, 0.42586207, 0.45862069,
       0.49137931, 0.52413793, 0.55689655, 0.58965517, 0.62241379,
       0.65517241, 0.68793103, 0.72068966, 0.75344828, 0.7862069 ,
       0.81896552, 0.85172414, 0.88448276, 0.91724138, 0.95      ]), 'x': array([0.    , 0.2375, 0.475 , 0.7125, 0.95  ])}}

Hakbang 2 (kung kinakailangan): I-optimize ang problema para sa quantum hardware execution

Sa default, ginagamit ng solver ang physically-informed parameters, na mga initial circuit parameters para sa isang ibinigay na qubit number at depth mula sa kung saan magsisimula ang aming solver.

Ang mga shots ay bahagi rin ng mga parameters na may default value, dahil mahalaga ang fine-tuning ng mga ito.

Depende sa configuration na sinusubukan ninyong lutasin, maaaring kailangang isaayos ang parameters ng algorithm upang makamit ang kasiya-siyang mga solusyon; halimbawa, maaaring mangailangan ito ng mas marami o mas kaunting qubits bawat variable $t$ at $x$ , depende sa $a$ at $b$ . Iniaayos ng sumusunod ang bilang ng mga qubits bawat function bawat variable, ang depth bawat function, at ang bilang ng mga shots.

Makikita rin ninyo kung paano tukuyin ang backend at ang execution mode.

Bilang karagdagan, ang physically-informed parameters ay maaaring mag-steer ng optimization process sa maling direksyon; sa kasong iyan, maaari ninyong i-disable ito sa pamamagitan ng pagtakda ng initialization strategy sa "RANDOM".

job_2 = quick.run(
    use_case="cfd",
    physical_parameters={"a": 0.5, "b": 0.25},
    nb_qubits={"u": {"t": 2, "x": 1}},
    depth={"u": 3},
    shots=[500, 2500, 5000, 10000],
    initialization="RANDOM",
    backend="ibm_kingston",
    mode="session",
)

print(job_2.result())

{'functions': {'u': array([[0.25      , 0.24856543, 0.24687708, 0.2449444 , 0.24277686,
        0.24038389, 0.23777496, 0.23495952, 0.23194702, 0.22874691,
        0.22536866, 0.22182171, 0.21811551, 0.21425952, 0.2102632 ,
        0.20613599, 0.20188736, 0.19752675, 0.19306361, 0.18850741,
        0.18386759, 0.1791536 , 0.17437491, 0.16954096, 0.16466122,
        0.15974512, 0.15480213, 0.1498417 , 0.14487328, 0.13990632],
       [0.36875   , 0.36681313, 0.36457201, 0.36203594, 0.35921422,
        0.35611615, 0.35275103, 0.34912817, 0.34525687, 0.34114643,
        0.33680614, 0.33224532, 0.32747327, 0.32249928, 0.31733266,
        0.31198271, 0.30645873, 0.30077002, 0.29492589, 0.28893564,
        0.28280857, 0.27655397, 0.27018116, 0.26369944, 0.2571181 ,
        0.25044645, 0.24369378, 0.23686941, 0.22998264, 0.22304275],
       [0.4875    , 0.48506084, 0.48226695, 0.47912748, 0.47565158,
        0.47184841, 0.46772711, 0.46329683, 0.45856672, 0.45354594,
        0.44824363, 0.44266894, 0.43683103, 0.43073904, 0.42440212,
        0.41782942, 0.4110301 , 0.4040133 , 0.39678818, 0.38936388,
        0.38174955, 0.37395435, 0.36598742, 0.35785791, 0.34957498,
        0.34114777, 0.33258544, 0.32389713, 0.315092  , 0.30617919],
       [0.60625   , 0.60330854, 0.59996188, 0.59621902, 0.59208895,
        0.58758067, 0.58270318, 0.57746549, 0.57187658, 0.56594545,
        0.55968112, 0.55309256, 0.54618879, 0.53897879, 0.53147158,
        0.52367614, 0.51560147, 0.50725658, 0.49865046, 0.48979211,
        0.48069053, 0.47135472, 0.46179367, 0.45201638, 0.44203186,
        0.4318491 , 0.42147709, 0.41092485, 0.40020136, 0.38931562],
       [0.725     , 0.72155625, 0.71765682, 0.71331056, 0.70852631,
        0.70331293, 0.69767926, 0.69163414, 0.68518643, 0.67834497,
        0.6711186 , 0.66351618, 0.65554655, 0.64721855, 0.63854104,
        0.62952285, 0.62017284, 0.61049986, 0.60051274, 0.59022035,
        0.57963151, 0.56875509, 0.55759992, 0.54617486, 0.53448874,
        0.52255042, 0.51036875, 0.49795257, 0.48531072, 0.47245205]])}, 'samples': {'t': array([0.        , 0.03275862, 0.06551724, 0.09827586, 0.13103448,
       0.1637931 , 0.19655172, 0.22931034, 0.26206897, 0.29482759,
       0.32758621, 0.36034483, 0.39310345, 0.42586207, 0.45862069,
       0.49137931, 0.52413793, 0.55689655, 0.58965517, 0.62241379,
       0.65517241, 0.68793103, 0.72068966, 0.75344828, 0.7862069 ,
       0.81896552, 0.85172414, 0.88448276, 0.91724138, 0.95      ]), 'x': array([0.    , 0.2375, 0.475 , 0.7125, 0.95  ])}}

Hakbang 3: Ikumpara ang mga pagganap ng algorithm

Maaari ninyong ikumpara ang convergence process ng aming solusyon (HDES) ng job_2 sa performance ng isang physics-informed neural networks (PINN) algorithm at solver (tingnan ang paper at ang kaugnay na GitHub repository).

Sa halimbawa ng output ng job_2 (quantum-based approach), tanging 13 parameters (12 circuit parameters plus 1 scaling parameter) ay na-optimize gamit ang classical solver. Ang convergence process ay sumusunod:

optimizers:
   CMA: {'ftarget': np.float64(0.1), 'verb_disp': 10, 'maxiter': 100}
   CMA: {'ftarget': np.float64(0.005), 'verb_disp': 10, 'maxiter': 20}
   CMA: {'ftarget': np.float64(0.0025), 'verb_disp': 10, 'maxiter': 30}
   CMA: {'ftarget': np.float64(0.0005), 'verb_disp': 10, 'maxiter': 10}

500 shots
================== CMA =================
option:  {'ftarget': np.float64(0.1), 'verb_disp': 10, 'maxiter': 100}
0/100, loss: 0.02456641

1000 shots
================== CMA =================
option:  {'ftarget': np.float64(0.005), 'verb_disp': 10, 'maxiter': 20}
0/20, loss: 0.03641833
1/20, loss: 0.02461719
2/20, loss: 0.0283689
3/20, loss: 0.009898383
4/20, loss: 0.04454522
5/20, loss: 0.007019971
6/20, loss: 0.00811147
7/20, loss: 0.01592619
8/20, loss: 0.00764708
9/20, loss: 0.01401516
10/20, loss: 0.01767467
11/20, loss: 0.01220387

5000 shots
================== CMA =================
option:  {'ftarget': np.float64(0.0025), 'verb_disp': 10, 'maxiter': 30}
0/30, loss: 0.01024792
1/30, loss: 0.004343748
2/30, loss: 0.01450951
3/30, loss: 0.008591284
4/30, loss: 0.00266414
5/30, loss: 0.007923613
6/30, loss: 0.02023853
7/30, loss: 0.01031438
8/30, loss: 0.009513116
9/30, loss: 0.008132266
10/30, loss: 0.005787766
11/30, loss: 0.00390582

10000 shots
================== CMA =================
option:  {'ftarget': np.float64(0.0005), 'verb_disp': 10, 'maxiter': 10}
0/10, loss: 0.002386168
1/10, loss: 0.004024823
2/10, loss: 0.001311999
3/10, loss: 0.003433991
4/10, loss: 0.002339664
5/10, loss: 0.002978438
6/10, loss: 0.005458391
7/10, loss: 0.002026701
8/10, loss: 0.00207467
9/10, loss: 0.001947627
final_loss: 0.00151994463476429

Nangangahulugan ito na ang isang loss na mas mababa sa 0.0015 ay maaabot pagkatapos ng 28 iterations, at sa pag-optimize lamang ng ilang classical parameters.

Ngayon ay maaari nating ikumpara ang pareho sa PINN solution gamit ang default configuration na iminungkahi ng paper gamit ang isang gradient-based optimizer. Ang katumbas ng aming circuit na may 13 parameters na i-optimize ay ang neural network, na nangangailangan ng hindi bababa sa walong layers ng 20 neurons, at samakatuwid ay nagsasangkot ng pag-optimize ng 3021 parameters. Pagkatapos, ang target loss ay naabot sa Step 315, loss: 0.0014988397.

Graph na nagpapakita ng PINN data kumpara sa HDES-Qiskit function.

Ngayon, dahil nais nating gumawa ng patas na pagkukumpara, dapat nating gamitin ang parehong optimizer sa parehong kaso. Ang pinakamababang bilang ng iterations na nahanap namin para sa 12 layers ng 20 neurons = 4701 parameters:

(10_w,20)-aCMA-ES (mu_w=5.9,w_1=27%) in dimension 4701 (seed=351961)
Iterat #Fevals   function value  axis ratio  sigma  min&max std  t[m:s]
   20 5.398521572351456e-02 1.0e+00 9.98e-03  1e-02  1e-02 0:02.3
   40 5.444650724530220e-02 1.0e+00 9.97e-03  1e-02  1e-02 0:05.1
   60 4.447407275438309e-02 1.0e+00 9.95e-03  1e-02  1e-02 0:08.2
   80 2.068969979882240e-02 1.0e+00 9.94e-03  1e-02  1e-02 0:11.7
  120 1.028892211616039e-02 1.0e+00 9.91e-03  1e-02  1e-02 0:20.1
  140 5.140972323715687e-03 1.0e+00 9.90e-03  1e-02  1e-02 0:25.4
  180 3.811701666563749e-03 1.0e+00 9.87e-03  1e-02  1e-02 0:37.4
  200 3.189878538250923e-03 1.0e+00 9.85e-03  1e-02  1e-02 0:44.2
  240 2.547040116041899e-03 1.0e+00 9.83e-03  1e-02  1e-02 0:59.7
  280 2.166548743844032e-03 1.0e+00 9.80e-03  1e-02  1e-02 1:18.0
  300 1.783065614290535e-03 1.0e+00 9.79e-03  1e-02  1e-02 1:28.4
  320 2.045844215899706e-03 1.0e+00 9.78e-03  1e-02  1e-02 1:39.8
Stopping early: loss 0.001405 <= target 0.0015
CMA-ES finished. Best loss: 0.001404788694344461

Maaari ninyong gawin ang pareho gamit ang inyong data mula sa job_2, at mag-plot ng pagkukumpara sa PINN solution.

# check the loss function and compare between the two approaches
print(job_2.logs())

Hakbang 4: Gamitin ang resulta

Sa inyong solusyon, maaari na ninyong piliin kung ano ang gagawin dito. Ipinapakita ng sumusunod kung paano i-plot ang resulta.

solution = job.result()

# Plot the solution of the second simulation job_2
_ = plt.figure()
ax = plt.axes(projection="3d")

# plot the solution using the 3d plotting capabilities of pyplot
t, x = np.meshgrid(solution["samples"]["t"], solution["samples"]["x"])
ax.plot_surface(
    t,
    x,
    solution["functions"]["u"],
    edgecolor="royalblue",
    lw=0.25,
    rstride=26,
    cstride=26,
    alpha=0.3,
)
ax.scatter(t, x, solution, marker=".")
ax.set(xlabel="t", ylabel="x", zlabel="u(t,x)")

plt.show()

Output ng nakaraang code cell

Pansinin ang pagkakaiba ng initial condition para sa pangalawang run at ang epekto nito sa resulta:

solution_2 = job_2.result()

# Plot the solution of the second simulation job_2
_ = plt.figure()
ax = plt.axes(projection="3d")

# plot the solution using the 3d plotting capabilities of pyplot
t, x = np.meshgrid(solution_2["samples"]["t"], solution_2["samples"]["x"])
ax.plot_surface(
    t,
    x,
    solution_2["functions"]["u"],
    edgecolor="royalblue",
    lw=0.25,
    rstride=26,
    cstride=26,
    alpha=0.3,
)
ax.scatter(t, x, solution_2, marker=".")
ax.set(xlabel="t", ylabel="x", zlabel="u(t,x)")

plt.show()

Output ng nakaraang code cell

Tutorial survey

Pakiusap po, maglaan ng isang minuto upang magbigay ng feedback sa tutorial na ito. Tutulungan kami ng inyong mga pananaw na pahusayin ang aming mga content offerings at user experience:

Link sa survey

Pinagmulan​

Mga Kinakailangan​

Pag-setup​

Hakbang 1: Itakda ang mga properties ng problemang sosolusyunan​

Hakbang 2 (kung kinakailangan): I-optimize ang problema para sa quantum hardware execution​

Hakbang 3: Ikumpara ang mga pagganap ng algorithm​

Hakbang 4: Gamitin ang resulta​

Tutorial survey​