Modelează un fluid non-vâscos în curgere folosind QUICK-PDE
Funcțiile Qiskit sunt o funcționalitate experimentală disponibilă doar utilizatorilor IBM Quantum® Premium Plan, Flex Plan și On-Prem (prin IBM Quantum Platform API) Plan. Acestea se află în stare de previzualizare și pot fi modificate.
Estimare de utilizare: 50 de minute pe un procesor Heron r2. (NOTĂ: Aceasta este doar o estimare. Timpul tău de execuție poate varia.)
Reține că timpul de execuție al acestei funcții este în general mai mare de 20 de minute, deci poate vrei să împarți acest tutorial în două secțiuni: prima, în care îl parcurgi și lansezi job-urile, și a doua, câteva ore mai târziu (oferind suficient timp pentru finalizarea job-urilor), pentru a lucra cu rezultatele.
Context
Acest tutorial urmărește să prezinte la nivel introductiv cum se folosește funcția QUICK-PDE pentru a rezolva probleme complexe de multi-fizică pe QPU-uri Heron R2 de 156Q, folosind H-DES (Hybrid Differential Equation Solver) de la ColibriTD. Algoritmul de bază este descris în articolul H-DES. Reține că acest solver poate rezolva și ecuații neliniare.
Problemele de multi-fizică — incluzând dinamica fluidelor, difuzia termică și deformarea materialelor, pentru a numi doar câteva — pot fi descrise ubic prin Ecuații cu Derivate Parțiale (EDP).
Astfel de probleme sunt extrem de relevante pentru diverse industrii și constituie o ramură importantă a matematicii aplicate. Cu toate acestea, rezolvarea EDP-urilor neliniare multivariabile cuplate cu instrumente clasice rămâne dificilă din cauza necesității unor resurse exponențial mari.
Această funcție este potrivită pentru ecuații cu complexitate și variabile în creștere și reprezintă primul pas pentru a debloca posibilități considerate odinioară inabordabile. Pentru a descrie complet o problemă modelată prin EDP, este necesar să cunoști condițiile inițiale și de frontieră. Acestea pot schimba semnificativ soluția EDP și calea de găsire a soluției.
Acest tutorial te învață cum să:
- Definești parametrii funcției de condiție inițială.
- Ajustezi numărul de qubiți (folosiți pentru a codifica funcția ecuației diferențiale), adâncimea și numărul de măsurători.
- Rulezi QUICK-PDE pentru a rezolva ecuația diferențială de bază.
Cerințe
Înainte de a începe acest tutorial, asigură-te că ai instalat următoarele:
- Qiskit SDK v2.0 sau mai nou (
pip install qiskit) - Qiskit Functions Catalog (
pip install qiskit-ibm-catalog) - Matplotlib (
pip install matplotlib) - Acces la funcția QUICK-PDE. Completează formularul pentru a solicita acces.
Configurare
Autentifică-te folosind cheia ta API și selectează funcția astfel:
# 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")
Pasul 1: Setează proprietățile problemei de rezolvat
Acest tutorial acoperă experiența utilizatorului din două perspective: problema fizică determinată de condițiile inițiale și componenta algoritmică pentru rezolvarea unui exemplu de dinamică a fluidelor pe un calculator cuantic.
Dinamica Computațională a Fluidelor (CFD) are o gamă largă de aplicații, deci este important să studiem și să rezolvăm EDP-urile de bază. O familie importantă de EDP-uri sunt ecuațiile Navier-Stokes, care constituie un sistem de ecuații diferențiale parțiale neliniare ce descriu mișcarea fluidelor. Acestea sunt extrem de relevante pentru problemele științifice și aplicațiile inginerești.
În anumite condiții, ecuațiile Navier-Stokes se reduc la ecuația Burgers', o ecuație de convecție-difuzie care descrie fenomene ce apar în dinamica fluidelor, dinamica gazelor și acustica neliniară, pentru a numi câteva, modelând sisteme disipative.
Versiunea unidimensională a ecuației depinde de două variabile: modelând dimensiunea temporală, reprezentând dimensiunea spațială. Forma generală a ecuației se numește ecuația Burgers' vâscoasă și se scrie:
unde este câmpul de viteză al fluidului la o poziție dată și momentul , iar este vâscozitatea fluidului. Vâscozitatea este o proprietate importantă a unui fluid care îi măsoară rezistența la mișcare sau deformare dependentă de rată, jucând astfel un rol crucial în determinarea dinamicii unui fluid. Când vâscozitatea fluidului este nulă ( = 0), ecuația devine o ecuație de conservare care poate dezvolta discontinuități (unde de șoc), din cauza lipsei de rezistență internă. În acest caz, ecuația se numește ecuația Burgers' non-vâscoasă și reprezintă un caz special de ecuație neliniară a undelor.
Strict vorbind, curgerea non-vâscoasă nu apare în natură, dar la modelarea fluxului aerodinamic, datorită efectului infinitezimal al transportului, utilizarea unei descrieri non-vâscoase a problemei poate fi utilă. Surprinzător, mai mult de 70% din teoria aerodinamică tratează curgeri non-vâscoase.
Acest tutorial folosește ecuația Burgers' non-vâscoasă ca exemplu CFD pentru a rezolva pe QPU-urile IBM® folosind QUICK-PDE, dată de ecuația:
Condiția inițială pentru această problemă este setată la o funcție liniară: unde și sunt constante arbitrare ce influențează forma soluției. Poți ajusta și și să observi cum influențează procesul de rezolvare și soluția.
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 ])}}
Pasul 2 (dacă este necesar): Optimizează problema pentru execuția pe hardware cuantic
În mod implicit, solverul folosește parametri informați fizic, care sunt parametrii inițiali ai circuitului pentru un număr dat de qubiți și adâncime, de la care solverul nostru va porni.
Numărul de măsurători face și el parte din parametri cu o valoare implicită, deoarece ajustarea fină a acestora este importantă.
În funcție de configurația pe care încerci să o rezolvi, parametrii algoritmului pentru a obține soluții satisfăcătoare ar putea trebui adaptați; de exemplu, poate necesita mai mulți sau mai puțini qubiți per variabilă și , în funcție de și . Următorul exemplu ajustează numărul de qubiți per funcție per variabilă, adâncimea per funcție și numărul de măsurători.
Poți vedea și cum să specifici Backend-ul și modul de execuție.
În plus, parametrii informați fizic pot conduce procesul de optimizare
într-o direcție greșită; în acest caz, îi poți dezactiva setând strategia
initialization la "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 ])}}
Pasul 3: Comparați performanțele algoritmului
Poți compara procesul de convergență al soluției noastre (HDES) din job_2 cu performanța unui algoritm și solver de rețele neuronale informate fizic (PINN) (vezi articolul și depozitul GitHub asociat).
În exemplul ieșirii din job_2 (abordare bazată pe cuantică), sunt optimizați doar 13 parametri (12 parametri de circuit plus 1 parametru de scalare) cu solverul clasic. Procesul de convergență este următorul:
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
Asta înseamnă că o pierdere sub 0.0015 poate fi atinsă după 28 de iterații, optimizând doar câțiva parametri clasici.
Acum putem compara același lucru cu soluția PINN cu configurația implicită sugerată de articol, folosind un optimizer bazat pe gradient. Echivalentul circuitului nostru cu 13 parametri de optimizat este rețeaua neuronală, care necesită cel puțin opt straturi de 20 de neuroni și, astfel, implică optimizarea a 3021 de parametri. Apoi, pierderea țintă este atinsă la Pasul 315, loss: 0.0014988397.

