البحث في الموقع

وصلنا الآن إلى نهاية سلسلة مقالات NumPy، ونأمل أن تكون أدركت أهمية مكتبة NumPy وأنها مكتبة قوية ومتعددة الاستخدامات. وبنفس الوقت تذكر أن لغة بايثون هي لغة قوية جدًا، وقد تكون في بعض الحالات المحددة أقوى من NumPy. لنفكر على سبيل المثال، في تمرين مثير للاهتمام على النحو التالي: جُمعت إجابات مختلفة من المجتمع، وأسفرت الحلول المقترحة عن نتائج مفاجئة. لكن لنبدأ بكتابة التعليمات البرمجية بلغة بايثون لحل هذا التمرين: def solution_1(): # Brute force # 14641 (=11*11*11*11) iterations & tests Z = [] for i in range(11): for j in range(11): for k in range(11): for l in range(11): if i+j+k+l == 10: Z.append((i,j,k,l)) return Z هذا الحل هو الأبطأ لأنه يتطلب 4 حلقات، كما أنه يختبر جميع المجموعات المختلفة والبالغ عددها 11641 والمكونة من 4 أعداد صحيحة بين 0 و 10 للاحتفاظ فقط بالمجموعات التي يكون مجموعها مساوٍ إلى 10. يمكننا طبعًا التخلص من الحلقات باستخدام أداة itertools، لكن التعليمات تظل بطيئة: import itertools as it def solution_2(): # Itertools # 14641 (=11*11*11*11) iterations & tests return [(i,j,k,l) for i,j,k,l in it.product(range(11),repeat=4) if i+j+k+l == 10] كان أحد أفضل الحلول المقترحة يستفيد من إمكانية الحصول على حلقات مخبأة ذكية تتيح لنا بناء كل صف tuple مباشرةً، دون أي اختبار كما هو موضح أدناه: def solution_3(): return [(a, b, c, (10 - a - b - c)) for a in range(11) for b in range(11 - a) for c in range(11 - a - b)] يستخدم كاتب هذا التمرين أفضل حل بالاعتماد على NumPy استراتيجية مختلفة مع مجموعة مقيدة من الاختبارات: def solution_4(): X123 = np.indices((11,11,11)).reshape(3,11*11*11) X4 = 10 - X123.sum(axis=0) return np.vstack((X123, X4)).T[X4 > -1] إذا قيّمنا مدة تنفيذ هذه التوابع نحصل على: >>> timeit("solution_1()", globals()) 100 loops, best of 3: 1.9 msec per loop >>> timeit("solution_2()", globals()) 100 loops, best of 3: 1.67 msec per loop >>> timeit("solution_3()", globals()) 1000 loops, best of 3: 60.4 usec per loop >>> timeit("solution_4()", globals()) 1000 loops, best of 3: 54.4 usec per loop كما تلاحظ حل Numpy هو الأسرع ولكن الحل باستخدام بايثون قابل للمقارنة. سنحاول الآن إضافة تعديل صغير على حل بايثون: def solution_3_bis(): return ((a, b, c, (10 - a - b - c)) for a in range(11) for b in range(11 - a) for c in range(11 - a - b)) وسنحصل على النتيجة: >>> timeit("solution_3_bis()", globals()) 10000 loops, best of 3: 0.643 usec per loop اكتسبنا هنا عاملًا قدره 100 فقط عن طريق استبدال قوسين مربعين بأقواس عادية. كيف يعقل ذلك؟ يمكن العثور على التفسير من خلال النظر في نوع الكائن المُعاد: >>> print(type(solution_3())) <class 'list'> >>> print(type(solution_3_bis())) <class 'generator'> يعيد التابع ()solution_3_bis مولدًا generator يمكن استخدامه لإنشاء القائمة الكاملة أو للتكرار على جميع العناصر المختلفة. على أي حال، تأتي عملية التسريع الضخمة من عدم الاستقرار في القائمة الكاملة، وبالتالي من المهم أن نتساءل عما إذا كنت بحاجة إلى مثيل فعلي لنتيجتك أو إذا كان المولد البسيط ينفذ المهمة. مكتبة Numpy وأخواتها هناك العديد من حزم بايثون الأخرى إلى جانب Numpy، والتي تستحق النظر لأنها تتناول أصناف مشابهة ولكن مختلفة من المشكلات باستخدام التقنيات المختلفة (التجميع، الجهاز الافتراضي، وحدة معالجة الرسوميات GPU، الضغط، إلخ). قد تكون حزمة واحدة أفضل من الأخرى اعتمادًا على مشكلتك المحددة وأجهزتك الخاصة. دعونا نوضح استخدامهم من خلال مثال بسيط جدًا، إذ نريد حساب تعبير بناءً على متجهين عددين float: import numpy as np a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = 2*a + 3*b حزمة NumExpr توفر حزمة numexpr إجراءات للتقييم السريع لتعبيرات المصفوفة من ناحية العناصر باستخدام جهاز افتراضي يعتمد على المتجهات vector-based، وهي مشابهة لحزمة SciPy، لكنها لا تتطلب خطوة ترجمة منفصلة لشيفرات C أو ++C. import numpy as np import numexpr as ne a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = ne.evaluate("2*a + 3*b") مصرف Cython Cython هو مصرّف ثابت محسن لكل من لغة برمجة بايثون ولغة كايثون Cython الموسعة (استنادًا إلى Pyrex). يجعل مصرّف كايثون كتابة امتدادات C لبايثون سهلة مثل بايثون نفسها. import numpy as np def evaluate(np.ndarray a, np.ndarray b): cdef int i cdef np.ndarray c = np.zeros_like(a) for i in range(a.size): c[i] = 2*a[i] + 3*b[i] return c a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = evaluate(a, b) مصرف Numba يمنحك مصرّف Numba القدرة على تسريع تطبيقاتك من خلال وظائف عالية الأداء مكتوبة مباشرةً بلغة بايثون. مع بعض التعليقات التوضيحية يمكن تجميع تعليمات بايثون array-oriented و math-heavy في الوقت المناسب لتعليمات الآلة الأصلية، على غرار C و ++C وفورتران Fortran، دون الحاجة إلى تبديل اللغات أو مترجم بايثون. from numba import jit import numpy as np @jit def evaluate(np.ndarray a, np.ndarray b): c = np.zeros_like(a) for i in range(a.size): c[i] = 2*a[i] + 3*b[i] return c a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = evaluate(a, b) مكتبة Theano مكتبة Theano هي مكتبة بايثون تتيح لك تحديد التعبيرات الرياضية التي تتضمن مصفوفات متعددة الأبعاد وتحسينها وتقييمها بكفاءة. تتميز Theano بالتكامل الوثيق مع الاستخدام غير المستقر والشفاف لوحدة معالجة الرسوميات GPU ، والتمايز الرمزي symbolic differentiation الفعال، وتحسينات السرعة والثبات، وإنشاء شيفرة C الديناميكية، واختبار الوحدة الشامل والتحقق الذاتي. import numpy as np import theano.tensor as T x = T.fvector('x') y = T.fvector('y') z = 2*x + 3*y f = function([x, y], z) a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = f(a, b) PyCUDA يتيح لك PyCUDA الوصول إلى واجهة برمجة تطبيقات الحساب المتوازي CUDA الخاصة بشركة Nvidia من بايثون. import numpy as np import pycuda.autoinit import pycuda.driver as drv from pycuda.compiler import SourceModule mod = SourceModule(""" __global__ void evaluate(float *c, float *a, float *b) { const int i = threadIdx.x; c[i] = 2*a[i] + 3*b[i]; } """) evaluate = mod.get_function("evaluate") a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = np.zeros_like(a) evaluate(drv.Out(c), drv.In(a), drv.In(b), block=(400,1,1), grid=(1,1)) مكتبة PyOpenCL تتيح لك PyOpenCL الوصول إلى وحدات معالجة الرسومات وأجهزة الحوسبة المتوازية الأخرى من بايثون. import numpy as np import pyopencl as cl a = np.random.uniform(0, 1, 1000).astype(np.float32) b = np.random.uniform(0, 1, 1000).astype(np.float32) c = np.empty_like(a) ctx = cl.create_some_context() queue = cl.CommandQueue(ctx) mf = cl.mem_flags gpu_a = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=a) gpu_b = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=b) evaluate = cl.Program(ctx, """ __kernel void evaluate(__global const float *gpu_a; __global const float *gpu_b; __global float *gpu_c) { int gid = get_global_id(0); gpu_c[gid] = 2*gpu_a[gid] + 3*gpu_b[gid]; } """).build() gpu_c = cl.Buffer(ctx, mf.WRITE_ONLY, a.nbytes) evaluate.evaluate(queue, a.shape, None, gpu_a, gpu_b, gpu_c) cl.enqueue_copy(queue, c, gpu_c) Scipy وأخواتها إذا كان هناك عدة حزم إضافية لمكتبة Numpy، فهناك الملايين من الحزم الإضافية لـ scipy؛ إذ ربما يكون لكل مجال من مجالات العلوم مجموعته الخاصة، ومعظم الأمثلة التي كنا ندرسها حتى الآن كان من الممكن حلها في استدعائين أو ثلاثة استدعاءات باستخدام الحزمة ذات الصلة. لكن طبعًا لم يكن هذا هو الهدف. تعد برمجة الوظائف بنفسك عمومًا تمرينًا جيدًا إذا كان لديك بعض وقت الفراغ، وتتمثل أكبر صعوبة في هذه المرحلة في العثور على هذه الحزم ذات الصلة. فيما يلي قائمة قصيرة جدًا من الحزم التي جرى صيانتها واختبارها جيدًا والتي قد تبسط حياتك العلمية (اعتمادًا على مجالك). هناك طبعًا الكثير من الخيارات وذلك حسب احتياجاتك الخاصة، ومن المحتمل ألا تكون مضطرًا إلى برمجة كل شيء بنفسك. مكتبة scikit-learn scikit-learn هي مكتبة تعلم آلي مجانية للغة برمجة بايثون، وتتميز بخوارزميات تصنيف وانحدار وتجميع مختلفة بما في ذلك دعم خوارزمية أجهزة المتجهات vector machines، والغابات العشوائية random forests وتعزيز التدرج gradient boosting و k-means و DBSCAN، وهذه الحزمة مصممة للتفاعل مع مكتبات بايثون الرقمية والعلمية numpy و SciPy. مكتبة scikit-image scikit-image عبارة عن حزمة بايثون مخصصة لمعالجة الصور، واستخدام المصفوفات المتداخلة مثل كائنات صور. يصف هذا المقال كيفية استخدام scikit-image في مهام معالجة الصور المختلفة، ويركز على الارتباط بوحدات بايثون النمطية العلمية الأخرى مثل numpy و SciPy. مكتبة SymPy SymPy هي مكتبة بايثون للرياضيات الرمزية، وتهدف إلى أن تصبح كل ميزات النظام الحاسوبي جبرية -أو اختصارًا CAS- مع الحفاظ على الشيفرة بسيطة قدر الإمكان حتى تكون مفهومة وقابلةً للتوسّع بسهولة. كُتبت كامل مكتبة SymPy بلغة بايثون. حزمة Astropy مشروع Astropy هو جهد مجتمعي لتطوير حزمة أساسية واحدة لعلم الفلك astronomy في بايثون وتعزيز إمكانية التشغيل البيني بين حزم علم الفلك في بايثون. حزمة Cartopy Cartopy هي حزمة بايثون مصممة لتسهيل رسم الخرائط لتحليل البيانات والتصوير. تستفيد Cartopy من المكتبات القوية PROJ.4 والمكتوبة والمحددة الشكل ولها واجهة رسم بسيطة وبديهية لحزمة matplotlib لإنشاء خرائط جاهزة النشر. محاكي Brian Brian هو محاكي مجاني ومفتوح المصدر للشبكات العصبية، مكتوب بلغة البرمجة بايثون وهو متاح على جميع الأنظمة الأساسية تقريبًا. يساعد المحاكي في توفير وقت المعالجات إضافةً إلى توفير وقت الباحثين. صُمم محاكي Brian ليكون سهل التعلم والاستخدام ومرن جدًا وقابلًا للتوسيع بسهولة. مكتبة Glumpy Glumpy هي مكتبة تصور تفاعلي مبنية على OpenGL في بايثون، وهدفها هو تسهيل إنشاء تمثيلات مرئية سريعة وقابلة للتطوير وجميلة وتفاعلية وديناميكية. الخلاصة تُعد Numpy مكتبةً متعددة الاستخدامات، ولكن هذا لا يعني أنه ينبغي عليك استخدامها في كل الحالات. عرضنا في هذا الدليل بعض البدائل (بما في ذلك تعليمات بايثون العادية) وهي تستحق البحث والتعب، ولكن الاختيار عائدٌ لك، إذ عليك أن تفكر بالحل الأفضل بالنسبة لك من حيث وقت التطوير ووقت الحساب والجهد الذي ستحتاجه للتعديل. في حال قررت تصميم حلك الخاص فيجب عليك اختباره والمحافظة عليه، ولكن في المقابل ستكون حرًا في تصميمه بالطريقة التي تريدها. من ناحية أخرى: إذا قررت الاعتماد على حزمة طرف ثالث، فستوفر الكثير من وقت التطوير ولكن ستضطر إلى تكييف الحزمة مع احتياجاتك الخاصة. ترجمة -وبتصرّف- للفصل Beyond NumPy من كتاب From Python to Numpy لصاحبه Nicolas P. Rougier. اقرأ أيضًا المقال السابق: تخصيص أسلوب المتجهات Custom vectorization عبر استعمال مكتبة NumPy مفاهيم متقدمة حول مكتبة NumPy في بايثون أهم 8 مكتبات بلغة البايثون تستخدم في المشاريع الصغيرة

