نمذجة سائل غير لزج متدفق باستخدام QUICK-PDE

ملاحظة

دوال Qiskit هي ميزة تجريبية متاحة فقط لمستخدمي خطة IBM Quantum® Premium Plan وFlex Plan وOn-Prem (عبر IBM Quantum Platform API) Plan. وهي في مرحلة الإصدار التجريبي وقابلة للتغيير.

تقدير الاستخدام: 50 دقيقة على معالج Heron r2. (ملاحظة: هذا تقدير فحسب. قد يختلف وقت التشغيل الفعلي.)

لاحظ أن وقت تنفيذ هذه الدالة يتجاوز عموماً 20 دقيقة، لذا قد ترغب في تقسيم هذا البرنامج التعليمي إلى قسمين: الأول تقرأ فيه المحتوى وتُطلق المهام، والثاني بعد بضع ساعات (لإتاحة وقت كافٍ لاكتمال المهام) للعمل مع نتائج المهام.

الخلفية النظرية

يهدف هذا البرنامج التعليمي إلى تعليم المستخدمين على مستوى تمهيدي كيفية استخدام دالة QUICK-PDE لحل مسائل فيزيائية متعددة التخصصات معقدة على وحدات المعالجة الكمومية 156Q Heron R2 باستخدام H-DES (حل المعادلات التفاضلية الهجين) من ColibriTD. الخوارزمية الأساسية موصوفة في ورقة H-DES. لاحظ أن هذا الحل يمكنه أيضاً حل المعادلات غير الخطية.

يمكن وصف مسائل الفيزياء المتعددة التخصصات - بما فيها ديناميكيات الموائع، وانتشار الحرارة، وتشوه المواد، على سبيل المثال - بشكل شامل من خلال المعادلات التفاضلية الجزئية (PDEs).

هذه المسائل ذات صلة بالغة الأهمية لمختلف الصناعات وتُشكّل فرعاً مهماً من الرياضيات التطبيقية. غير أن حل المعادلات التفاضلية الجزئية غير الخطية متعددة المتغيرات والمترابطة بالأدوات الكلاسيكية يظل تحدياً كبيراً نظراً للحاجة إلى موارد هائلة تتضاعف بشكل أسي.

هذه الدالة مناسبة للمعادلات ذات التعقيد والمتغيرات المتزايدة، وهي الخطوة الأولى لفتح إمكانيات كانت تُعدّ حتى وقت قريب مستعصية على الحل. لوصف مسألة مُنمذَجة بالمعادلات التفاضلية الجزئية وصفاً كاملاً، من الضروري معرفة الشروط الابتدائية وشروط الحدود. يمكن أن تُغيّر هذه الشروط تغييراً جوهرياً حل المعادلة التفاضلية الجزئية والمسار المتبع للوصول إلى ذلك الحل.

يُعلّمك هذا البرنامج التعليمي كيفية:

1. تحديد معاملات دالة الشرط الابتدائي.
2. ضبط عدد الكيوبتات (المستخدمة لترميز دالة المعادلة التفاضلية) والعمق وعدد اللقطات.
3. تشغيل QUICK-PDE لحل المعادلة التفاضلية الأساسية.

المتطلبات

قبل البدء في هذا البرنامج التعليمي، تأكد من تثبيت ما يلي:

• Qiskit SDK v2.0 أو أحدث (pip install qiskit)
• Qiskit Functions Catalog (pip install qiskit-ibm-catalog)
• Matplotlib (pip install matplotlib)
• الوصول إلى دالة QUICK-PDE. املأ النموذج لطلب الوصول.

الإعداد

قم بالمصادقة باستخدام مفتاح API الخاص بك واختر الدالة كما يلي:

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")

الخطوة 1: تعيين خصائص المسألة المراد حلها

يتناول هذا البرنامج التعليمي تجربة المستخدم من منظورين: المسألة الفيزيائية المحددة بالشروط الابتدائية، والمكوّن الخوارزمي لحل مثال على ديناميكيات الموائع على حاسوب كمومي.

لديناميكيات الموائع الحسابية (CFD) نطاق واسع من التطبيقات، ومن ثم فإنه من المهم دراسة وحل المعادلات التفاضلية الجزئية الأساسية. ومن أهم عائلات هذه المعادلات معادلات Navier-Stokes، وهي منظومة من المعادلات التفاضلية الجزئية غير الخطية التي تصف حركة الموائع. وهي ذات أهمية بالغة للمسائل العلمية والتطبيقات الهندسية.

في ظروف معينة، تُختزل معادلات Navier-Stokes إلى معادلة Burgers'، وهي معادلة حمل-انتشار تصف الظواهر الحادثة في ديناميكيات الموائع وديناميكيات الغازات والصوتيات غير الخطية، على سبيل المثال، من خلال نمذجة الأنظمة التبديدية.

يعتمد الشكل أحادي البُعد للمعادلة على متغيرين: $t \in \mathbb{R}_{\geq 0}$ يُمثّل البُعد الزمني، $x \in \mathbb{R}$ يُمثّل البُعد المكاني. يُسمى الشكل العام للمعادلة معادلة Burgers' اللزجة وصيغتها:

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