Acum, deoarece dorim o comparație echitabilă, ar trebui să folosim același optimizer în ambele cazuri. Numărul minim de iterații găsit pentru 12 straturi de 20 de neuroni = 4701 parametri:
(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]
1 20 5.398521572351456e-02 1.0e+00 9.98e-03 1e-02 1e-02 0:02.3
2 40 5.444650724530220e-02 1.0e+00 9.97e-03 1e-02 1e-02 0:05.1
3 60 4.447407275438309e-02 1.0e+00 9.95e-03 1e-02 1e-02 0:08.2
4 80 2.068969979882240e-02 1.0e+00 9.94e-03 1e-02 1e-02 0:11.7
6 120 1.028892211616039e-02 1.0e+00 9.91e-03 1e-02 1e-02 0:20.1
7 140 5.140972323715687e-03 1.0e+00 9.90e-03 1e-02 1e-02 0:25.4
9 180 3.811701666563749e-03 1.0e+00 9.87e-03 1e-02 1e-02 0:37.4
10 200 3.189878538250923e-03 1.0e+00 9.85e-03 1e-02 1e-02 0:44.2
12 240 2.547040116041899e-03 1.0e+00 9.83e-03 1e-02 1e-02 0:59.7
14 280 2.166548743844032e-03 1.0e+00 9.80e-03 1e-02 1e-02 1:18.0
15 300 1.783065614290535e-03 1.0e+00 9.79e-03 1e-02 1e-02 1:28.4
16 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
Poți face același lucru cu datele tale din job_2 și să trasezi o comparație cu soluția PINN.
# check the loss function and compare between the two approaches
print(job_2.logs())
Pasul 4: Folosește rezultatul
Cu soluția ta, poți acum alege ce să faci cu ea. Cele de mai jos demonstrează cum să trasezi rezultatul.
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()

Observă diferența condiției inițiale pentru al doilea rulaj și efectul acesteia asupra rezultatului:
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()

Sondaj tutorial
Te rog să acorzi un minut pentru a oferi feedback cu privire la acest tutorial. Părerile tale ne vor ajuta să îmbunătățim oferta de conținut și experiența utilizatorilor:
Note: This survey is provided by IBM Quantum and relates to the original English content. To give feedback on doQumentation's website, translations, or code execution, please open a GitHub issue.