تتمثل إحدى نقاط قوة مكتبة NumPy في إمكانية استخدامها لبناء كائنات جديدة أو لتصنيف كائن ndarray إلى صنف فرعي. هذه العملية السابقة مملة بعض الشيء ولكنها تستحق الجهد لأنها تسمح لك بتحسين كائن ndarray ليناسب مشكلتك. سندرس في هذا المقال حالتين من العالم الحقيقي (قائمة مقيدة بنوع typed list ومصفوفة مراعية للذاكرة memory-aware array) اللتان تُستخدمان على نطاق واسع في مشروع glumpy، بينما تُعد الحالة الثانية -وهي مصفوفة مضاعفة الدقة double precision array- شائعة في الدراسة الأكاديمية. القائمة المقيدة بنوع Typed list تسمى أيضًا المصفوفة غير المنتظمة، وهي قائمة بالعناصر التي تحتوي جميعها على نفس نوع البيانات (بمعنى NumPy). وتوفر كلًا من القائمة وواجهة برمجة تطبيقات ndarray. ولكن نظرًا لأن واجهات برمجة التطبيقات APIs المخصصة قد لا تكون متوافقة في بعض الحالات، علينا تحديد الخيارات. على سبيل المثال فيما يتعلق بالعامل "+" سنختار استخدام NumPy API، إذ تُضاف القيمة إلى كل عنصر على حدة بدلًا من توسيع القائمة عن طريق إلحاق عنصر جديد بها. >>> l = TypedList([[1,2], [3]]) >>> print(l) [1, 2], [3] >>> print(l+1) [2, 3], [4] نريد أن نوفر للكائن الجديد إمكانية إدراج العناصر وإنشائها وإزالتها بسلاسة من قائمة API، وهذا ما سنوضحه في الفقرات التالية. الإنشاء Creation نظرًا لأن الكائن ديناميكي بحكم التعريف، فمن المهم تقديم تابع إنشاء للأغراض العامة general-purpose قوي بما يكفي لتجنب الاضطرار إلى معالجات لاحقة، إذ تكلف معالجات مثل الإدراج، أو الحذف كثيرًا من العمليات ونحن نسعى إلى تجنبها. إليك اقتراح لإنشاء كائن TypedList: def __init__(self, data=None, sizes=None, dtype=float) """ Parameters ---------- data : array_like An array, any object exposing the array interface, an object whose __array__ method returns an array, or any (nested) sequence. sizes: int or 1-D array If `itemsize is an integer, N, the array will be divided into elements of size N. If such partition is not possible, an error is raised. If `itemsize` is 1-D array, the array will be divided into elements whose successive sizes will be picked from itemsize. If the sum of itemsize values is different from array size, an error is raised. dtype: np.dtype Any object that can be interpreted as a NumPy data type. """ تسمح واجهة برمجة التطبيقات API هذه بإنشاء قائمة فارغة أو إنشاء قائمة من بعض البيانات الخارجية. لاحظ أنه في الحالة الأخيرة نحتاج إلى تحديد كيفية تقسيم البيانات إلى عدة عناصر أو تقسيمها إلى عناصر بحجم 1، ويمكن أن يكون قسمًا عاديًا (أي أن كل عنصر بطول بيانات 2)، أو قسمًا مخصصًا (أي يجب تقسيم البيانات إلى عناصر بحجم 1 و2 و3 و4 عنصر). >>> L = TypedList([[0], [1,2], [3,4,5], [6,7,8,9]]) >>> print(L) [ [0] [1 2] [3 4 5] [6 7 8] ] >>> L = TypedList(np.arange(10), [1,2,3,4]) [ [0] [1 2] [3 4 5] [6 7 8] ] في هذه المرحلة السؤال هو ما إذا كان يجب إنشاء صنف فرعي مشتق من الصنف ndarray، أو استخدام مصفوفة ndarray داخلية لتخزين بياناتنا. في هذه الحالة، ليس من المنطقي استخدام صنف فرعي ndarray لأننا لا نريد تقديم واجهة ndarray. بدلًا من ذلك، سنستخدم مصفوفة ndarray لتخزين بيانات القائمة وسيوفر لنا اختيار هذا التصميم مزيدًا من المرونة. سنستخدم لتخزين حد limit لكل عنصر مصفوفة عناصر items، ستهتم بتخزين الموضع (البداية والنهاية) لكل عنصر. توجد حالتان منفصلتان لإنشاء قائمة، هما: لم تُقدم أية بيانات، أو قُدمت بعض البيانات. الحالة الأولى سهلة وتتطلب فقط إنشاء مصفوفتي data_ و items_. لاحظ أن حجم المصفوفتين غير فارغ لأن تغيير حجم المصفوفة سيكون مكلفًا جدًا في كل مرة ندخل عنصر جديد فيها. لذلك من الأفضل حجز بعض المساحة. الحالة الأولى: لم تُقدّم أية بيانات، فقط dtype. self._data = np.zeros(512, dtype=dtype) self._items = np.zeros((64,2), dtype=int) self._size = 0 self._count = 0 الحالة الثانية: قُدمت بعض البيانات بالإضافة إلى قائمة بأحجام العناصر. self._data = np.array(data, copy=False) self._size = data.size self._count = len(sizes) indices = sizes.cumsum() self._items = np.zeros((len(sizes),2),int) self._items[1:,0] += indices[:-1] self._items[0:,1] += indices الوصول Access بمجرد الانتهاء من عملية الإنشاء، لا بد من إجراء قليل من الحسابات في كل التوابع وتغيير مفتاح key مختلف عند الحصول على عنصر أو إدراجه أو ضبطه. إليك تعليمات تابع __getitem__، وهو تابع سهل ولكن هناك خطوة سلبية محتملة فيه: def __getitem__(self, key): if type(key) is int: if key < 0: key += len(self) if key < 0 or key >= len(self): raise IndexError("Tuple index out of range") dstart = self._items[key][0] dstop = self._items[key][1] return self._data[dstart:dstop] elif type(key) is slice: istart, istop, step = key.indices(len(self)) if istart > istop: istart,istop = istop,istart dstart = self._items[istart][0] if istart == istop: dstop = dstart else: dstop = self._items[istop-1][1] return self._data[dstart:dstop] elif isinstance(key,str): return self._data[key][:self._size] elif key is Ellipsis: return self.data else: raise TypeError("List indices must be integers") تمرين يُعد تعديل القائمة أكثر تعقيدًا، لأنه يتطلب إدارة الذاكرة بطريقة صحيحة، وبما أن ذلك لا يمثل صعوبة حقيقة، فقد تركنا هذا تمرينًا لك. في حال احتجت إلى مساعدة، يمكنك إلقاء نظرة على التعليمات البرمجية أدناه. كن حذرًا مع الخطوات السلبية ونطاق المفتاح وتوسع المصفوفة. عندما تحتاج المصفوفة الأساسية إلى التوسيع، فمن الأفضل توسيعها أكثر من اللازم لتجنب الحاجة إلى التوسع المستقبلي. ضبط عنصر Setitem يمكننا ضبط عنصر كما يلي: L = TypedList([[0,0], [1,1], [0,0]]) L[1] = 1,1,1 حذف عنصر Delitem يجري حذف عنصر على النحو التالي: L = TypedList([[0,0], [1,1], [0,0]]) del L[1] إدراج عنصر insert يمكننا إضافة عنصر كما يلي: L = TypedList([[0,0], [1,1], [0,0]]) L.insert(1, [3,3]) مصفوفة مراعية للذاكرة Memory aware array مكتبة Glumpy Glumpy هي مكتبة تمثيل مرئي تفاعلي مبنية على OpenGL في بايثون هدفها تسهيل إنشاء عروض مرئية سريعة وقابلة للتطوير وجميلة وتفاعلية وديناميكية. يوضح الشكل التالي محاكاة مجرة حلزونية باستخدام نظرية موجة الكثافة density wave theory: محاكاة مجرة حلزونية باستخدام نظرية موجة الكثافة. عرض لنمر باستخدام التجميعات واستدعائين GL. تعتمد Glumpy على تكامل محكم وسلس مع المصفوفات المعقدة، وهذا يعني أنه يمكنك معالجة بيانات وحدة معالجة الرسومات GPU كما تفعل مع مصفوفات numpy العادية، وستتولى مكتبة Glumpy بقية المهام، وأسرع طريقة لفهم ذلك هي بالمثال التالي: from glumpy import gloo dtype = [("position", np.float32, 2), # x,y ("color", np.float32, 3)] # r,g,b V = np.zeros((3,3),dtype).view(gloo.VertexBuffer) V["position"][0,0] = 0.0, 0.0 V["position"][1,1] = 0.0, 0.0 V هو VertexBuffer وهو يمثّل GPUData ومصفوفة numpy. عندما يُعدل V، تعتني Glumpy بحساب أصغر كتلة متجاورة من ذاكرة dirty حتى تحميلها على ذاكرة GPU. عندما يُستخدم هذا المخزن المؤقت على وحدة معالجة الرسومات GPU، تتولى Glumpy تحميل منطقة dirty في اللحظات الأخيرة، وهذا يعني أنه في حال عدم استخدامك V فلن يُحمّل أي شيء على وحدة معالجة الرسومات. في الحالة المذكورة أعلاه، تتكون آخر منطقة ذاكرة dirty محسوبة من 88 بايت بدءًا من الإزاحة 0 كما هو موضح أدناه: وبالتالي سينتهي الأمر بتحميل Glumpy 88 بايت بينما عُدّل 16 بايت فقط فعليًا. قد تتساءل عما إذا كان هذا هو الأمثل. في الواقع يكون الأمر كذلك في معظم الأوقات، لأن تحميل بعض البيانات إلى مخزن مؤقت يتطلب الكثير من العمليات على جانب GL ولكل استدعاء تكلفة ثابتة. صنف فرعي من Array كما هو موضح في توثيق Subclassing ndarray، فإن التصنيف الفرعي ndarray معقد بسبب حقيقة أن الحالات الجديدة لأصناف ndarray يمكن أن تظهر بثلاث طرق مختلفة: استدعاء صريح للباني. تحويل العرض view casting. قالب جديد. ومع ذلك فإن حالتنا أبسط لأننا مهتمون فقط بتحويل العرض، لذلك نحتاج فقط إلى تعريف التابع __new__ الذي سيُستدعى عند إنشاء كل نسخة. لذلك ستجهز صنف GPUData بخاصيتين: النطاقات extents: يمثل هذا النطاق الكامل للعرض نسبيًا بالنسبة للمصفوفة الأساسية، ويُخزّن مثل إزاحة بايت وحجم بايت. البيانات المنتظرة pending_data: يمثل هذا تجاور لمنطقة dirty ويُخزّن مثل (إزاحة البايت، حجم البايت) نسبيًا بالنسبة إلى خاصية extents. class GPUData(np.ndarray): def __new__(cls, *args, **kwargs): return np.ndarray.__new__(cls, *args, **kwargs) def __init__(self, *args, **kwargs): pass def __array_finalize__(self, obj): if not isinstance(obj, GPUData): self._extents = 0, self.size*self.itemsize self.__class__.__init__(self) self._pending_data = self._extents else: self._extents = obj._extents نطاقات الحوسبة نحتاج في كل مرة يُطلب فيها عرض جزئي للمصفوفة إلى حساب نطاقات هذا العرض الجزئي حتى نتمكن من الوصول إلى المصفوفة الأساسية. def __getitem__(self, key): Z = np.ndarray.__getitem__(self, key) if not hasattr(Z,'shape') or Z.shape == (): return Z Z._extents = self._compute_extents(Z) return Z def _compute_extents(self, Z): if self.base is not None: base = self.base.__array_interface__['data'][0] view = Z.__array_interface__['data'][0] offset = view - base shape = np.array(Z.shape) - 1 strides = np.array(Z.strides) size = (shape*strides).sum() + Z.itemsize return offset, offset+size else: return 0, self.size*self.itemsize تعقب البيانات المنتظرة pending data باستمرار تتمثل إحدى الصعوبات الإضافية في أننا لا نريد أن تتعقب جميع طرق العرض منطقة dirty المُستخدمة ولكن نريد تعقب فقط المصفوفة الأساسية. هذا هو سبب عدم إنشاء نسخة منself._pending_data في الحالة الثانية من تابع __array_finalize__. سنتعامل مع هذا عندما نحتاج إلى تحديث بعض البيانات أثناء استدعاء __setitem__ على سبيل المثال: def __setitem__(self, key, value): Z = np.ndarray.__getitem__(self, key) if Z.shape == (): key = np.mod(np.array(key)+self.shape, self.shape) offset = self._extents[0]+(key * self.strides).sum() size = Z.itemsize self._add_pending_data(offset, offset+size) key = tuple(key) else: Z._extents = self._compute_extents(Z) self._add_pending_data(Z._extents[0], Z._extents[1]) np.ndarray.__setitem__(self, key, value) def _add_pending_data(self, start, stop): base = self.base if isinstance(base, GPUData): base._add_pending_data(start, stop) else: if self._pending_data is None: self._pending_data = start, stop else: start = min(self._pending_data[0], start) stop = max(self._pending_data[1], stop) self._pending_data = start, stop الخلاصة كما هو موضح على موقع ويب NumPy، فإن مكتبة NumPy هي الحزمة الأساسية للحوسبة العلمية باستخدام بايثون، ولكن وكما هو موضح في هذا المقال يتجاوز استخدام نقاط قوة NumPy مجرد حاوية متعددة الأبعاد multi-dimensional container للبيانات العامة. يمكن استخدام ndarray مثل خاصية خاصة في حالة واحدة (أسلوب TypedList)، أو تصنيفًا فرعيًا مباشرًا لصنف ndarray (أسلوب GPUData) لتتبع الذاكرة في حالات أخرى، وقد وضحنا كيف يمكن أن تتوسع قدرات NumPy لتناسب احتياجات محددة للغاية. ترجمة -وبتصرّف- للفصل Custom vectorization من كتاب From Python to Numpy لصاحبه Nicolas P. Rougier. اقرأ أيضًا المقال السابق: الاعتماد على المتجهات في حل المشاكل باستعمال مكتبة NumPy في بايثون استخدام المتجهات vectorization في لغة بايثون مع مكتبة NumPy أهم 8 مكتبات بلغة البايثون تستخدم في المشاريع الصغيرة

يعد الاعتماد على المتجهات في المشاكل أصعب بكثير من متجهات الشيفرات لأنه يعني أنه يتعين عليك إعادة التفكير في مشكلتك من أجل جعلها معتمدة على المتجهات vectorizable، وهذا يعني غالبًا أن عليك استخدام خوارزمية مختلفة لحل مشكلتك، أو حتى الأسوأ من ذلك، فقد تحتاج إلى ابتكار خوارزمية جديدة. وبالتالي تكمن الصعوبة في التفكير خارج الصندوق. لتوضيح ذلك دعنا نفكر في مشكلة بسيطة، ولنفرض لدينا متجهين X و Y، ونريد حساب مجموع [X * Y [j لجميع أزواج المؤشرات i و j. أحد الحلول البسيطة والواضحة هو كتابة: def compute_python(X, Y): result = 0 for i in range(len(X)): for j in range(len(Y)): result += X[i] * Y[j] return result ولكن تتضمن الشيفرات السابقة حلقتين وبالتالي فإن التنفيذ سيكون بطيئًا: >>> X = np.arange(1000) >>> timeit("compute_python(X,X)") 1 loops, best of 3: 0.274481 sec per loop كيف يمكن جعل المشكلة معتمدة على المتجهات إذًا؟ إذا كنت تتذكر مقرر الجبر الخطي، فربما تكون قد حددت التعبير [X * Y [j ليكون مشابهًا جدًا لتعبير المصفوفة، لذلك ربما يمكننا الاستفادة من سرعة NumPy. أحد الحلول الخاطئة هو بكتابة: def compute_numpy_wrong(X, Y): return (X*Y).sum() هذا الحل خاطئ لأن التعبير X * Y يحسب في الواقع متجهًا جديدًا Z وفقًا لما يلي: [Z = X * Y [i، وليس هذا ما نريده، لذلك سنجرب استغلال ميزة البث broadcasting لمكتبة NumPy وإعادة تشكيل المتجهين أولًا ثم ضربهما: def compute_numpy(X, Y): Z = X.reshape(len(X),1) * Y.reshape(1,len(Y)) return Z.sum() أصبح لدينا [Z[i,j] == X[i,0]*Y[0,j، وإذا حسبنا مجموع كل عناصر Z، سنحصل على النتيجة المتوقعة. سنتأكد الآن من مقدار التحسن في السرعة الذي اكتسبناه في هذه الطريقة: >>> X = np.arange(1000) >>> timeit("compute_numpy(X,X)") 10 loops, best of 3: 0.00157926 sec per loop هذا أفضل، فقد حصلنا على عامل سرعة أكبر بحوالي 150 مرة، ولكن لا يزال بإمكاننا تحسين هذا العامل. بالعودة إلى الطريقة الأولى، يمكنك أن تلاحظ أن الحلقة الداخلية تستخدم [X [i الذي لا يعتمد على المؤشر j، مما يعني أنه يمكن إزالته من الحلقة الداخلية. يمكن إعادة كتابة الشيفرة على النحو التالي: def compute_numpy_better_1(X, Y): result = 0 for i in range(len(X)): Ysum = 0 for j in range(len(Y)): Ysum += Y[j] result += X[i]*Ysum return result ونظرًا لأن الحلقة الداخلية لا تعتمد على المؤشر i، فيمكننا أن نحسبها مرةً واحدةً فقط: def compute_numpy_better_2(X, Y): result = 0 Ysum = 0 for j in range(len(Y)): Ysum += Y[j] for i in range(len(X)): result += X[i]*Ysum return result خفضنا التعقيد بفضل إزالة الحلقة الداخلية من "(O(n2" إلى "(O(n". وبنفس الطريقة يمكننا كتابة: def compute_numpy_better_3(x, y): Ysum = 0 for j in range(len(Y)): Ysum += Y[j] Xsum = 0 for i in range(len(X)): Xsum += X[i] return Xsum*Ysum أخيرًا، بعد أن أدركنا أننا نحتاج فقط إلى حاصل ضرب المجموع لكلٍ من X و Y على التوالي، يمكننا الاستفادة من دالة np.sum وكتابة: def compute_numpy_better(x, y): return np.sum(y) * np.sum(x) وبالتالي حصلنا على شيفرات برمجية أقصر وأرتب وأسرع: >>> X = np.arange(1000) >>> timeit("compute_numpy_better(X,X)") 1000 loops, best of 3: 3.97208e-06 sec per loop بهذه الطريقة نكون قد أعدنا صياغة المشكلة، مستغلين حقيقة أن: "‎∑ij XiYj = ∑iXi∑jYj". وتعلمنا الآن أن هناك نوعين من التحويل المعتمد على المتجهات، وهما: متجهات الشيفرات ومتجهات المشكلة. يُعد النوع الثاني من المتجهات هو الأصعب والأكثر أهمية لأنك من خلاله يمكنك توقع مكاسب هائلة في تحسين السرعة. اكتسبنا في هذا المثال البسيط عامل سرعة أكبر بمقدار 150 مرة عند استخدام متجهات الشيفرة، واكتسبنا عامل قدره 70000 باستخدام متجهات المشكلة، أي من خلال كتابة مشكلتنا بطريقة مختلفة فقط. إذا أعدنا كتابة الحل الأخير بطريقة بايثون فإن التحسين سيكون جيدًا ولكن ليس مثل استخدام numpy: def compute_python_better(x, y): return sum(x)*sum(y) هذا الإصدار الجديد من بايثون أسرع بكثير من إصدار بايثون السابق، لكنه لا يزال أبطأ 50 مرةً من الإصدار numpy: >>> X = np.arange(1000) >>> timeit("compute_python_better(X,X)") 1000 loops, best of 3: 0.000155677 sec per loop إيجاد المسار يتعلق اكتشاف المسار بإيجاد أقصر مسار في الرسم البياني. يمكن تقسيم هذا إلى مشكلتين مختلفتين: العثور على مسار بين كل عقدتين في الرسم البياني والعثور على أقصر مسار، وسنوضح ذلك من خلال اكتشاف المسار في متاهة، لذلك ستكون المهمة الأولى هي بناء المتاهة. متاهة سياج أخضر تحيط بمنزل Longleat الفخم في إنجلترا. صورة Prince Rurik، 2005. بناء متاهة توجد عدة خوارزميات لإنشاء متاهة، لكننا سنستخدم هنا خوارزمية معينة تعمل عن طريق إنشاء "n" (الكثافة) جزيرة بطول "p" (التعقيد). تنشأ الجزيرة من خلال اختيار نقطة بداية عشوائية بإحداثيات فردية، ثم نختار اتجاهًا عشوائيًا. إذا كانت الخلية التي تملك خطوتين في اتجاه معين حرة سيُضاف جدار في أول وثاني خطوة في هذا الاتجاه، وتُكرر العملية n خطوة لهذه الجزيرة، وبالتالي ننشيء p جزيرة. يكون كل من n و p أعداد عشرية float لتتكيف مع حجم المتاهة. تكون الجزر صغيرة جدًا مع هذا التعقيد المنخفض ويكون حل المتاهة سهلًا، وتحتوي المتاهة على مزيدٍ من "الغرف الفارغة الكبيرة" مع الكثافة المنخفضة. def build_maze(shape=(65, 65), complexity=0.75, density=0.50): # Only odd shapes shape = ((shape[0]//2)*2+1, (shape[1]//2)*2+1) # ضبط التعقيد والكثافة نسبةً إلى حجم المتاهة n_complexity = int(complexity*(shape[0]+shape[1])) n_density = int(density*(shape[0]*shape[1])) # بناء المتاهة الفعلية Z = np.zeros(shape, dtype=bool) # ملء الحدود Z[0, :] = Z[-1, :] = Z[:, 0] = Z[:, -1] = 1 # تبدأ الجزر بنقطة مع انحياز من طرف الحدود P = np.random.normal(0, 0.5, (n_density, 2)) P = 0.5 - np.maximum(-0.5, np.minimum(P, +0.5)) P = (P*[shape[1], shape[0]]).astype(int) P = 2*(P//2) # إنشاء الجزر for i in range(n_density): # Test for early stop: if all starting point are busy, this means we # won't be able to connect any island, so we stop. T = Z[2:-2:2, 2:-2:2] if T.sum() == T.size: break x, y = P[i] Z[y, x] = 1 for j in range(n_complexity): neighbours = [] if x > 1: neighbours.append([(y, x-1), (y, x-2)]) if x < shape[1]-2: neighbours.append([(y, x+1), (y, x+2)]) if y > 1: neighbours.append([(y-1, x), (y-2, x)]) if y < shape[0]-2: neighbours.append([(y+1, x), (y+2, x)]) if len(neighbours): choice = np.random.randint(len(neighbours)) next_1, next_2 = neighbours[choice] if Z[next_2] == 0: Z[next_1] = 1 Z[next_2] = 1 y, x = next_2 else: break return Z هنا رسم متحرك يوضح عملية التوليد. بناء متاهة متدرجة مع التحكم في التعقيد والكثافة. خوارزمية الاتساع أولا Breadth-first تعالج خوارزمية البحث ذات الاتساع أولًا Breadth-first (وكذلك خوارزمية العمق أولًا depth-first) مشكلة إيجاد مسار بين عقدتين من خلال فحص جميع الاحتمالات بدءًا من عقدة الجذر والتوقف بمجرد العثور على حل (الوصول إلى العقدة الوجهة). تعمل هذه الخوارزمية في الوقت الخطي بالتعقيد (|O(|V| + |E، إذ يمثل V عدد الرؤوس vertices و E هو عدد الحواف. كتابة مثل هذه الخوارزمية ليس صعبًا بشرط أن يكون لديك بنية البيانات الصحيحة. في حالتنا فإن تمثيل المتاهة على شكل مصفوفة ليس هو الحل الأنسب، ونحتاج إلى تحويله إلى رسم بياني حقيقي وفق اقتراح فالنتين بريوخانوف Valentin Bryukhanov. def build_graph(maze): height, width = maze.shape graph = {(i, j): [] for j in range(width) for i in range(height) if not maze[i][j]} for row, col in graph.keys(): if row < height - 1 and not maze[row + 1][col]: graph[(row, col)].append(("S", (row + 1, col))) graph[(row + 1, col)].append(("N", (row, col))) if col < width - 1 and not maze[row][col + 1]: graph[(row, col)].append(("E", (row, col + 1))) graph[(row, col + 1)].append(("W", (row, col))) return graph ملاحظة: إذا استخدمنا خوارزمية العمق أولاً فلا يوجد ضمان للعثور على أقصر مسار، ويمكننا فقط للعثور على مسار إذا كان موجودًا. بعد الانتهاء من ذلك، تكون كتابة خوارزمية الاتساع أولًا أمرًا سهلًا. نبدأ من عقدة البداية ونزور العقد في العمق الحالي ونكرر العملية حتى الوصول إلى العقدة النهائية إن أمكن. والسؤال هنا: هل نحصل بهذه الطريقة على أقصر طريق لاستكشاف الرسم البياني؟ في هذه الحالة المحددة "نعم"، لأننا لا نملك رسمًا بيانيًا مرجحًا للحافة، أي أن جميع الحواف لها نفس الوزن (أو التكلفة). def breadth_first(maze, start, goal): queue = deque([([start], start)]) visited = set() graph = build_graph(maze) while queue: path, current = queue.popleft() if current == goal: return np.array(path) if current in visited: continue visited.add(current) for direction, neighbour in graph[current]: p = list(path) p.append(neighbour) queue.append((p, neighbour)) return None طريقة بيلمان فورد Bellman-Ford خوارزمية بيلمان فورد Bellman-Ford قادرة على إيجاد المسار الأمثل في الرسم البياني باستخدام عملية الانتشار، إذ تعثر على المسار الأفضل عن طريق تصاعد التدرج الناتج. تعمل هذه الخوارزمية في الوقت التربيعي (|O(|V||E، إذ يمثّل V عدد الرؤوس و E عدد الحواف. في حالتنا البسيطة لن نواجه الكثير من المشاكل. الخوارزمية موضحة أدناه في الشكل (القراءة من اليسار إلى اليمين ومن أعلى إلى أسفل). بعد فهم الخوارزمية يمكننا صعود التدرج اللوني من عقدة البداية، ويمكنك التحقق من الشكل الذي يؤدي إلى أقصر طريق. قيم خوارزمية التكرار في متاهة بسيطة. بمجرد الوصول إلى المدخل، يكون من السهل العثور على أقصر طريق عن طريق صعود تدرج القيمة. نبدأ بتعيين عقدة الخروج على القيمة 1، بينما تُعيّن كل عقدة أخرى بالقيمة 0 باستثناء الجدران، ثم نكرر العملية بحيث تحسب القيمة الجديدة لكل خلية على أنها القيمة القصوى بين قيمة الخلية الحالية وقيم الخلايا الأربعة المُخفّضة المتجاورة (gamma=0.9 في الحالة أدناه). تبدأ العملية بمجرد أن تصبح قيمة عقدة البداية موجبة تمامًا. يُعد التنفيذ باستخدام NumPy أبسط إذا استفدنا من ميزات عامل الترشيح العام generic_filter (من scipy.ndimage) لعملية الانتشار: def diffuse(Z): # North, West, Center, East, South return max(gamma*Z[0], gamma*Z[1], Z[2], gamma*Z[3], gamma*Z[4]) # Build gradient array G = np.zeros(Z.shape) # Initialize gradient at the entrance with value 1 G[start] = 1 # Discount factor gamma = 0.99 # We iterate until value at exit is > 0. This requires the maze # to have a solution or it will be stuck in the loop. while G[goal] == 0.0: G = Z * generic_filter(G, diffuse, footprint=[[0, 1, 0], [1, 1, 1], [0, 1, 0]]) لكن في هذه الحالة بالتحديد يكون الأمر بطيئًا إلى حد ما، ومن الأفضل إعداد حلنا الخاص، وإعادة استخدام جزء من تشفير لعبة الحياة: # Build gradient array G = np.zeros(Z.shape) # Initialize gradient at the entrance with value 1 G[start] = 1 # Discount factor gamma = 0.99 # We iterate until value at exit is > 0. This requires the maze # to have a solution or it will be stuck in the loop. G_gamma = np.empty_like(G) while G[goal] == 0.0: np.multiply(G, gamma, out=G_gamma) N = G_gamma[0:-2,1:-1] W = G_gamma[1:-1,0:-2] C = G[1:-1,1:-1] E = G_gamma[1:-1,2:] S = G_gamma[2:,1:-1] G[1:-1,1:-1] = Z[1:-1,1:-1]*np.maximum(N,np.maximum(W, np.maximum(C,np.maximum(E,S)))) بمجرد الانتهاء من ذلك، يمكننا صعود التدرج للعثور على أقصر طريق كما هو موضح في الشكل أدناه: إيجاد المسار باستخدام خوارزمية بيلمان-فورد. تشير الألوان المتدرجة إلى القيم المنتشرة من نقطة نهاية المتاهة (أسفل اليمين). عُثر على المسار من خلال التدرج التصاعدي من الهدف. خوارزمية تدفقات السوائل Fluid Dynamics فيما يلي سنوضح طريقة عمل خوارزمية تدفقات السوائل، وكيفية تنفيذها برمجيًا. التدفق الهيدروديناميكي بمستويين مختلفين من التقريب، نهر نيكار، هايدلبرغ، ألمانيا. صورة Steven Mathey، 2012. طريقة لاغرانج Lagrangian مقابل أويلريان Eulerian في نظرية الحقل الكلاسيكية، فإن مواصفات لاغرانج للحقل هي طريقة للنظر إلى حركة السوائل، إذ يتبّع المراقب حركة السائل عبر المكان والزمان، ويعطينا رسم موقع كل جزء عبر الزمن مسار ذلك الجزء، ويمكن تخيل هذا على أنه الجلوس في قارب والانجراف في النهر. تعد مواصفات أويلريان لحقل التدفق طريقةً للنظر إلى حركة السوائل التي تركز على مواقع محددة في الفضاء الذي يتدفق من خلاله السائل مع مرور الوقت. يمكن تصور ذلك من خلال الجلوس على ضفة النهر ومشاهدة الماء يمر من الموقع الثابت. أي بعبارة أخرى، في حالة أويلريان يُقسم جزءٌ من الفضاء إلى خلايا وتحتوي كل خلية على متجه السرعة ومعلومات أخرى، مثل الكثافة ودرجة الحرارة. نحتاج في حالة لاغرانج إلى فيزياء تعتمد على الجسيمات مع تفاعلات ديناميكية وعمومًا نحتاج إلى عدد كبير من الجسيمات. لكلتا الطريقتين مزايا وعيوب ويعتمد الاختيار بين الطريقتين على طبيعة مشكلتك، ويمكنك أيضًا مزج الطريقتين في طريقة هجينة. مع ذلك، فإن أكبر مشكلة في المحاكاة القائمة على الجسيمات هي أن تفاعل الجسيمات يتطلب إيجاد الجسيمات المجاورة وهذا له تكلفة كما رأينا في حالة أسراب الطيور Boids في المقال السابق. إذا استهدفنا بايثونو NumPy فقط، فمن الأفضل اختيار طريقة أويلريان لأن الاعتماد على المتجهات سيكون عديم الأهمية تقريبًا مقارنةً بطريقة لاغرانج. كتابة الشيفرة باستخدام مكتبة NumPy لن نشرح كل النظرية الكامنة وراء ديناميك السوائل الحسابية لوجود العديد من الموارد عبر الإنترنت التي تشرح هذا بطريقة جيدة. لماذا نختار حسابات السائل إذًا؟ لأن النتائج تكون غالبًا جيدة ورائعة. سنعمل على تبسيط المشكلة أكثر من خلال تنفيذ طريقة من رسوميات الحاسوب بحيث يكون الهدف هو السلوك المقنع. سنعمل على تبسيط المشكلة أكثر من خلال تنفيذ طريقة من رسوميات الحاسوب، إذ لا يكون الهدف هو التصحيح ولكن السلوك المقنع. كتب جو ستام Jos Stam مقالًا لطيفًا جدًا لصالح SIGGRAPH عام 1999 يصف تقنية الحصول على سوائل مستقرة بمرور الوقت (أي لا يتباعد حلها على المدى الطويل). وكتب ألبيرتو سانتينيAlberto Santini تعليمات باستخدام لغة البرمجة بايثون (استخدم حينها numarray). الآن، علينا فقط تكييفها مع مكتبة NumPy الحديثة وتسريعها قليلًا باستخدام الحيل المعقدة. لن نشرح الشيفرات لأنها طويلة جدًا. محاكاة الدخان باستخدام خوارزمية السوائل المستقرة بواسطة جو ستام. تأتي معظم مقاطع الفيديو الصحيحة من حزمة Glumpy وتستخدم GPU (عمليات الإطارات الاحتياطية، أي عدم وجود OpenCL أو CUDA) لإجراء عمليات حسابية أسرع. أخذ عينات الضجيج الأزرق يشير الضجيج الأزرق blue noise إلى مجموعات العينات التي لها توزيعات عشوائية وموحدة مع عدم وجود أي انحياز طيفي. هذا الضجيج مفيد جدًا في مجموعة متنوعة من تطبيقات الرسوم، مثل التصيير rendering والتردد dithering وما إلى ذلك، وقد اقترحت عدة طرق لتحقيق مثل هذا الضجيج، ولكن الطريقة الأكثر بساطة هي بالتأكيد طريقة دارت DART. طريقة دارت DART تُعد طريقة دارت DART واحدةً من أقدم وأبسط الطرق، وهي تعمل عن طريق الرسم المتسلسل لنقاط عشوائية موحدة وقبول فقط تلك التي تقع على مسافة لا تقل عن كل عينة مقبولة سابقة، ولذلك فهي طريقة بطيئة جدًا لأن كل مرشح جديد يحتاج إلى اختباره مقابل المرشحين المقبولين السابقين؛ وكلما زادت النقاط التي تقبلها كانت الطريقة أبطأ. لنفكر في وحدة السطح والحد الأدنى لنصف القطر r لفرضه في كل نقطة. التعبئة الأكثر كثافة للدوائر في المستوى هي الشبكة السداسية لشكل خلية النحل، وتحسب الكثافة بالصيغة: بالنظر إلى الدوائر ذات نصف القطر r، يمكننا أن نكتب: نحن نعلم الحد الأعلى النظري لعدد الأقراص الممكن وضعها على السطح ولكن من المحتمل ألا نصل إلى هذا الحد الأقصى بسبب المواضع العشوائية. علاوةً على ذلك، نظرًا لرفض الكثير من النقاط بعد قبول القليل منها، سنحتاج إلى وضع حد لعدد المحاكمات الفاشلة المتتالية قبل أن نوقف العملية كاملةً. import math import random def DART_sampling(width=1.0, height=1.0, r = 0.025, k=100): def distance(p0, p1): dx, dy = p0[0]-p1[0], p0[1]-p1[1] return math.hypot(dx, dy) points = [] i = 0 last_success = 0 while True: x = random.uniform(0, width) y = random.uniform(0, height) accept = True for p in points: if distance(p, (x, y)) < r: accept = False break if accept is True: points.append((x, y)) if i-last_success > k: break last_success = i i += 1 return points تركنا لك تمرين إيجاد اتجاه توجيه لطريقة دارت DART. إذ تكمن الفكرة في إجراء حساب مسبق لعدد منتظم كافٍ من العينات العشوائية، إضافةً إلى المسافات المقترنة واختبار تضمينها التسلسلي. طريقة بريدسون Bridson إذا لم تسبب متجهات الطريقة السابقة أي صعوبة حقيقية، فإن تحسين السرعة ليس جيدًا وتظل الجودة منخفضة وتعتمد على المعامل k؛ فكلما كان هذا المعامل أعلى، كان ذلك أفضل لأنه يتحكم بصورةٍ أساسية في مدى صعوبة محاولة إدخال عينة جديدة، ولكن عندما يكون هناك عدد كبير فعلًا من العينات المقبولة، فإن فرصة واحد فقط ستسمح لنا بالعثور على موقع لإدخال عينة جديدة. يمكننا زيادة قيمة k ولكن هذا سيجعل الطريقة أبطأ دون أي ضمان في الجودة. حان الوقت للتفكير خارج الصندوق ولحسن الحظ، فعل روبرت بريدسون ذلك لنا واقترح طريقة بسيطة لكنها فعالة: الخطوة 0: تهيئة شبكة خلفية ذات n بُعد لتخزين العينات وتسريع البحث المكاني. نختار حجم الخلية المراد تقييدها من خلال بحيث تحتوي كل خلية شبكية على عينة واحدة على الأكثر، وبالتالي يمكن تنفيذ الشبكة مثل مصفوفة بسيطة من الأعداد الصحيحة ذات البعد n: يشير الرقم الافتراضي "1-" إلى عدم وجود عينة، بينما يعطي العدد الصحيح غير السالب مؤشر العينة الموجود في الخلية. الخطوة 1. حدد العينة الأولية x0، التي تُختار عشوائيًا من المجال، وأدخلها في شبكة الخلفية وهيّئ "القائمة النشطة" (مجموعة من مؤشرات العينة) باستخدام هذا المؤشر (صفر). الخطوة 2: عندما لا تكون القائمة النشطة فارغة، اختر منها مؤشرًا عشوائيًا وليكن i. عليك توليد ما يصل إلى "k" من النقاط المختارة بطريقة موحدة من الحلقة الكروية بين نصف القطر "r" و "2r" وحول xi. تحقق من أجل كل نقطة مما إذا كانت على مسافة "r" من العينات الموجودة (باستخدام شبكة الخلفية لاختبار العينات القريبة فقط). إذا كانت نقطة ما بعيدة بحد كافٍ عن العينات الموجودة فأرسلها على أنها عينة تالية وأضفها إلى القائمة النشطة. إذا لم يعثر على مثل هذه النقط بعد "k" محاولة، عليك إزالة i من القائمة النشطة. لا يطرح التنفيذ مشكلة حقيقية ويُترك بمثابة تمرين للقارئ. لاحظ أن هذه الطريقة ليست سريعة فحسب ولكنها تقدم أيضًا جودة أفضل (المزيد من العينات) من طريقة دارت DART ومع وجود عدد كبير من المحاولات. الخلاصة يعد المثال الأخير الذي كنا ندرسه هو مثال جيد، لأن متجهات المشكلة أكثر أهمية من متجهات الشيفرات، وقد كنا محظوظين في هذه الحالة المحددة بما يكفي لإنجاز العمل ولكن لن يكون الأمر كذلك دائمًا. نأمل أن تكون الآن مقتنعًا أنه من الجيد البحث عن حلول بديلة بمجرد العثور على حل. ستعمل دائمًا على تحسين السرعة من خلال متجهات الشيفرات والتعليمات البرمجية ولكن في هذه العملية فقد تفوتك تحسينات كبيرة. ترجمة -وبتصرّف- للفصل Problem vectorization من كتاب From Python to Numpy لصاحبه Nicolas P. Rougier. اقرأ أيضًا المقال السابق: استخدام المتجهات vectorization في لغة بايثون مع مكتبة NumPy مفاهيم متقدمة حول مكتبة NumPy في بايثون أهم 8 مكتبات بلغة البايثون تستخدم في المشاريع الصغيرة

يُقصد باستخدام المتجهات vectorization في الشيفرة أن المشكلة التي تحاول حلها قابلةٌ للتحليل بطبيعتها ولا تتطلب سوى استخدام بعض الحيل لجعلها أسرع، وهذا لا يعني أن الأمر سهل أو مباشر، ولكنه على الأقل لا يستلزم إعادة التفكير في كامل مشكلتك. يتطلب الأمر بعض الخبرة لمعرفة أين يمكن استخدام المتجهات في الشيفرة. سنوضح ذلك من خلال مثال بسيط، إذ نريد تلخيص قائمتين من الأعداد الصحيحة. سنستخدم الدالة البسيطة التالية: def add_python(Z1,Z2): return [z1+z2 for (z1,z2) in zip(Z1,Z2)] إعادة كتابة الشيفرة السابقة عبر المتجهات من مكتبة NumPy بسهولة كبيرة بالشكل التالي: def add_numpy(Z1,Z2): return np.add(Z1,Z2) وكما هو متوقع تظهر المقارنة أن الطريقة الثانية هي الأسرع: >>> Z1 = random.sample(range(1000), 100) >>> Z2 = random.sample(range(1000), 100) >>> timeit("add_python(Z1, Z2)", globals()) 1000 loops, best of 3: 68 usec per loop >>> timeit("add_numpy(Z1, Z2)", globals()) 10000 loops, best of 3: 1.14 usec per loop كما أنها تتكيف أيضًا مع شكل Z1 و Z2، وهذا هو سبب عدم كتابة Z1 + Z2 لأن هذه الطريقة لن تعمل إذا كان كل من Z1 و Z2 قائمتين. تُفسّر عملية الجمع الداخلية في الطريقة الأولى المعتمدة على تابع لغة بايثون بطريقة مختلفة اعتمادًا على طبيعة الكائنين؛ وإذا أخذنا في الحسبان قائمتين متداخلتين، نحصل على المخرجات التالية: >>> Z1 = [[1, 2], [3, 4]] >>> Z2 = [[5, 6], [7, 8]] >>> Z1 + Z2 [[1, 2], [3, 4], [5, 6], [7, 8]] >>> add_python(Z1, Z2) [[1, 2, 5, 6], [3, 4, 7, 8]] >>> add_numpy(Z1, Z2) [[ 6 8] [10 12]] تربط الطريقة الأولى القائمتين معًا، بينما تربط الطريقة الثانية القوائم الداخلية معًا، وتحسب الطريقة الأخيرة ما هو متوقع. يمكنك إعادة كتابة تعليمات بايثون بحيث تضم قوائم متداخلة مختلفة الحجم. المتجه الموحد Uniform vectorization يُعد المتجه الموحد أبسط شكل من أشكال المتجهات، إذ تشترك جميع العناصر في نفس العملية الحسابية في كل خطوة زمنية مع عدم وجود معالجة محددة لأي عنصر. لعبة الحياة هي إحدى الحالات النمطية التي اخترعها جون كونواي John Conway وهي إحدى أقدم الأمثلة على المشغلات الآلية الخلوية cellular automata، التي تمثّل مجموعةً من الخلايا التي ترتبط معًا بمفهوم الجيران وتكون متجّهاتها بسيطة. سنحدّد بدايةً قواعد اللعبة ثم سندرس كيفية جعلها معتمدة على المتجهات: يظهر نسيج صدفة الحلزون نمط مشغل آلي خلوي على غلافها الخارجي. الصورة لريتشارد لينج 2005. لعبة الحياة يمكنك الحصول على فكرة حول لعبة الحياة في موقع ويكيبيديا؛ وهي لعبة معتمدة على المشغلات الآلية الخلوية، وابتكرها عالم الرياضيات البريطاني جون هورتون كونواي John Horton Conway في عام 1970. اللعبة هي لعبة بدون لاعب، مما يعني أن تطورها يُحدّد من خلال حالتها الأولية، ولا تحتاج اللعبة إلى مدخلات من اللاعبين، ويتفاعل المراقب مع لعبة الحياة من خلال إنشاء تكوين أولي ومراقبة كيفية تطور اللعبة. عالم لعبة الحياة هي شبكة متعامدة ثنائية الأبعاد لا نهائية من الخلايا المربعة، لكل منها حالتين محتملتين (حياة أو موت)، وترتبط حالة كل خلية بجاراتها الثمانية، وهي الخلايا المتجاورة أفقيًا، أو رأسيًا، أو قطريًا، وتحدث في كل خطوة في الوقت المناسب التحولات التالية: تموت أي خلية حية لها أقل من اثنين من الجيران الأحياء (كما لو كانت بسبب الاحتياجات الناجمة عن نقص السكان). تموت أي خلية لها أكثر من ثلاثة جيران أحياء (كما لو كان السبب هو الاكتظاظ). تعيش أي خلية حية لها اثنان أو ثلاثة من الجيران الأحياء للجيل القادم. تصبح أي خلية ميتة مع ثلاثة جيران أحياء بالضبط خليةً حية. يشكل النمط الأولي "أصل seed" النظام. يُنشأ الجيل الأول من خلال القواعد السابقة في وقت واحد على كل خلية (تحدث الولادات والوفيات في وقت واحد)، أي كل جيل هو نتيجة خالصة للجيل السابق، ويستمر تطبيق القواعد بطريقة مستمرة لخلق أجيال أخرى. كتابة شيفرة بلغة بايثون يمكننا كتابة شيفرة لعبة Game of Life في لغة بايثون باستخدام قائمة list تمثّل اللوحة، إذ من المفترض أن تتطور الخلايا. تُجهز هذه اللوحة بحدود صفرية "0" مما يسمح بتسريع الأمور قليلًا عن طريق تجنب إجراء اختبارات محددة للحدود عند حساب الجيران. Z = [[0,0,0,0,0,0], [0,0,0,1,0,0], [0,1,0,1,0,0], [0,0,1,1,0,0], [0,0,0,0,0,0], [0,0,0,0,0,0]] يكون عدد الجيران واضحًا عند أخذ الحدود بالحسبان: def compute_neighbours(Z): shape = len(Z), len(Z[0]) N = [[0,]*(shape[0]) for i in range(shape[1])] for x in range(1,shape[0]-1): for y in range(1,shape[1]-1): N[x][y] = Z[x-1][y-1]+Z[x][y-1]+Z[x+1][y-1] \ + Z[x-1][y] +Z[x+1][y] \ + Z[x-1][y+1]+Z[x][y+1]+Z[x+1][y+1] return N نكرر الخطوة في الوقت المناسب، ولذلك نحسب عدد الجيران لكل خلية داخلية ثم نحدث كامل اللوحة وفقًا للقواعد الأربعة المذكورة أعلاه: def iterate(Z): N = compute_neighbours(Z) for x in range(1,shape[0]-1): for y in range(1,shape[1]-1): if Z[x][y] == 1 and (N[x][y] < 2 or N[x][y] > 3): Z[x][y] = 0 elif Z[x][y] == 0 and N[x][y] == 3: Z[x][y] = 1 return Z يوضح الشكل أدناه أربعة تكرارات على منطقة 4×4، إذ أن الحالة الأولية هي طائرة شراعية glider، وهي بنية اكتُشفت عام 1970. من المعروف أن نمط الطائرة الشراعية يكرر نفسه خطوةً واحدةً قطريًا خلال أربعة تكرارات. كتابة شيفرة Numpy بناءً على التعليمات السابقة المكتوبة بلغة بايثون، نلاحظ أن شيفرة لعبة الحياة ينقسم إلى قسمين، أحدهما مسؤول عن حساب عدد الجيران والآخر مسؤول عن تطبيق القواعد. يُعد حساب عدد الجيران أمرًا سهلًا نسبيًا، لأننا حرصنا على إضافة حد فارغ حول النطاق arena. ومن خلال النظر في العروض الجزئية للنطاق يمكننا الوصول إلى الجيران بطريقة حدسية كما توضح هذه الحالة أحادية البعد: يتطلب الانتقال إلى الحالة ثنائية البعد حسابات إضافية للتأكد من مراعاة جميع الجيران الثمانية. N = np.zeros(Z.shape, dtype=int) N[1:-1,1:-1] += (Z[ :-2, :-2] + Z[ :-2,1:-1] + Z[ :-2,2:] + Z[1:-1, :-2] + Z[1:-1,2:] + Z[2: , :-2] + Z[2: ,1:-1] + Z[2: ,2:]) يمكننا لتطبيق القاعدة كتابة التعليمات البرمجية باستخدام تابع argwhere الذي سينتج عنه رقم المؤشر عندما يكون الشرط المحدد صحيحًا: # المصفوفات المسطحة N_ = N.ravel() Z_ = Z.ravel() # تطبيق القواعد R1 = np.argwhere( (Z_==1) & (N_ < 2) ) R2 = np.argwhere( (Z_==1) & (N_ > 3) ) R3 = np.argwhere( (Z_==1) & ((N_==2) | (N_==3)) ) R4 = np.argwhere( (Z_==0) & (N_==3) ) # ضبط قيم جديدة Z_[R1] = 0 Z_[R2] = 0 Z_[R3] = Z_[R3] Z_[R4] = 1 # كن متأكدًا من بقاء الحدود فارغة Z[0,:] = Z[-1,:] = Z[:,0] = Z[:,-1] = 0 رغم عدم استخدام الإصدار الأول للحلقات المتداخلة لكنه يبقى بعيدًا على أن يكون أمثليًا بسبب استخدام استدعاءات argwhere التي قد تكون بطيئة جدًا، ويمكننا بدلًا من ذلك تحليل القواعد إلى الخلايا التي ستبقى على قيد الحياة (أي التي ستبقى عند قيمة 1) والخلايا التي ستلد وذلك بالاستفادة من القدرة المنطقية لمكتبة Numpy: birth = (N==3)[1:-1,1:-1] & (Z[1:-1,1:-1]==0) survive = ((N==2) | (N==3))[1:-1,1:-1] & (Z[1:-1,1:-1]==1) Z[...] = 0 Z[1:-1,1:-1][birth | survive] = 1 إذا دققت في أول سطرين ستلاحظ أن المتغيرين birth و survive مصفوفتان يمكن استخدامهما لتعيين قسم المصفوفة Z إلى القيمة "1" بعد محوها. تمرين يمكن أن ينتج عن تفاعل الأصناف الكيميائية وانتشارها مجموعةً متنوعةً من الأنماط، تذكرنا غالبًا بتلك التي نراها في الطبيعة، ومن أمثلة هذا التفاعل نموذج معادلات جراي سكوت Gray-Scott، ولمزيد من المعلومات حول هذا النظام الكيميائي راجع مقالة الأنماط المعقدة في نظام بسيط. لنفكر في نوعين كيميائيين U و V، ولكل منهما تركيز u و v ومعدل انتشاء Du و Dv. يُحوّل V إلى P بمعدل تحويل k. يمثل f معدل العملية التي تشبع U وتجعل كل من U و V و P يتلاشى. ويمكن كتابة ذلك على النحو التالي: استنادًا إلى مثال لعبة الحياة سنجرب تنفيذ نظام نشر التفاعل. فيما يلي مجموعة من المعاملات التي سنحتاج إلى اختبارها: table { width: 100%; } thead { vertical-align: middle; text-align: center; } td, th { border: 1px solid #dddddd; text-align: right; padding: 8px; text-align: inherit; } tr:nth-child(even) { background-color: #dddddd; } الاسم Du Dv f k Bacteria 1 0.16 0.08 0.035 0.065 Bacteria 2 0.14 0.06 0.035 0.065 Coral 0.16 0.08 0.060 0.062 Fingerprint 0.19 0.05 0.060 0.062 Spirals 0.10 0.10 0.018 0.050 Spirals Dense 0.12 0.08 0.020 0.050 Spirals Fast 0.10 0.16 0.020 0.050 Unstable 0.16 0.08 0.020 0.055 Worms 1 0.16 0.08 0.050 0.065 Worms 2 0.16 0.08 0.054 0.063 Zebrafish 0.16 0.08 0.035 0.060 والنتائج هي كالآتي: والآتي: إضافةً إلى ما يلي: المتجهات الزمنية Temporal vectorization مجموعة Mandelbrot هي مجموعةٌ من الأعداد العقدية c التي لا تتباعد diverge فيها الدالة "fc(z) = z2 + c" عند تكرارها بدءًا من z = 0، أي من أجل التسلسل "… ,((fc(0), fc(fc(0" وتظل محدودة بالقيمة المطلقة. من السهل حسابها ولكنها قد تستغرق وقتًا طويلَا جدًا لأنك تحتاج إلى التأكد من أن رقمًا معينًا لا يتباعد، ويُنفّذ ذلك من خلال تكرار الحساب حتى أقصى عدد من التكرارات، وبعد ذلك إذا كان الرقم لا يزال ضمن الحدود يكون غير متباعد، وبالتأكيد كلما زاد عدد مرات التكرار زادت الدقة التي تحصل عليها. بروكلي رومانسكو، يُظهر شكلًا قريبًا من كسورية fractal طبيعية. صورة جون سوليفان، 2004. كتابة شيفرة بايثون تكتب التعليمات البرمجية على النحو التالي: def mandelbrot_python(xmin, xmax, ymin, ymax, xn, yn, maxiter, horizon=2.0): def mandelbrot(z, maxiter): c = z for n in range(maxiter): if abs(z) > horizon: return n z = z*z + c return maxiter r1 = [xmin+i*(xmax-xmin)/xn for i in range(xn)] r2 = [ymin+i*(ymax-ymin)/yn for i in range(yn)] return [mandelbrot(complex(r, i),maxiter) for r in r1 for i in r2] الجزء المثير للاهتمام والبطيء في الشيفرات السابقة هو دالة mandelbrot التي تحسب التسلسل "(((…fc(fc(fc". جعل هذه الشيفرات معتمدة على المتجهات ليس أمرًا واضحًا تمامًا بسبب تعليمة return الداخلية التي توحي بعملية فاضلية للعنصر. بمجرد أن يتباعد لا نحتاج إلى التكرار مرةً أخرى ويمكننا إرجاع عدد التكرارات في حال عدم التباعد، ولكن كيف يمكن أن تفعل Numpy الشيء نفسه؟ كتابة الشيفرة باستخدام مكتبة Numpy الحيلة هي البحث في كل قيم التكرار التي لم تتباعد بعد وتحديث المعلومات المتعلقة بهذه القيم فقط. نظرًا لأننا نبدأ من "Z = 0"، فإننا نعلم أن كل قيمة ستُحدث مرةً واحدةً على الأقل (عندما تكون مساويةً للصفر فإنها لم تتباعد بعد)، وستتوقف عن التحديث بمجرد تباعدها. ولفعل ذلك سنستخدم فهرسة fancy indexing مع الدالة (less(x1,x2 التي ترجع القيمة الحقيقية للعنصر المحقق (x1<x2). def mandelbrot_numpy(xmin, xmax, ymin, ymax, xn, yn, maxiter, horizon=2.0): X = np.linspace(xmin, xmax, xn, dtype=np.float32) Y = np.linspace(ymin, ymax, yn, dtype=np.float32) C = X + Y[:,None]*1j N = np.zeros(C.shape, dtype=int) Z = np.zeros(C.shape, np.complex64) for n in range(maxiter): I = np.less(abs(Z), horizon) N[I] = n Z[I] = Z[I]**2 + C[I] N[N == maxiter-1] = 0 return Z, N وإليك النتيجة: >>> xmin, xmax, xn = -2.25, +0.75, int(3000/3) >>> ymin, ymax, yn = -1.25, +1.25, int(2500/3) >>> maxiter = 200 >>> timeit("mandelbrot_python(xmin, xmax, ymin, ymax, xn, yn, maxiter)", globals()) 1 loops, best of 3: 6.1 sec per loop >>> timeit("mandelbrot_numpy(xmin, xmax, ymin, ymax, xn, yn, maxiter)", globals()) 1 loops, best of 3: 1.15 sec per loop كتابة الشيفرة الأسرع باستخدام Numpy سيكون الربح بزيادة السرعة خمسة أضعاف تقريبًا وليس بالقدر الذي نتوقعه، ومن أسباب هذه المشكلة هو أن الدالة np.less تشير إلى اختبارات xn × yn في كل تكرار، بينما نعلم أن بعض القيم قد تباعدت فعلًا، وحتى إذا أجريت هذه الاختبارات على المستوى C (من خلال numpy)، ستبقى التكلفة كبيرة. توجد طريقةٌ أخرى اقترحها دان جودمان Dan Goodman وذلك بالعمل على مصفوفة ديناميكية في كل تكرار بحيث تخزن فقط النقاط التي لم تتباعد بعد، وهذا يتطلب إضافة أسطر برمجية إضافية ولكن ستكون النتيجة أسرع وتؤدي إلى تحسين عامل السرعة عشر مرات مقارنةً بالطريقة الأولى. def mandelbrot_numpy_2(xmin, xmax, ymin, ymax, xn, yn, itermax, horizon=2.0): Xi, Yi = np.mgrid[0:xn, 0:yn] Xi, Yi = Xi.astype(np.uint32), Yi.astype(np.uint32) X = np.linspace(xmin, xmax, xn, dtype=np.float32)[Xi] Y = np.linspace(ymin, ymax, yn, dtype=np.float32)[Yi] C = X + Y*1j N_ = np.zeros(C.shape, dtype=np.uint32) Z_ = np.zeros(C.shape, dtype=np.complex64) Xi.shape = Yi.shape = C.shape = xn*yn Z = np.zeros(C.shape, np.complex64) for i in range(itermax): if not len(Z): break # Compute for relevant points only np.multiply(Z, Z, Z) np.add(Z, C, Z) # Failed convergence I = abs(Z) > horizon N_[Xi[I], Yi[I]] = i+1 Z_[Xi[I], Yi[I]] = Z[I] # Keep going with those who have not diverged yet np.negative(I,I) Z = Z[I] Xi, Yi = Xi[I], Yi[I] C = C[I] return Z_.T, N_.T وتصبح النتيجة كما يلي: >>> timeit("mandelbrot_numpy_2(xmin, xmax, ymin, ymax, xn, yn, maxiter)", globals()) 1 loops, best of 3: 510 msec per loop إظهار النتائج من أجل تصور نتائجنا يمكننا عرض المصفوفة N مباشرةً باستخدام الأمر imshow من مكتبة matplotlib ولكن هذا سينتج صورة مخططة banded وهي نتيجة معروفة لخوارزمية عد الهروب escape count المستخدمة. يمكن التخلص من هذا الربط من خلال حساب الهروب الجزئي وذلك عن طريق قياس مدى هبوط النقطة المتكررة خارج حدود قطع الهروب. فيما يلي صورة للنتيجة، إذ استخدمنا التسوية الطبيعية recount normalization، وأضفنا خريطة ألوان طبيعية (gamma=0.3) إضافةً إلى تظليل فاتح. تمرين نريد الآن قياس البعد الكسري لمجموعة ماندلبورت Mandelbrot باستخدام بُعد Minkowski–Bouligand dimension، ولفعل ذلك نحتاج إلى حسابات لمربع بأبعاد متناقص، كما هو موضح في لقطة الشاشة التالية، وكما تتوقع لا يمكننا استخدام تعليمات بايثون العادية لأنها ستكون بطيئة جدًا. الهدف من التمرين هو كتابة دالة باستخدام NumPy تأخذ مصفوفة عشرية ثنائية الأبعاد وتُرجع البعد. ينبغي توحيد القيم الموجودة في المصفوفة بين 0 و1. المتجهات المكانية تشير المتجهات المكانية إلى الحالة التي تشارك فيها العناصر نفس الحساب ولكنها تتفاعل مع مجموعة فرعية فقط من العناصر الأخرى، وكانت هذه الحالة موجودةٌ فعلًا في لعبة الحياة، ولكن في بعض المواقف هناك صعوبة إضافية لأن المجموعة الفرعية ديناميكية وتحتاج إلى التحديث في كل تكرار؛ فمثلًا في أنظمة الجسيمات تتفاعل الجسيمات غالبًا مع الجيران المحليين، وينطبق هذا أيضًا على حالة "boids" التي تحاكي سلوكيات السرب flocking. برنامج Boids الطيور المتدفقة هي مثال على التنظيم الذاتي في علم الأحياء. صورة Christoffer A Rasmussen، 2012. يُعد Boids برنامج حياة اصطناعية طُور بواسطة كريج رينولدز Craig Reynolds في عام 1986 وهو يحاكي سلوك سرب الطيور. اسم Boid هو اختصار إلى " bird-oid object" والذي يشير إلى كائن يمثل الطيور. كما هو الحال مع معظم برامج محاكاة الحياة الاصطناعية، فإن Boids هو مثال على السلوك الناشئ، أي أن تعقيد Boids ينشأ من تفاعل العوامل الفردية والالتزام بمجموعة من القواعد البسيطة. القواعد المطبقة في عالم Boids هي كما يلي: الانفصال separation: إرشاد لتجنب ازدحام رفاق السرب المحليين. المحاذاة alignment: إرشاد نحو متوسط اتجاه رفاق السرب المحليين. التماسك cohesion: الإرشاد للتحرك نحو متوسط الموضع (مركز الكتلة) للرفاق الحاليين. تنفيذ الشيفرة باستخدام تعليمات بايثون العادية نظرًا لأن كل boid هو كيان مستقل له العديد من الخاصيات، مثل الموضع والسرعة فمن الطبيعي أن نبدأ بكتابة الصنف Boid: import math import random from vec2 import vec2 class Boid: def __init__(self, x=0, y=0): self.position = vec2(x, y) angle = random.uniform(0, 2*math.pi) self.velocity = vec2(math.cos(angle), math.sin(angle)) self.acceleration = vec2(0, 0) الكائن vec2 هو صنف بسيط جدًا يتعامل مع جميع عمليات المتجه الشائعة بمكونين، وهذا سيوفر علينا بعض الكتابة في صنف Boid الرئيسي. يوجد بعض حزم المتجهات في فهرس حزمة Python، ولكن استخدامها سيكون مبالغة في مثل هذا المثال البسيط. يُعد السرب حالةً صعبة لبايثون العادية لأنه يتفاعل مع الجيران المحليين، ومع ذلك ونظرًا لأن أسراب الطيور تتحرك فإن العثور على مثل هؤلاء الجيران المحليين يتطلب حسابات في كل مرة تقطع المسافة إلى سرب طيور آخر لفرز تلك الطيور الموجودة في نصف قطر تفاعل معين، وبالتالي تكون الطريقة النموذجية لكتابة القواعد الثلاثة: def separation(self, boids): count = 0 for other in boids: d = (self.position - other.position).length() if 0 < d < desired_separation: count += 1 ... if count > 0: ... def alignment(self, boids): ... def cohesion(self, boids): … لإكمال الصورة، يمكننا أيضًا إنشاء كائن Flock: class Flock: def __init__(self, count=150): self.boids = [] for i in range(count): boid = Boid() self.boids.append(boid) def run(self): for boid in self.boids: boid.run(self.boids) يمكن أن نحصل باستخدام هذا النهج على ما يصل إلى 50 سربًا حتى يصبح وقت الحساب بطيئًا جدًا للحصول على رسم متحرك سلس. سنحصل بالاعتماد على مكتبة Numpy على أداء أفضل بكثير، ولكن قبل الانتقال إلى استخدام هذه المكتبة سنشير إلى المشكلة الرئيسية الموجودة في شيفرات بايثون السابقة؛ فكما تلاحظ يوجد في الشيفرات الكثير من التكرار، أي بتعبير أدق نحن لا نستغل حقيقة أن المسافة الإقليدية انعكاسية أي: |x − y| = |y − x|. والحقيقة أن كل دالة تحسب "n2" مسافة بينما حساب مسافة "n2/2" كافي، إذا خُزّنت تخزينًا مؤقتًا بطريقة صحيحة. وأيضًا تعيد كل قاعدة حساب كل مسافة دون تخزين نتيجة الدوال الأخرى تخزينًا مؤقتًا، وفي النهاية، نكون قد حسبنا "3n2" مسافة بدلًا من "n2/2". كتابة الشيفرة باستخدام مكتبة NumPy تتخذ كتابة الشيفرة باستخدام مكتبة Numpy نهجًا مختلفًا، سنجمع كل كائنات أسراب الطيور boids في مصفوفة الموضع ومصفوفة السرعة: n = 500 velocity = np.zeros((n, 2), dtype=np.float32) position = np.zeros((n, 2), dtype=np.float32) الخطوة الأولى هي حساب الجيران لكل كائنات boids، ولذلك نحتاج إلى حساب المسافات المزدوجة: dx = np.subtract.outer(position[:, 0], position[:, 0]) dy = np.subtract.outer(position[:, 1], position[:, 1]) distance = np.hypot(dx, dy) كان بإمكاننا استخدام دالة scipy cdist لكننا سنحتاج إلى المصفوفتين dx و dy لاحقًا، وبمجرد حسابها سيكون من الأسرع استخدام تابع hypot. لاحظ شكل المسافة (n,n) وكل خط يرتبط بكائن boid واحد، أي أن كل سطر يعطي المسافة إلى جميع كائنات boids الأخرى. يمكننا الآن حساب الجيران من هذه المسافات وفقًا للقواعد الثلاثة، مستفيدين من إمكانية دمجها معًا. يمكننا حساب قناع للمسافة والذي يكون موجبًا تمامًا (أي ليس له ارتباط ذاتي)، ثم جداء هذا القناع مع أقنعة المسافات الأخرى. mask_0 = (distance > 0) mask_1 = (distance < 25) mask_2 = (distance < 50) mask_1 *= mask_0 mask_2 *= mask_0 mask_3 = mask_2 ملاحظة: إذا افترضنا أن أسراب الطيور boids لا يمكن أن تحتل نفس الموضع، فكيف يمكن حساب mask_0 بفعالية أكثر؟ ثم نحسب عدد الجيران داخل نصف القطر المحدد ونتأكد من أنه 1 على الأقل لتجنب القسمة على الصفر. mask_1_count = np.maximum(mask_1.sum(axis=1), 1) mask_2_count = np.maximum(mask_2.sum(axis=1), 1) mask_3_count = mask_2_count الآن، سنكتب القواعد الثلاث: المحاذاة # حساب متوسط السرعة للجيران المحليين target = np.dot(mask, velocity)/count.reshape(n, 1) # توحيد النتيجة norm = np.sqrt((target*target).sum(axis=1)).reshape(n, 1) target *= np.divide(target, norm, out=target, where=norm != 0) # المحاذاة بسرعة ثابتة target *= max_velocity # الناتج steering حساب مسار الإرشاد alignment = target - velocity التماسك # حساب مركز الجاذبية للجيران المحليين center = np.dot(mask, position)/count.reshape(n, 1) # حساب الاتجاه إلى المركز target = center - position # توحيد النتيجة norm = np.sqrt((target*target).sum(axis=1)).reshape(n, 1) target *= np.divide(target, norm, out=target, where=norm != 0) # (max_velocity) التماسك بسرعة ثابتة target *= max_velocity # الناتج steering حساب مسار الإرشاد cohesion = target - velocity الانفصال # حساب قوة التنافر من الجيران المحليين repulsion = np.dstack((dx, dy)) # تتناسب القوة عكسيًا مع المسافة repulsion = np.divide(repulsion, distance.reshape(n, n, 1)**2, out=repulsion, where=distance.reshape(n, n, 1) != 0) # حساب الاتجاه بعيدًا عن الآخرين target = (repulsion*mask.reshape(n, n, 1)).sum(axis=1)/count.reshape(n, 1) # توحيد النتيجة norm = np.sqrt((target*target).sum(axis=1)).reshape(n, 1) target *= np.divide(target, norm, out=target, where=norm != 0) # (max_velocity) الانفصال بسرعة ثابتة target *= max_velocity # الناتج steering حساب مسار الإرشاد separation = target - velocity يجب أن تكون مسارات الإرشاد الثلاث الناتجة (الانفصال والمحاذاة والتماسك) محدودةً من حيث الحجم. وسنترك لك الدمج بين هذه القواعد. التحديث الناتج للسرعة والموضع: acceleration = 1.5 * separation + alignment + cohesion velocity += acceleration position += velocity أخيرًا نتخيل النتيجة باستخدام مخطط مبعثر scatter موجه مخصص. Boids هو برنامج حياة اصطناعية، طُوِّر بواسطة Craig Reynolds في عام 1986، والذي يحاكي سلوك الطيور المتدفقة. تمرين نحن الآن جاهزون لعرض أسراب الطيور boids، وأسهل طريقة لذلك هي استخدام دالة الرسوم المتحركة matplotlib ومخطط مبعثر. لسوء الحظ لا يمكننا توجيه التبعثر بطريقة فردية ونحتاج إلى صنع كائنات خاصة باستخدام توابع PathCollection العائد للدالة matplotlib. يمكن تعريف مسار المثلث البسيط كما يلي: v= np.array([(-0.25, -0.25), ( 0.00, 0.50), ( 0.25, -0.25), ( 0.00, 0.00)]) c = np.array([Path.MOVETO, Path.LINETO, Path.LINETO, Path.CLOSEPOLY]) يمكن تكرار هذا المسار عدة مرات داخل المصفوفة ويمكن جعل كل مثلث مستقلاً. n = 500 vertices = np.tile(v.reshape(-1), n).reshape(n, len(v), 2) codes = np.tile(c.reshape(-1), n) لدينا الآن مصفوفة (n, 4, 2) ومجموعة (n, 4) للشيفرات التي تمثل عدد n من boids. نحن مهتمون بمعالجة مصفوفة vertices لتعكس الترجمة والقياس والدوران لكل كائنات boids. كيف يمكن أن تكتب وظائف الترجمة والقياس والتدوير؟ الخلاصة درسنا من خلال هذه الأمثلة الثلاثة ثلاثة أشكال لاستخدام المتجهات في الشيفرات: المتجهات الموحدة، إذ تشترك العناصر في نفس الحساب دون قيد أو شرط وللمدة نفسها. المتجهات الزمانية، إذ تشترك العناصر في نفس الحساب ولكنها تتطلب عددًا مختلفًا من التكرارات. المتجهات المكانية، إذ تشترك العناصر في نفس الحساب ولكن على الوسطاء المكانية الديناميكية. ومن المحتمل أن يكون هناك العديد من أشكال متجهات الشيفرات. كما أوضحنا سابقًا يعد هذا النوع من المتجهات من أبسط الأنواع على الرغم من أننا رأينا أنه من الصعب حقًا تنفيذه ويتطلب بعض الخبرة أو بعض المساعدة أو كليهما. ترجمة -وبتصرف- للفصل Code vectorization من كتاب From Python to Numpy لصاحبه Nicolas P. Rougier. اقرأ أيضًا المقال السابق: مفاهيم متقدمة حول مكتبة NumPy في بايثون أهم 8 مكتبات بلغة البايثون تستخدم في المشاريع الصغيرة

تُعد مكتبة NumPy إحدى مكتبات لغة بايثون Python وتُستخدم للتعامل مع المصفوفات، وتهدف إلى توفير كائن مصفوفة أسرع بما يصل إلى 50 مرة من قوائم بايثون التقليدية. يطلق على كائن المصفوفة في مكتبة NumPy اسم ndarray، ويوفر العديد من الوظائف الداعمة التي تجعل التعامل مع المصفوفات أمرًا سهلًا جدًا. تعتمد مكتبة NumPy على العمليات على المتجهات vectorization، فإذا كنت معتادًا على التعامل مع لغة بايثون، فهذه هي الصعوبة الرئيسية التي ستواجهها لأنك ستحتاج إلى تغيير طريقة تفكيرك، كما ستتغير نوع العناصر التي ستستخدمها. أكثر العناصر شيوعًا في مكتبة NumPy هي المتجهات vectors والمصفوفات arrays والعروض views والدوال العمومية ufuns. هذا المقال جزء من سلسلة متقدمة حول NumPy وإليك كامل روابط السلسلة: مفاهيم متقدمة حول مكتبة NumPy في بايثون استخدام المتجهات vectorization في لغة بايثون مع مكتبة NumPy الاعتماد على المتجهات في حل المشاكل باستعمال مكتبة NumPy في بايثون تخصيص أسلوب المتجهات Custom vectorization عبر استعمال مكتبة NumPy ما بعد مكتبة NumPy في بايثون مثال بسيط سنأخذ مثالًا بسيطًا وهو المشي العشوائي random walk، إذ يتوجب علينا في الطريقة التقليدية المعتمدة على البرمجة كائنية التوجه OOP تعريف صنف class يأخذ الاسم RandomWalker، ثم إنشاء التابع walk لإرجاع الموقع الحالي بعد كل خطوة عشوائية. يمكن تنفيذ هذا المثال باستخدام التعليمات البرمجية التالية: class RandomWalker: def __init__(self): self.position = 0 def walk(self, n): self.position = 0 for i in range(n): yield self.position self.position += 2*random.randint(0, 1) - 1 walker = RandomWalker() walk = [position for position in walker.walk(1000)] وتعطي تجربة قياس السرعة النتيجة التالية: >>> from tools import timeit >>> walker = RandomWalker() >>> timeit("[position for position in walker.walk(n=10000)]", globals()) 10 loops, best of 3: 15.7 msec per loop النهج الإجرائي Procedural approach يمكننا لحل مثل هذه المشكلة البسيطة حفظ تعريف الصنف والتركيز على دالة المشي walk التي تحسب المواقع المتتالية بعد كل خطوة عشوائية. def random_walk(n): position = 0 walk = [position] for i in range(n): position += 2*random.randint(0, 1)-1 walk.append(position) return walk walk = random_walk(1000) توفر هذه الطريقة من استهلاك وحدة المعالجة المركزية CPU ولكن ليس لحدٍ كبير لأنها تشبه الطريقة السابقة المعتمدة على البرمجة كائنية التوجه. >>> from tools import timeit >>> timeit("random_walk(n=10000)", globals()) 10 loops, best of 3: 15.6 msec per loop النهج المعتمد على المتجهات Vectorized approach يمكننا باستخدام وحدة itertools تحسين الأداء بطريقة أفضل، إذ تقدم هذه الوحدة مجموعةً من الدوال لإنشاء تكرارات لحلقات فعالة. نلاحظ أن مثال المشي العشوائي ما هو إلا تراكم خطوات، لذلك سنعيد كتابة الدالة بإنشاء جميع الخطوات أولًا، ثم تجميعها بدون استخدام أي حلقة، كما توضح التعليمات التالية: def random_walk_faster(n=1000): from itertools import accumulate # Only available from Python 3.6 steps = random.choices([-1,+1], k=n) return [0]+list(accumulate(steps)) walk = random_walk_faster(1000) اعتمدنا في هذا المثال على المتجهات في الدالة، فبدلًا من استخدام التكرار لاختيار الخطوات وإضافتها إلى الموضع الحالي، أنشأنا أولًا جميع الخطوات دفعةً واحدةً، ثم استخدمنا دالة التجميع accumulate لحساب جميع المواضع. تخلصنا من الحلقة وبالتالي أصبح تنفيذ التعليمات أسرع: >>> from tools import timeit >>> timeit("random_walk_faster(n=10000)", globals()) 10 loops, best of 3: 2.21 msec per loop لاحظ أننا وفرنا 85% من وقت الحساب مقارنةً بالإصدار السابق، لكن الإيجابية التي تميز الإصدار الجديد هي تبسيط العمليات المعقدة على المتجهات. سنحتاج إلى استبدال أدوات itertools بأدوات numpy: def random_walk_fastest(n=1000): # No 's' in numpy choice (Python offers choice & choices) steps = np.random.choice([-1,+1], n) return np.cumsum(steps) walk = random_walk_fastest(1000) هذه الطريقة سهلة وتوفر الكثير من الوقت كما توضح النتائج التالية: >>> from tools import timeit >>> timeit("random_walk_fastest(n=10000)", globals()) 1000 loops, best of 3: 14 usec per loop يركز هذا الدليل يركز على الاعتماد على المتجهات سواءٌ كان ذلك على مستوى المشكلة أو الشيفرات البرمجية، وسنوضح أهمية هذا الاختلاف قبل الانتقال إلى دراسة المتجهات المخصصة. سهولة القراءة مقابل السرعة قبل الانتقال إلى القسم التالي، سنلقي نظرةً على مشكلة محتملة قد تواجهها بعد أن تصبح على دراية بمكتبة NumPy؛ فمكتبة NumPy مكتبة قوية جدًا ويمكنك أن تصنع العجائب بها، ولكن هذا على حساب سهولة القراءة، لذلك إذا أهملت إضافة التعليقات أثناء كتابة التعليمات البرمجية، فلن تتمكن من معرفة وظيفة بعض التعليمات بعد بضع أسابيع أو ربما أيام. على سبيل المثال هل يمكنك معرفة ما هي نتيجة تنفيذ الدالتين التاليتين؟ ربما ستتمكن من فهم الدالة الأولى ولكن من غير المرجح أن تفهم الثانية: def function_1(seq, sub): return [i for i in range(len(seq) - len(sub)) if seq[i:i+len(sub)] == sub] def function_2(seq, sub): target = np.dot(sub, sub) candidates = np.where(np.correlate(seq, sub, mode='valid') == target)[0] check = candidates[:, np.newaxis] + np.arange(len(sub)) mask = np.all((np.take(seq, check) == sub), axis=-1) return candidates[mask] استخدمت الدالة الثانية مكتبة NumPy وهي أسرع بعشر مرات، ولكن المشكلة هي بصعوبة قراءة وفهم التعليمات البرمجية التي تستخدمها. تشريح المصفوفة كما أوضحنا سابقًا ينبغي أن تكون لديك خبرةً أساسيةً في التعامل مع مكتبة NumPy قبل المتابعة في هذا الدليل، وإذا لم يكن لديك هذه الخبرة فمن الأفضل أن تتعلم عنها قبل العودة إلى هنا، وبناءً على ذلك سنبدأ هذا الدليل بتذكير سريع على الهيكلية الأساسية للمصفوفات المعقدة، خاصةً بما يتعلق بتخطيط الذاكرة وعرضها ونسخ البيانات وتحديد نوعها، وهذه مفاهيم أساسية يجب أن تفهمها إذا كنت تريد أن تستفيد من فلسفة مكتبة numpy. لنفكر في مثال بسيط نريد من خلاله مسح جميع القيم في مصفوفة من النوع np.float32. كيف يمكن كتابة تعليمات تنفذ المثال السابق بسرعة؟ الصيغة أدناه واضحة إلى حد ما (على الأقل لأولئك الذين هم على دراية بمكتبة NumPy) ولكن السؤال أعلاه يطلب العثور على أسرع عملية. >>> Z = np.ones(4*1000000, np.float32) >>> Z[...] = 0 إذا دققت بنوع عناصر المصفوفة dtype وحجمها، يمكنك ملاحظة أنه يمكن تحويل هذه المصفوفة إلى العديد من أنواع البيانات "المتوافقة" الأخرى. نقصد هنا أنه يمكن تقسيم Z.size * Z.itemsize إلى أنواع عناصر dtype جديدة. >>> timeit("Z.view(np.float16)[...] = 0", globals()) 100 loops, best of 3: 2.72 msec per loop >>> timeit("Z.view(np.int16)[...] = 0", globals()) 100 loops, best of 3: 2.77 msec per loop >>> timeit("Z.view(np.int32)[...] = 0", globals()) 100 loops, best of 3: 1.29 msec per loop >>> timeit("Z.view(np.float32)[...] = 0", globals()) 100 loops, best of 3: 1.33 msec per loop >>> timeit("Z.view(np.int64)[...] = 0", globals()) 100 loops, best of 3: 874 usec per loop >>> timeit("Z.view(np.float64)[...] = 0", globals()) 100 loops, best of 3: 865 usec per loop >>> timeit("Z.view(np.complex128)[...] = 0", globals()) 100 loops, best of 3: 841 usec per loop >>> timeit("Z.view(np.int8)[...] = 0", globals()) 100 loops, best of 3: 630 usec per loop ولكن في الواقع فإن طريقة تصفية جميع القيم ليست هي الأسرع، فقد زاد معدل السرعة بنسبة 25% من خلال تحويل بعض عناصر المصفوفة إلى نوع البيانات np.float64، بينما زاد معدل السرعة بنسبة 50% من خلال عرض المصفوفة واستخدام نوع البيانات np.int8. تتعلق أسباب التسريع بهيكلية عمل مكتبة numpy. تخطيط الذاكرة Memory layout يُعرّف الصنف ndarray وفقًا لتوثيق Numpy كما يلي: المصفوفة هي كتلة مجاورة من الذاكرة يمكن الوصول إلى عناصرها باستخدام مخطط الفهرسة indexing scheme، الذي يحدد شكل ونوع البيانات، وهذا هو المطلوب لتعريف مصفوفة جديدة. Z = np.arange(9).reshape(3,3).astype(np.int16) حجم عناصر المصفوفة Z في التعليمة السابقة هو 2 بايت (int16)، وعدد أبعاد (len(Z.shape هو 2. >>> print(Z.itemsize) 2 >>> print(Z.shape) (3, 3) >>> print(Z.ndim) 2 نظرًا لأن Z ليست طريقة عرض، يمكننا استنتاج خطوات المصفوفة التي تحدد عدد البايتات التي يجب أن تخطوها لاجتياز كل بعد من أبعاد المصفوفة. >>> strides = Z.shape[1]*Z.itemsize, Z.itemsize >>> print(strides) (6, 2) >>> print(Z.strides) (6, 2) لنتمكن من الوصول إلى عنصر معين (مصمم بواسطة فهرس index tuple)، أي كيفية حساب إزاحة البداية والنهاية: offset_start = 0 for i in range(ndim): offset_start += strides[i]*index[i] offset_end = offset_start + Z.itemsize سنستخدم الآن طريقة التحويل tobytes للتأكد من المعلومة السابقة: >>> Z = np.arange(9).reshape(3,3).astype(np.int16) >>> index = 1,1 >>> print(Z[index].tobytes()) b'\x04\x00' >>> offset = 0 >>> for i in range(Z.ndim): ... offset + = Z.strides[i]*index[i] >>> print(Z.tobytes()[offset_start:offset_end] b'\x04\x00' يمكن تمثيل هذه المصفوفة بطرق (تخطيطات) مختلفة: تخطيط العنصر item layout تخطيط العنصر المسطّح Flattened item layout تخطيط الذاكرة Memory layout (C order, big endian)‎ إذا أخذنا الآن شريحة من Z، فإن النتيجة هي عرض المصفوفة الأساسية Z: V = Z[::2,::2] يُحدّد هذا العرض باستخدام نوع dtype وشكل وخطوات، إذ أصبح من غير الممكن استنتاج الخطوات strides من النوع والشكل فقط: تخطيط العنصر item layout تخطيط العنصر المسطّح Flattened item layout تخطيط الذاكرة Memory layout (C order, big endian)‎ العروض والنسخ تُعد العروض Views والنسخ copies مفاهيمًا مهمة لتحسين حساباتك الرقمية. بالرغم من تعاملنا مع هذه المفاهيم في الفقرات السابقة إلا أن الموضوع أكثر تعقيدًا. الوصول المباشر وغير المباشر أولًا علينا التمييز بين نوعي الفهرسة indexing و fancy indexing؛ إذ يعيد النوع الأول عرضًا view دائمًا، بينما يعيد النوع الثاني نسخةً copy. يؤدي تعديل العرض في الحالة الأولى إلى تعديل المصفوفة الأساسية، بينما لن يحدث ذلك في الحالة الثانية: >>> Z = np.zeros(9) >>> Z_view = Z[:3] >>> Z_view[...] = 1 >>> print(Z) [ 1. 1. 1. 0. 0. 0. 0. 0. 0.] >>> Z = np.zeros(9) >>> Z_copy = Z[[0,1,2]] >>> Z_copy[...] = 1 >>> print(Z) [ 0. 0. 0. 0. 0. 0. 0. 0. 0.] إذا كنت تريد استخدام نوع الفهرسة fancy indexing فمن الأفضل الاحتفاظ بنسخة الفهرس fancy index ثم التعامل معه: >>> Z = np.zeros(9) >>> index = [0,1,2] >>> Z[index] = 1 >>> print(Z) [ 1. 1. 1. 0. 0. 0. 0. 0. 0.] إذا لم تكن متأكدًا مما إذا كانت نتيجة الفهرسة هي عرض أو نسخة، فيمكنك التحقق باستخدام الأداة base؛ فإذا كانت النتيجة "None"، فسيكون النوع هو نسخة: >>> Z = np.random.uniform(0,1,(5,5)) >>> Z1 = Z[:3,:] >>> Z2 = Z[[0,1,2], :] >>> print(np.allclose(Z1,Z2)) True >>> print(Z1.base is Z) True >>> print(Z2.base is Z) False >>> print(Z2.base is None) True لاحظ أن بعض وظائف numpy تُعيد عرضًا view، مثل ravel عندما يكون ذلك ممكنًا، بينما يعيد البعض الآخر نسخةً copy، مثل flatten: >>> Z = np.zeros((5,5)) >>> Z.ravel().base is Z True >>> Z[::2,::2].ravel().base is Z False >>> Z.flatten().base is Z False نسخة مؤقتة Temporary copy يمكن إنشاء النسخ بطريقة صريحة كما في الأمثلة السابقة، ولكن الطريقة الأكثر استخدمًا هي الإنشاء الضمني للنسخ الوسيطة. يوضح المثال التالي إجراء العمليات الحسابية باستخدام المصفوفات: >>> X = np.ones(10, dtype=np.int) >>> Y = np.ones(10, dtype=np.int) >>> A = 2*X + 2*Y أنشأنا في هذا المثال ثلاث مصفوفات وسيطة: مصفوفة للاحتفاظ بنتيجة X*‏2 والثانية للاحتفاظ بنتيجة Y*‏2، والأخيرة لإيجاد نتيجة Y*‏2 + X*‏2، وتكون المصفوفات في هذه الحالة صغيرة الحجم ولا تحدث هذه الطريقة أي فرق واضح؛ ولكن في حال كانت المصفوفات كبيرة الحجم فعليك أن تكون حذرًا مع هذه التعبيرات وتفكر بطريقة مختلفة لإيجاد الحل، فعلى سبيل المثال إذا كانت النتيجة النهائية هي فقط المهمة ولم تكن بحاجة إلى X أو Y، فسيكون الحل البديل: >>> X = np.ones(10, dtype=np.int) >>> Y = np.ones(10, dtype=np.int) >>> np.multiply(X, 2, out=X) >>> np.multiply(Y, 2, out=Y) >>> np.add(X, Y, out=X) استغنينا باستخدام هذا الحل البديل عن فكرة إنشاء مصفوفات مؤقتة، لكن المشكلة أن هناك بعض الحالات التي تتطلب إنشاء هذه النسخ المؤقتة وهذا يؤثر على الأداء كما هو موضح في المثال التالي: >>> X = np.ones(1000000000, dtype=np.int) >>> Y = np.ones(1000000000, dtype=np.int) >>> timeit("X = X + 2.0*Y", globals()) 100 loops, best of 3: 3.61 ms per loop >>> timeit("X = X + 2*Y", globals()) 100 loops, best of 3: 3.47 ms per loop >>> timeit("X += 2*Y", globals()) 100 loops, best of 3: 2.79 ms per loop >>> timeit("np.add(X, Y, out=X); np.add(X, Y, out=X)", globals()) 1000 loops, best of 3: 1.57 ms per loop الخلاصة في الختام سندرس المثال التالي: لنفرض لدينا متجهين Z1 و Z2، ونرغب في معرفة ما إذا كان Z2 هو عرض view للمتجه Z1، وإذا كانت الإجابة نعم فما هو هذا العرض view؟ >>> Z1 = np.arange(10) >>> Z2 = Z1[1:-1:2] أولاً سنحدد نوع الفهرسة وهل Z1 هي أساس Z2: >>> print(Z2.base is Z1) True حتى الآن تأكدنا أن Z2 هي عرض view للمتغير Z1، مما يعني أنه يمكن التعبير عن Z2على النحو التالي: [Z1[start:stop:step. وتكمن الصعوبة في التعرف على البداية والتوقف والخطوة؛ إذ يمكننا بالنسبة للخطوة استخدام خاصية strides لأي مصفوفة تعطي عدد البايتات للانتقال من عنصر إلى آخر في كل بُعد. في حالتنا ولأن كلا المصفوفتين أحادي البعد، يمكننا مقارنة الخطوة الأولى فقط. >>> step = Z2.strides[0] // Z1.strides[0] >>> print(step) 2 الصعوبة التالية هي العثور على مؤشرات البداية والتوقف، إذ يمكننا فعل ذلك بالاستفادة من طريقة byte_bounds التي تُرجع مؤشرًا إلى نقاط نهاية المصفوفة. >>> offset_start = np.byte_bounds(Z2)[0] - np.byte_bounds(Z1)[0] >>> print(offset_start) # bytes 8 >>> offset_stop = np.byte_bounds(Z2)[-1] - np.byte_bounds(Z1)[-1] >>> print(offset_stop) # bytes -16 يعد تحويل هذه الإزاحات offset إلى مؤشرات أمرًا سهلاً باستخدام حجم العناصر مع الأخذ بالحسبان أن offset_stop سلبي (حد نهاية Z2 أصغر منطقيًا من حد النهاية لمصفوفة Z1)، لذلك نحتاج إلى إضافة حجم عناصر Z1 للحصول على فهرس النهاية الصحيحة: >>> start = offset_start // Z1.itemsize >>> stop = Z1.size + offset_stop // Z1.itemsize >>> print(start, stop, step) 1, 8, 2 وعند اختبار النتائج نحصل على ما يلي: >>> print(np.allclose(Z1[start:stop:step], Z2)) True يمكنك التمرين من خلال تحسين هذا التنفيذ الأول والبسيط جدًا من خلال مراعاة: خطوات سلبية. مصفوفة متعددة الأبعاد. ترجمة -وبتصرف- للفصلين Introduction و Anatomy of an array من كتاب From Python to Numpy لصاحبه Nicolas P. Rougier. اقرأ أيضًا المرجع الشامل إلى تعلم لغة بايثون أهم 8 مكتبات بلغة البايثون تستخدم في المشاريع الصغيرة