حيث $u(x,t)$ هو حقل سرعة السائل عند موضع $x$ وزمن $t$، و$\nu$ هي لزوجة السائل. اللزوجة خاصية مهمة للسائل تقيس مقاومته للحركة أو الشوه المعتمدة على معدل الإجهاد، وبالتالي فهي تؤدي دوراً محورياً في تحديد ديناميكيات السائل. عندما تكون لزوجة السائل معدومة ($\nu$ = 0)، تصبح المعادلة معادلة حفظ يمكن أن تُنشئ تقطعات (موجات صدمية)، بسبب غياب مقاومتها الداخلية. في هذه الحالة، تُسمى المعادلة معادلة Burgers' غير اللزجة وهي حالة خاصة من معادلة الموجة غير الخطية.

من الناحية الدقيقة، لا توجد تدفقات غير لزجة في الطبيعة، لكن عند نمذجة التدفق الهوائي، ونظراً للأثر المتناهي الصغر للنقل، يمكن أن يكون استخدام وصف غير لزج للمسألة مفيداً. ومما يثير الدهشة أن أكثر من 70% من نظرية الديناميكيات الهوائية تتعامل مع التدفقات غير اللزجة.

يستخدم هذا البرنامج التعليمي معادلة Burgers' غير اللزجة كمثال على CFD لحله على وحدات معالجة كمومية IBM® باستخدام QUICK-PDE، وصيغتها:

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

يُعيَّن الشرط الابتدائي لهذه المسألة كدالة خطية: $u(t=0,x) = ax + b,\text{ with }a,b\in\mathbb{R}$ حيث $a$ و$b$ ثابتان اختياريان يؤثران في شكل الحل. يمكنك ضبط $a$ و$b$ ومشاهدة تأثيرهما على عملية الحل والنتيجة.

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  ])}}

الخطوة 2 (عند الحاجة): تحسين المسألة لتنفيذها على الأجهزة الكمومية

بشكل افتراضي، يستخدم الحل معاملات مستنيرة فيزيائياً، وهي معاملات دائرة ابتدائية لعدد كيوبتات وعمق معين، يبدأ منها الحل.

اللقطات هي أيضاً جزء من المعاملات بقيمة افتراضية، إذ يكون ضبطها الدقيق أمراً مهماً.

بناءً على التهيئة التي تحاول حلها، قد تحتاج معاملات الخوارزمية للوصول إلى حلول مُرضية إلى التكيف؛ على سبيل المثال، قد تتطلب كيوبتات أكثر أو أقل لكل متغير $t$ و$x$، بحسب قيم $a$ و $b$. يضبط ما يلي عدد الكيوبتات لكل دالة لكل متغير، والعمق لكل دالة، وعدد اللقطات.

يمكنك أيضاً الاطلاع على كيفية تحديد الخادم الخلفي ووضع التنفيذ.

بالإضافة إلى ذلك، قد تُوجّه المعاملات المستنيرة فيزيائياً عملية التحسين في اتجاه خاطئ؛ في تلك الحالة، يمكنك تعطيلها بتعيين استراتيجية initialization إلى "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  ])}}

الخطوة 3: مقارنة أداء الخوارزميات

يمكنك مقارنة عملية التقارب لحلنا (HDES) الخاص بـ job_2 مع أداء خوارزمية ومحلل الشبكات العصبية المستنيرة فيزيائياً (PINN) (راجع الورقة ومستودع GitHub المرتبط بها).

في مثال مخرجات job_2 (النهج الكمومي)، يتم تحسين 13 معاملاً فحسب (12 معاملاً للدائرة بالإضافة إلى معامل قياس واحد) باستخدام المحلل الكلاسيكي. تسير عملية التقارب على النحو التالي:

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

هذا يعني أنه يمكن الوصول إلى خسارة أقل من 0.0015 بعد 28 تكراراً، مع تحسين عدد قليل فحسب من المعاملات الكلاسيكية.

يمكننا الآن مقارنة ذلك بحل PINN مع التهيئة الافتراضية المقترحة في الورقة باستخدام محسِّن يعتمد على التدرج. ما يعادل دائرتنا ذات 13 معاملاً للتحسين هو الشبكة العصبية، التي تتطلب ثماني طبقات على الأقل من 20 خلية عصبية، وبالتالي يستلزم تحسين 3021 معاملاً. بعد ذلك، يُبلَغ عن بلوغ خسارة الهدف عند الخطوة 315، الخسارة: 0.0014988397.

الآن، وبما أننا نريد إجراء مقارنة عادلة، ينبغي استخدام المحسِّن ذاته في كلتا الحالتين. أدنى عدد من التكرارات وجدناه لـ 12 طبقة من 20 خلية عصبية = 4701 معاملاً:

(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

يمكنك فعل الأمر ذاته مع بياناتك من job_2، ورسم مقارنة مع حل PINN.

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

الخطوة 4: استخدام النتيجة

بعد الحصول على حلك، يمكنك الآن اختيار ما تفعله به. يوضح ما يلي كيفية رسم النتيجة.

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()

لاحظ الفرق في الشرط الابتدائي للتشغيل الثاني وأثره على النتيجة:

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()

استطلاع البرنامج التعليمي

يُرجى تخصيص دقيقة لتقديم ملاحظاتك حول هذا البرنامج التعليمي. ستساعدنا آراؤك في تحسين محتوياتنا وتجربة المستخدم:

رابط الاستطلاع